Optimization of nested "for" loops

8 次查看(过去 30 天)
Hi.
I would like to optimize my code. I examined it and found a serious bottle neck. Below is the chunk of code with nested "for" loops which in my opinion are very inefficient.
PHI = zeros(length(X),length(X));
PHI_X = zeros(length(X),length(X));
PHI_XX = zeros(length(X),length(X));
PHI_Y = zeros(length(X),length(X));
PHI_YY = zeros(length(X),length(X));
for i=1:length(X)
for j=1:length(Y)
x = X(i)-X(j);
y = Y(i)-Y(j);
PHI (i,j) = (x^2 + y^2 + e^2) ^ p;
PHI_X (i,j) = 2 * p * x * (x^2 + y^2 + e^2) ^ (p-1);
PHI_Y (i,j) = 2 * p * y * (x^2 + y^2 + e^2) ^ (p-1);
PHI_XX (i,j) = 2 * p * ((2*p-1)*x^2 + e^2) * (x^2 + y^2 + e^2) ^ (p-2);
PHI_YY (i,j) = 2 * p * ((2*p-1)*y^2 + e^2) * (x^2 + y^2 + e^2) ^ (p-2);
end
end
X and Y are vectors while e and p are real-valued coefficients. This part of the code is very important because it is executed even several hundred times (for different values of e and p) during one code execution. I would like to optimize it somehow but I don't know how because I'm quite a new user of MATLAB.
Kind regards

采纳的回答

DGM
DGM 2021-5-3
Should be able to do something like this without the loop (don't need to preallocate either):
x = X'-X;
y = Y'-Y;
PHI = (x.^2 + y.^2 + e^2) .^ p;
PHI_X = 2 * p * x .* (x.^2 + y.^2 + e^2) .^ (p-1);
PHI_Y = 2 * p * y .* (x.^2 + y.^2 + e^2) .^ (p-1);
PHI_XX = 2 * p * ((2*p-1)*x.^2 + e^2) .* (x.^2 + y.^2 + e^2) .^ (p-2);
PHI_YY = 2 * p * ((2*p-1)*y.^2 + e^2) .* (x.^2 + y.^2 + e^2) .^ (p-2);
Those first couple lines will only work in R2016b or later, otherwise, you'd need to do:
x = bsxfun(@minus,X,X');
etc

更多回答(1 个)

Matt J
Matt J 2021-5-3
编辑:Matt J 2021-5-3
Replace the loop with,
x=X(:)-X(:).';
y=Y(:)-Y(:).';
xsq=x.^2;
ysq=y.^2;
xysq= xsq+ysq;
tmp1=xysq + e^2;
PHI = tmp1 .^ p;
tmp2=tmp1 .^ (p-1);
PHI_X = 2 * p * x .*tmp2;
PHI_Y = 2 * p * y .*tmp2;
clear tmp2
tmp3=tmp1 .^ (p-2);
PHI_XX = (2*p*(2*p-1)*xsq + 2*p*e^2) .* tmp3;
PHI_YY = (2*p*(2*p-1)*ysq + 2*p*e^2) .* tmp3;

类别

Help CenterFile Exchange 中查找有关 Solver Outputs and Iterative Display 的更多信息

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by