2

我有一个 5D 矩阵 Cij(3,3,Nx,Ny,Nz)，其中 Nx,Ny 和 Nz 作为输入给出。

我需要执行以下操作：

for ikx=1:Nx,
    for iky=1:Ny,
        for ikz=1:Nz,

            %Random simulation of fourier components
            n=zeros((3),'double');
            for j=1:9,
                ncomponent=randn(2);
                n(j)=complex(ncomponent(1),ncomponent(2));
                %Calculation of H
                H(:,ikx,iky,ikz)=dot(Cij(:,:,ikx,iky,ikz),n);
            end;
        end;
    end;
end;

问题是增加 Nx,Ny,Nz 循环需要非常长的时间来计算 H 矩阵。

有人知道获得 H 矩阵的更快方法吗？

4

2 回答 2

3

首先应该注意的是，在最里面的循环中，您执行了 9 次点积，H(:,ikx,iky,ikz)每次都覆盖。这没有任何意义。您应该只填写循环内的随机值，并n 在该循环后计算H(:,ikx,iky,ikz)一次。

但是，所有循环都是不必要的，因为您可以利用函数DOT是矢量化的并且可以处理 5-D 数组的事实（即它会自动在第一个非单维上执行点操作）。您所要做的就是制作n一个 3×3×Nx×Ny×Nz 的复数值矩阵。这两行应该给你与上面的代码相同的结果：

n = complex(rand([3 3 Nx Ny Nz]), rand([3 3 Nx Ny Nz]));
H = squeeze(dot(Cij, n));

函数SQUEEZE用于从 中删除单个维度H，这将使其成为 3×Nx×Ny×Nz 矩阵。

于 2011-12-19T17:23:02.583 回答
1

您可以使用 some reshape(and permute)来做到这一点

C=rand(3,3,Nx,Ny,Nz);
n=rand(3,3);

如果你想在 n 和 C 的每个元素之间进行矩阵乘法：

H=reshape(n*reshape(C,3,3*Nx*Ny*Nz),[3,3,Nx,Ny,Nz])

如果你想要 n 和 C 的每个元素之间的点积：

H=reshape(reshape(n,1,[])*reshape(C,3*3,Nx*Ny*Nz),[Nx,Ny,Nz])
于 2011-12-19T17:25:23.210 回答