Classifying sparse matrix in parfor loop

5 visualizzazioni (ultimi 30 giorni)
Bjorn
Bjorn il 11 Nov 2021
Commentato: Matt J il 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 Commento
Matt J
Matt J il 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.

Accedi per commentare.

Risposta accettata

Matt J
Matt J il 11 Nov 2021
Modificato: Matt J il 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 Commenti
Bjorn
Bjorn il 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 il 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.

Accedi per commentare.

Più risposte (0)

Categorie

Scopri di più su Loops and Conditional Statements in Help Center e File Exchange

Prodotti


Release

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by