我有一个从截断的正态分布生成一个值的函数,并带有一个while循环,该循环可确保舍弃位于截断之外的任何生成的值,并将其替换为另一个生成值,直到它位于该范围内。
def gen_truncated(minimum, maximum, ave, sigma):
# min=0.9, max=1,
x = 0.
while x < minimum or x > maximum:
x = np.random.normal(0,1)*sigma+ave
return x
我如何以如下方式矢量化此函数:x
现在是由许多x
值组成的数组,其生成方式是始终有while循环确保只要条件满足,就可以重新生成数组元素达到x < minimum
和x > maximum
?是否有一种向量化的方式将x
的每个元素与一个数字进行比较,即minimum
或maximum
?
编辑:如果我还有更多需要满足的约束怎么办?最终,我希望对通过多个约束生成的4x4矩阵的生成进行矢量化处理,gen_truncated()
中的约束只是众多约束之一。我有一个gen_sigma()
,它首先会生成3个值lambda1, lambda2, lambda3
,现在lambda3
再次需要满足几个条件,而lambda1
和lambda2
的值会重新绘制。一旦它们正确,就将所有三个值馈入get_tau()
以生成3个值。同样,这些tau值需要满足更多约束,否则它们将被丢弃并再次生成,直到正确为止。最终,它们形成一个称为sigma_gen
的4x4矩阵,将其左右乘以create_rotation()
至gen_channel
,以创建单个4x4矩阵channel
。
import numpy as np from numpy.linalg import norm def gen_sigma(minimum, maximum, ave, sigma): lambda1 = gen_truncated(minimum, maximum, ave, sigma) lambda2 = gen_truncated(minimum, maximum, ave, sigma) lambda3 = gen_truncated(minimum, maximum, ave, sigma) while 1+lambda3 < abs(lambda1+lambda2) or 1-lambda3 < abs(lambda2-lambda1): lambda3 = gen_truncated(minimum, maximum, ave, sigma) tau = get_tau(lambda1, lambda2, lambda3) lambdas = [lambda1, lambda2, lambda3] while (norm(tau)**2 > 1-sum([x**2 for x in [lambda1, lambda2, lambda3]]) + 2*lambda1*lambda2*lambda3) or (z_eta(tau, lambdas) < 0): tau = get_tau(lambda1, lambda2, lambda3) sigma_gen = np.array([[ 1, 0, 0, 0], [tau[0], lambda1, 0, 0], [tau[1], 0, lambda2, 0], [tau[2], 0, 0, lambda3]]) return sigma_gen def get_tau(einval1, einval2, einval3): max_tau1 = 1 - abs(einval1) max_tau2 = 1 - abs(einval2) max_tau3 = 1 - abs(einval3) tau1 = max_tau1*(2*np.random.uniform(0,1)-1) tau2 = max_tau2*(2*np.random.uniform(0,1)-1) tau3 = max_tau3*(2*np.random.uniform(0,1)-1) return [tau1, tau2, tau3] def z_eta(t: np.ndarray, l: np.ndarray): condition = (norm(t)**4 - 2*norm(t)**2 - 2*sum([(l[i]**2)*(2*(t[i]**2-norm(t)**2)) for i in range(3)])+ q(l)) return condition def q(e: np.ndarray): # e are the eigenvalues return (1+e[0]+e[1]+e[2])*(1+e[0]-e[1]-e[2])*(1-e[0]+e[1]-e[2])*(1-e[0]-e[1]+e[2]) def create_rotation(angles: np.ndarray) -> np.ndarray: "random rotation in PL form" # input np.random.normal(0,1,3)*0.06 rotation = np.eye(4, dtype=complex) left = np.array([[ np.cos(angles[0]), np.sin(angles[0]), 0], [-np.sin(angles[0]), np.cos(angles[0]), 0], [ 0, 0, 1]]) mid = np.array([[1, 0, 0], [0, np.cos(angles[1]), np.sin(angles[1])], [0, -np.sin(angles[1]), np.cos(angles[1])]]) right = np.array([[ np.cos(angles[2]), np.sin(angles[2]), 0], [-np.sin(angles[2]), np.cos(angles[2]), 0], [ 0, 0, 1]]) rotation[1:4,1:4] = left@mid@right return rotation def gen_channel(r1, r2, ave, sigma): rand1 = np.random.normal(0,1,3) rand2 = np.random.normal(0,1,3) channel = create_rotation(rand1*r1)@gen_sigma(0.9, 1, ave, sigma)@\ create_rotation(rand2*r2) return channel
频道的示例运行
gen_channel(0.05, 0.05, 0.98, 0.15)
例如]
Out[140]: array([[ 1. +0.j, 0. +0.j, 0. +0.j, 0. +0.j], [-0.05828008+0.j, 0.91805971+0.j, 0.14291751+0.j, -0.00946994+0.j], [-0.00509449+0.j, -0.14170308+0.j, 0.90034613+0.j, -0.11548884+0.j], [ 0.0467522 +0.j, -0.00851749+0.j, 0.11450963+0.j, 0.90259637+0.j]])
现在,如果我想创建100个4x4矩阵,我将不得不使用列表理解,即
np.array([gen_channel(0.05, 0.05, 0.98, 0.15) for i in range(100)])
将进行所有约束比较,并一一创建4x4矩阵。现在,我最初的问题是由以下事实引起的:我想对它们进行矢量化处理,而不是一次比较一个值,而只是使用numpy广播生成一个值数组并检查约束,以便获得[C0生成100个这样的4x4矩阵,而无需列表理解。列表理解方法包含重复使用生成单个随机数的方法,这会导致其运行速度出现瓶颈。我要做的只是生成随机数数组,进行检查,然后生成4x4通道数组,以减少瓶颈。
我有一个从截断的正态分布生成一个值的函数,并带有一个while循环,该循环可确保舍弃位于截断之外的任何生成的值,并将其替换为...
您可以从原始分布中抽取大量样本,然后确定哪些条目位于正确的范围内,然后从中抽取:
gen_channel