Classifying sparse matrix in parfor loop

4 vues (au cours des 30 derniers jours)
Bjorn
Bjorn le 11 Nov 2021
Commenté : Matt J le 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 commentaire
Matt J
Matt J le 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.

Connectez-vous pour commenter.

Réponse acceptée

Matt J
Matt J le 11 Nov 2021
Modifié(e) : Matt J le 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 commentaires
Bjorn
Bjorn le 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 le 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.

Connectez-vous pour commenter.

Plus de réponses (0)

Catégories

En savoir plus sur Loops and Conditional Statements dans Help Center et File Exchange

Produits


Version

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by