时间戳

首页 > 解决方案 > 是否可以沿轴为 ndarray 执行 linalg.multi_dot ？

numpy - 是否可以沿轴为 ndarray 执行 linalg.multi_dot ？

问题描述

首先，我有一组 12 (2x2) 个矩阵。

II = np.identity(2, dtype=complex)
X = np.array([[0, 1], [1, 0]], dtype=complex)
Y = np.array([[0, -1j], [1j, 0]], dtype=complex)
Z = np.array([[1, 0], [0, -1]], dtype=complex)
PPP = (-II + 1j*X + 1j*Y + 1j*Z)/2
PPM = (-II + 1j*X + 1j*Y - 1j*Z)/2
PMM = (-II + 1j*X - 1j*Y - 1j*Z)/2
MMM = (-II - 1j*X - 1j*Y - 1j*Z)/2
MMP = (-II - 1j*X - 1j*Y + 1j*Z)/2
MPP = (-II - 1j*X + 1j*Y + 1j*Z)/2
PMP = (-II + 1j*X - 1j*Y + 1j*Z)/2
MPM = (-II - 1j*X + 1j*Y - 1j*Z)/2

目前我有一个函数operator_groups，它从这个组中为每个j循环绘制一个随机矩阵，并将它附加到一个列表中sequence。然后使用在所有单独循环之间绘制的随机矩阵j进行一些计算，与我们在这里的讨论无关。在j循环结束时，列表元素的序列sequence被反转，然后linalg.multi_dot执行，然后取其厄米特共轭（因此.conj().T）

def operator_groups():
    return random.choice([II, X, Y, Z, PPP, PPM, PMM, MMM, MMP, MPP, PMP, MPM])


for i in range(1, sample_size+1, 1):
    sequence = []
    for j in range(1, some_number, 1):
        noise = operator_groups()
        """some matrix calculations here"""
        sequence.append(noise)
    sequence_inverse = np.linalg.multi_dot(sequence[::-1]).conj().T

现在我希望ij通过在一个大矩阵中执行循环来对循环进行矢量化。现在noise是从组中随机采样的 N 个矩阵（而不仅仅是 1 个矩阵）的 ndarray，每个矩阵代表的迭代j，但并行化。代码现在看起来像这样。

def operator_groups(sample_size):
    return random.sample([II, X, Y, Z], sample_size)


sequence = []
for j in range(1, some_number, 1):
    noise = operator_groups(sample_size)
    sequence.append(noise)
sequence_inverse = np.linalg.multi_dot(sequence[::-1]).conj().T

现在这sequence是一个多维数组，我在将多维附加noise到 内的正确顺序时遇到了麻烦sequence，然后在执行linalg.multidot逆向sequence并取其 Hermitian 共轭时也遇到了问题。在这种情况下，我想要对与每个循环对应的每一行multi_dot存储的所有内容进行反转。如何才能做到这一点？noisejj

我将在下面提供一些“伪示例”来进一步演示我的问题，使用j= 3。为简单起见，这里我将仅“随机绘制” X, Y, Z。

非向量化案例：

i = 1
sequence = []
    j = 1
    noise = X (randomised)
    sequence.append(noise)
    sequence = [X]
    j = 2
    noise = Y (randomised)
    sequence.append(noise)
    sequence = [X, Y]
    j = 3
    noise = Z (randomised)
    sequence.append(noise)
    sequence = [X, Y, Z]

    end of j loop
take reverse order: [Z, Y, X]
do multi_dot: [ZYX] (Note: dot products, not element-wise multiplication)
take conjugate and tranpose(to get Hermitian): [ZYX].conj().T = [ZYX.conj().T]

矢量化案例（假设我在做sample_size= 3）：

sequence = []
    j = 1
    noise = [X,Z,Y](randomised)
    sequence.append(noise)
    sequence = [[X,Z,Y]]
    j = 2
    noise = [Z,Y,X] (randomised)
    sequence.append(noise)
    sequence = [[X,Z,Y],
                [Z,Y,X]]
    j = 3
    noise = [Z,Z,X] (randomised)
    sequence.append(noise)
    sequence = [[X,Z,Y],
                [Z,Y,X],
                [Z,Z,X]]
    end of j loop
take reverse order: [[Z,Z,X],
                     [Z,Y,X],
                     [X,Z,Y]]
do multi_dot(along an axis,
which is what I have trouble with): [ZZX,ZYZ,XXY]
take conjugate and tranpose(to get Hermitian): 
[ZZX,ZYZ,XXY].conj().T = [ZZX.conj().T, ZYZ.conj().T, XXY.conj().T]

我希望这些例子能说明我的问题

标签： numpymultidimensional-arraynumpy-ndarrayarray-broadcastingdot-product

解决方案

使用您的两个随机选择器：

In [13]: operator_groups()        # returns one (2,2) array                                                                                   
Out[13]: 
array([[-0.5+0.5j,  0.5-0.5j],
       [-0.5-0.5j, -0.5-0.5j]])
In [14]: operator_groups1(4)      # returns a list of (2,2) arrays                                                                           
Out[14]: 
[array([[0.+0.j, 1.+0.j],
        [1.+0.j, 0.+0.j]]), array([[ 0.+0.j, -0.-1.j],
        [ 0.+1.j,  0.+0.j]]), array([[ 1.+0.j,  0.+0.j],
        [ 0.+0.j, -1.+0.j]]), array([[1.+0.j, 0.+0.j],
        [0.+0.j, 1.+0.j]])]

您的循环创建一个数组列表：

In [15]: seq=[] 
    ...: for j in range(4): 
    ...:     seq.append(operator_groups()) 
    ...:                                                                                                     
In [16]: seq                                                                                                 
Out[16]: 
[array([[-0.5-0.5j, -0.5+0.5j],
        [ 0.5+0.5j, -0.5+0.5j]]), array([[1.+0.j, 0.+0.j],
        [0.+0.j, 1.+0.j]]), array([[-0.5+0.5j, -0.5-0.5j],
        [ 0.5-0.5j, -0.5-0.5j]]), array([[-0.5-0.5j,  0.5-0.5j],
        [-0.5-0.5j, -0.5+0.5j]])]

可以给multi_dot顺序打点：

In [17]: np.linalg.multi_dot(seq)                                                                            
Out[17]: 
array([[0.-1.j, 0.+0.j],
       [0.+0.j, 0.+1.j]])

如果我们使用组选择器构建序列，我们会得到一个列表列表：

In [18]: seq=[] 
    ...: for j in range(4): 
    ...:     seq.append(operator_groups1(3)) 
    ...:                                                                                                     
In [19]: seq                                                                                                 
Out[19]: 
[[array([[ 0.+0.j, -0.-1.j],
         [ 0.+1.j,  0.+0.j]]), array([[ 1.+0.j,  0.+0.j],
         [ 0.+0.j, -1.+0.j]]), array([[0.+0.j, 1.+0.j],
         [1.+0.j, 0.+0.j]])], [array([[ 0.+0.j, -0.-1.j],
         [ 0.+1.j,  0.+0.j]]), array([[ 1.+0.j,  0.+0.j],
         [ 0.+0.j, -1.+0.j]]), array([[0.+0.j, 1.+0.j],
         [1.+0.j, 0.+0.j]])], [array([[1.+0.j, 0.+0.j],
         [0.+0.j, 1.+0.j]]), array([[ 1.+0.j,  0.+0.j],
         [ 0.+0.j, -1.+0.j]]), array([[ 0.+0.j, -0.-1.j],
         [ 0.+1.j,  0.+0.j]])], [array([[1.+0.j, 0.+0.j],
         [0.+0.j, 1.+0.j]]), array([[ 1.+0.j,  0.+0.j],
         [ 0.+0.j, -1.+0.j]]), array([[ 0.+0.j, -0.-1.j],
         [ 0.+1.j,  0.+0.j]])]]
In [20]: len(seq)                                                                                            
Out[20]: 4
In [21]: len(seq[0])                                                                                         
Out[21]: 3

我们可以“堆叠”内部列表，创建一个 (n,2,2) 数组列表：

In [22]: seq1 = [np.stack(el) for el in seq]                                                                 
In [23]: seq1                                                                                                
Out[23]: 
[array([[[ 0.+0.j, -0.-1.j],
         [ 0.+1.j,  0.+0.j]],

        [[ 1.+0.j,  0.+0.j],
         [ 0.+0.j, -1.+0.j]],

        [[ 0.+0.j,  1.+0.j],
         [ 1.+0.j,  0.+0.j]]]), array([[[ 0.+0.j, -0.-1.j],
         [ 0.+1.j,  0.+0.j]],

        [[ 1.+0.j,  0.+0.j],
         [ 0.+0.j, -1.+0.j]],

        [[ 0.+0.j,  1.+0.j],
         [ 1.+0.j,  0.+0.j]]]), array([[[ 1.+0.j,  0.+0.j],
         [ 0.+0.j,  1.+0.j]],

        [[ 1.+0.j,  0.+0.j],
         [ 0.+0.j, -1.+0.j]],

        [[ 0.+0.j, -0.-1.j],
         [ 0.+1.j,  0.+0.j]]]), array([[[ 1.+0.j,  0.+0.j],
         [ 0.+0.j,  1.+0.j]],

        [[ 1.+0.j,  0.+0.j],
         [ 0.+0.j, -1.+0.j]],

        [[ 0.+0.j, -0.-1.j],
         [ 0.+1.j,  0.+0.j]]])]

然后我们可以matmul在这个列表上重复申请：

In [25]: res = seq1[0] 
    ...: for el in seq1[1:]: 
    ...:     res = res@el 
    ...:      
    ...:                                                                                                     
In [26]: res                                                                                                 
Out[26]: 
array([[[1.+0.j, 0.+0.j],
        [0.+0.j, 1.+0.j]],

       [[1.+0.j, 0.+0.j],
        [0.+0.j, 1.+0.j]],

       [[1.+0.j, 0.+0.j],
        [0.+0.j, 1.+0.j]]])

实际上matmul就像dot，但它将前导维度视为“批量”维度。

使用随机选择比较不同的结果很痛苦（除非我设置了种子），所以我把验证留给你。

推荐阅读