Nested for loops optimisation

Question

Virgile Favre 2016 年 4 月 25 日

0
リンク

この質問への直接リンク

https://jp.mathworks.com/matlabcentral/answers/280924-nested-for-loops-optimisation

コメント済み: Virgile Favre 2016 年 5 月 3 日

Hi there,

I am currently working on a project that involves computation. I need to determine the values of a Ns * Ns symmetric, real matrix.

At the moment I use two for loops to create this matrix, one (let us say with variable a) from 1 to Ns, and the other one ( with j) from j to Ns, in order to use the symmetry of the matrix. the code looks like this :

MATLAB code :
Ns = 30; nc = 1; m = 1; k = 6;
H = rand(Ns, Ns);
H = H + H'; % symmetric, thus will have real eigenvalues
[V, D] = eig(H);
T = Ns * nc / k;
lambda = diag(D);
Y = zeros(Ns, Ns);
for i = 1 : Ns
    for j = i : Ns
        for l = 1 : T
            Y(i, j) = Y(i, j) + real(conj(V(j, l)) * V(i, l) *...
                sum(conj(V(i, T + 1 : Ns)) .* V(j, T + 1 : Ns)...
                ./ (lambda(l) - lambda(T + 1 : Ns))'));
        end
        if i == j
            Y(i, i) = Y(i, i) / 2; % Useful to perform Y = (Y + Y') / 2
        end
    end
end 
Y = -k * (Y + Y'); % This is the actual Y in which I'm interested
% code
V is a Ns * Ns matrix and lambda is some vector of length Ns.
The problem I have is that T is never equal to Ns, thus I cannot write these operations as simple matrix elements, which can be very time-expensive, for example when Ns = 400, 500, ...

I would like to know if anyone could help me vectorize this code.

Edit : clarified somethings for this code to run : example values

2 件のコメント
なしを表示なしを非表示

John BG 2016 年 4 月 26 日

MATLAB Online で開く

can you get here a working script to start from? I got the following, please amend as appropriate so I can try translating your loops into bsxfun lines, if possible, if i know how to:

Ns=48;nc=10;k=1;V=randi(10,Ns,Ns);D=V;
T = Ns * nc / k;
lambda = diag(D);
Y = zeros(Ns, Ns);
for i = 1 : Ns
  for j = i : Ns
      for l = 1 : T
          Y(i, j) = Y(i, j) + real(conj(V(j, l)) * V(i, l) *sum(conj(V(i, T + 1 : Ns)) .* V(j, T + 1 : Ns)./ (lambda(l) - lambda(T + 1 : Ns))'));
      end
%       if i == j
%           Y(i, i) = Y(i, i) / 2; % Useful to perform X = (X + X') / 2
%       end
  end
end