Classifying sparse matrix in parfor loop

5 Ansichten (letzte 30 Tage)
Bjorn
Bjorn am 11 Nov. 2021
Kommentiert: Matt J am 13 Nov. 2021
Dear all,
The following for-loop is defined:
iend = nz - 1;
jend = nr - 1;
for i = 1 : iend
for j = 1 : jend
k = i * Nr + j;
dz = z_grid(i+1) - z_grid(i);
dr = r_grid(j+1) - r_grid(j);
A_sol(k,k) = - dz * (r_grid(j+1) * a(i+1,j+1) + r_grid(j) * b(i+1,j)) - dr * r_grid(j+1) * (c(i+1,j+1) + d(i,j+1));
end
end
with z_grid and r_grid being vectors and a, b, c and d matrices. The loop takes long to compute and I would like to use parloop to speed up the process.
From MATLAB Help I know that the code should be written in the form of:
iend = nz - 1;
jend = nr - 1;
for i = 1 : iend
for j = 1 : jend
k(j) = i * Nr + j;
dr = r_grid(j+1) - r_grid(j);
end
dz = z_grid(i+1) - z_grid(i);
A_sol(k,k) = -dz* (r_grid(j+1) * a(i+1,j+1) + r_grid(j) * b(i+1,j)) - dr * r_grid(j+1) * (c(i+1,j+1) + d(i,j+1));
end
And write A_sol(k,k) in the form A_sol(k,:). However, the rows and columns have the same indexation, and I cannot come up with a solution to still make this work. Could someone help me out?
Thank you in advance,
WIth kind regards,
Bjorn
  1 Kommentar
Matt J
Matt J am 11 Nov. 2021
Shouldn't k be given by,
k = (i-1) * Nr + j;
As you have it now, k will not start at k=1, but instead it will run from Nr+1 to iend*Nr+jend.

Melden Sie sich an, um zu kommentieren.

Akzeptierte Antwort

Matt J
Matt J am 11 Nov. 2021
Bearbeitet: Matt J am 11 Nov. 2021
The loop takes long to compute and I would like to use parloop to speed up the process.
I doubt it will, but the loop can be further optimized by pre-computing some things. Also, it will be much faster if A_sol is computed using the sparse() command.
r_grid=r_grid(:).'; %make sure these are row vectors
z_grid=z_grid(:).';
DZ=[0,diff(z_grid)].';
DR=[0,diff(r_grid)];
Term1=(-DZ .* a - DR .* c ).*r_grid;
Term2=(-DZ .* b).*r_grid;
Term3=(-DR .* d).*r_grid;
N=length(A_sol);
A_diag=zeros(N,1);
iend = nz - 1;
jend = nr - 1;
for i = 1 : iend
for j = 1 : jend
k = i * Nr + j;
A_diag(k) = Term1(i+1,j+1) +Term2(i+1,j) + Term3(i,j+1);
end
end
A_sol=spdiags(A_diag(:),0,N,N);
  2 Kommentare
Bjorn
Bjorn am 13 Nov. 2021
Dear Matt,
Thank you for your promt answer! Defining the matrix after the double for-loop and defining the terms before the double for loop decreased the computational time by up to a factor of ~10 for 2^10 x 2^10 matrices.
Thank you again,
WIth kind regards,
Bjorn
Matt J
Matt J am 13 Nov. 2021
Glad it helped but 2^10 is pretty small. You could probably make things faster just by using a non-sparse matrix.

Melden Sie sich an, um zu kommentieren.

Weitere Antworten (0)

Kategorien

Mehr zu Loops and Conditional Statements finden Sie in Help Center und File Exchange

Produkte


Version

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by