如何向量化在Matlab

Question

天线arrayfactor表达我有天线阵列因子表达这里：

我已编码的阵列因子表达下面给出：

lambda = 1; 
M = 100;N = 200; %an M x N array 
dx = 0.3*lambda; %inter-element spacing in x direction 
m = 1:M; 
xm = (m - 0.5*(M+1))*dx; %element positions in x direction 

dy = 0.4*lambda; 
n = 1:N; 
yn = (n - 0.5*(N+1))*dy; 

thetaCount = 360; % no of theta values 

thetaRes = 2*pi/thetaCount; % theta resolution 

thetas = 0:thetaRes:2*pi-thetaRes; % theta values 

phiCount = 180; 

phiRes = pi/phiCount; 

phis = -pi/2:phiRes:pi/2-phiRes; 

cmpWeights = rand(N,M); %complex Weights 

AF = zeros(phiCount,thetaCount); %Array factor 

tic 
for i = 1:phiCount 
    for j = 1:thetaCount 

     for p = 1:M 
      for q = 1:N 

       AF(i,j) = AF(i,j) + cmpWeights(q,p)*exp((2*pi*1j/lambda)*(xm(p)*sin(thetas(j))*cos(phis(i)) + yn(q)*sin(thetas(j))*sin(phis(i)))); 

      end 
     end 
    end 
end 

如何能我矢量化了计算阵列因子（AF）的代码。

欲行：

AF(i,j) = AF(i,j) + cmpWeights(q,p)*exp((2*pi*1j/lambda)*(xm(p)*sin(thetas(j))*cos(phis(i)) + yn(q)*sin(thetas(j))*sin(phis(i)))); 

要写入量化形式（通过修改for循环）。

Answer 1

方法一：全油门

内部嵌套循环产生这种每次迭代 - cmpWeights(q,p)*exp((2*pi*1j/lambda)*(xm(p)*sin(thetas(j))*cos(phis(i)) + yn(q)*sin(thetas(j))*sin(phis(i))))，这是总结了反复给我们AF最终输出。

让我们打电话exp(....部分为B。现在，B主要有两个部分，一个是标量(2*pi*1j/lambda)和从依赖于 原糊涂版本中使用的四个迭代变量构成的另一部分 (xm(p)*sin(thetas(j))*cos(phis(i)) + yn(q)*sin(thetas(j))*sin(phis(i))) - ​​。为便于以后参考，我们将其他部分称为C。

让我们把所有的进入角度：

• 糊涂的版本有AF(i,j) = AF(i,j) + cmpWeights(q,p)*exp((2*pi*1j/lambda)*(xm(p)*sin(thetas(j))*cos(phis(i)) + yn(q)*sin(thetas(j))*sin(phis(i))))，这是相当于现在AF(i,j) = AF(i,j) + cmpWeights(q,p)*B，其中B = exp((2*pi*1j/lambda)*(xm(p)*sin(thetas(j))*cos(phis(i)) + yn(q)*sin(thetas(j))*sin(phis(i))))。

• B可以简化为B = exp((2*pi*1j/lambda)* C)，其中C = (xm(p)*sin(thetas(j))*cos(phis(i)) + yn(q)*sin(thetas(j))*sin(phis(i)))。

• C将取决于迭代器 -​​。

因此，移植到一个量化的方式之后，它最终会因为这 -

%// 1) Define vectors corresponding to iterators used in the loopy version 
I = 1:phiCount; 
J = 1:thetaCount; 
P = 1:M; 
Q = 1:N; 

%// 2) Create vectorized version of C using all four vector iterators 
mult1 = bsxfun(@times,sin(thetas(J)),cos(phis(I)).'); %//' 
mult2 = bsxfun(@times,sin(thetas(J)),sin(phis(I)).'); %//' 

mult1_xm = bsxfun(@times,mult1(:),permute(xm,[1 3 2])); 
mult2_yn = bsxfun(@times,mult2(:),yn); 
C_vect = bsxfun(@plus,mult1_xm,mult2_yn); 

%// 3) Create vectorized version of B using vectorized C 
B_vect = reshape(exp((2*pi*1j/lambda)*C_vect),phiCount*thetaCount,[]); 

%// 4) Final output as matrix multiplication between vectorized versions of B and C 
AF_vect = reshape(B_vect*cmpWeights(:),phiCount,thetaCount); 

---

方法2：更少的内存密集型

这第二种方法会减少内存流量，并使用指数分布式属性 - exp(A+B) = exp(A)*exp(B)。

现在，原来糊涂的版本是这样的 -

AF(i,j) = AF(i,j) + cmpWeights(q,p)*exp((2*pi*1j/lambda)*... 
    (xm(p)*sin(thetas(j))*cos(phis(i)) + yn(q)*sin(thetas(j))*sin(phis(i)))) 

因此，使用分配律后，我们会像这样的东西endup -

K = (2*pi*1j/lambda) 
part1 = K*xm(p)*sin(thetas(j))*cos(phis(i)); 
part2 = K*yn(q)*sin(thetas(j))*sin(phis(i)); 
AF(i,j) = AF(i,j) + cmpWeights(q,p)*exp(part1)*exp(part2); 

因此，相关的量化方法将成为这样的东西 -

%// 1) Define vectors corresponding to iterators used in the loopy version 
I = 1:phiCount; 
J = 1:thetaCount; 
P = 1:M; 
Q = 1:N; 

%// 2) Define the constant used at the start of EXP() call 
K = (2*pi*1j/lambda); 

%// 3) Perform the sine-cosine operations part1 & part2 in vectorized manners 
mult1 = K*bsxfun(@times,sin(thetas(J)),cos(phis(I)).'); %//' 
mult2 = K*bsxfun(@times,sin(thetas(J)),sin(phis(I)).'); %//' 

%// Perform exp(part1) & exp(part2) in vectorized manners 
part1_vect = exp(bsxfun(@times,mult1(:),xm)); 
part2_vect = exp(bsxfun(@times,mult2(:),yn)); 

%// Perform multiplications with cmpWeights for final output 
AF = reshape(sum((part1_vect*cmpWeights.').*part2_vect,2),phiCount,[]) 

---

快速标杆

下面是与原始糊涂的做法的问题列出的输入数据运行时和提出的方法＃2 -

---------------------------- With Original Approach 
Elapsed time is 358.081507 seconds. 

---------------------------- With Proposed Approach #2 
Elapsed time is 0.405038 seconds. 

的运行时建议用一个疯狂的性能改进Approach #2！

Answer 2

基本的技巧是找出什么东西是恒定的，什么东西取决于下标词 - 因此是矩阵项。

在的总和：

• C(n,m)是矩阵
• 2π/λ是常数
• sin(θ)cos(φ)是常数
• x(m)和y(n)是矢量

这样的两件事我会做的是：

1. 使用meshgrid()
2. 采取一切外循环常数项的东西展开​​和ym成矩阵。

像这样：

... 

piFactor = 2 * pi * 1j/lambda; 

[xgrid, ygrid] = meshgrid(xm, ym);          % xgrid and ygrid will be size (N, M) 

for i = 1:phiCount 
    for j = 1:thetaCount 

     xFactor = sin(thetas(j)) * cos(phis(i)); 
     yFactor = sin(thetas(j)) * sin(phis(i)); 

     expFactor = exp(piFactor * (xgrid * xFactor + ygrid * yFactor)); % expFactor is size (N, M) 

     elements = cmpWeights .* expFactor;        % elements of sum, size (N, M) 

     AF(i, j) = AF(i, j) + sum(elements(:));       % sum and then integrate. 

    end 
end 

你也许可以弄清楚如何vectorise外环太多，但希望这给你一个起点。