Re: [deal.II] Efficient Matrix-based Calculation

[Jean-Paul Pelteret]

Dear Masoud,

> As you mentioned, we expect this to be slower than the equivalent code in 
> step-8; however, it is very slower than the same code that I wrote before in 
> Matlab.

A quick question: are you sure that you’re running your problem in release 
mode, rather than debug mode? I had a colleague who adopted this BDB approach 
that you’ve described, and he never mentioned performance issues to the degree 
that you’re seeing here.

Best,
Jean-Paul

> On 20. Oct 2021, at 14:08, Masoud Ahmadi  wrote:
> 
> Dear Prof. Bangerth,
> 
> Thanks for your answer. 
> Indeed, I took the declaration of matrices out of the loop; but still, as I 
> mentioned, the part that is substantially reducing the speed is the product 
> between the matrices. 
> As you mentioned, we expect this to be slower than the equivalent code in 
> step-8; however, it is very slower than the same code that I wrote before in 
> Matlab. 
> My question is, I want to be sure that "FullMatrix" is the most efficient and 
> fast way in dealii to evaluate such like matrix calculations (the products). 
> Or should I be looking for alternatives?
> 
> Best regards
> 
> On Tuesday, 19 October 2021 at 19:24:11 UTC+1 Wolfgang Bangerth wrote:
> 
> Masoud, 
> 
> > Here is a part of the tutorial code for example 8 to calculate stiffness 
> > matrix of an element: 
> > Screenshot 2021-10-15 at 12.48.39.png 
> > 
> > And here is my equivalent code for the same purpose with a matrix-based 
> > viewpoint: 
> > Screenshot 2021-10-15 at 23.12.13.png 
> > 
> > it is based on the following formulation for a 3D elastic problem: 
> > thumbnail_image.png 
> > By "slow", I mean in comparison with the same code in tutorial (example 8). 
> > I 
> > measured the lines of the code, and I am sure the difference is for the 
> > mentioned part of the code. To be more accurate: the operations on 
> > matrices, 
> > as follow in the second image: 
> > B[I].Tmmult(tmpM,D); 
> > tmpM.mmult(BDB,B[J]); 
> 
> You could probably make the code substantially faster already if you moved 
> the 
> declaration of the B, tmpM, and BDB matrices out of the loop over the 
> quadrature points, and simply set them to zero inside the loop. Creating 
> these 
> variables requires allocating memory every time you do one iteration over the 
> quadrature points, which is expensive. 
> 
> But in the end, compare what work you are doing in your loop with what the 
> equivalent code in step-8 is doing, and you shouldn't be surprised by your 
> observation. step-8 only does 2 (and if i==j, then 3) tensor contractions for 
> each index i,j,q. You are allocating and releasing memory, and doing two 
> products between 6x3 and 3x3 matrices, plus a lot of index work. It does not 
> surprise me that this is slow. 
> 
> Best 
> W. 
> 
> -- 
>  
> Wolfgang Bangerth email: bang...@colostate.edu 
>  
> www: http://www.math.colostate.edu/~bangerth/ 
>  
> 
> 
> -- 
> The deal.II project is located at http://www.dealii.org/ 
> 
> For mailing list/forum options, see 
> https://groups.google.com/d/forum/dealii?hl=en 
> 
> --- 
> You received this message because you are subscribed to the Google Groups 
> "deal.II User Group" group.
> To unsubscribe from this group and stop receiving emails from it, send an 
> email to dealii+unsubscr...@googlegroups.com 
> .
> To view this discussion on the web visit 
> https://groups.google.com/d/msgid/dealii/00c2608c-7489-4538-8bd7-09c989d844c4n%40googlegroups.com
>  
> .

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/E3BA93D2-831E-453F-9CD6-1818587CC00D%40gmail.com.

Re: [deal.II] Efficient Matrix-based Calculation

[Masoud Ahmadi]

Dear Prof. Bangerth,

Thanks for your answer. 
Indeed, I took the declaration of matrices out of the loop; but still, as I 
mentioned, the part that is substantially reducing the speed is the product 
between the matrices. 
As you mentioned, we expect this to be slower than the equivalent code in 
step-8; however, it is very slower than the same code that I wrote before 
in Matlab. 
My question is, I want to be sure that "FullMatrix" is the most efficient 
and fast way in dealii to evaluate such like matrix calculations (the 
products). Or should I be looking for alternatives?

Best regards

On Tuesday, 19 October 2021 at 19:24:11 UTC+1 Wolfgang Bangerth wrote:

>
> Masoud,
>
> > Here is a part of the tutorial code for example 8 to calculate stiffness 
> > matrix of an element:
> > Screenshot 2021-10-15 at 12.48.39.png
> > 
> > And here is my equivalent code for the same purpose with a matrix-based 
> viewpoint:
> > Screenshot 2021-10-15 at 23.12.13.png
> > 
> > it is based on the following formulation for a 3D elastic problem:
> > thumbnail_image.png
> > By "slow", I mean in comparison with the same code in tutorial (example 
> 8). I 
> > measured the lines of the code, and I am sure the difference is for the 
> > mentioned part of the code. To be more accurate: the operations on 
> matrices, 
> > as follow in the second image:
> > B[I].Tmmult(tmpM,D);
> > tmpM.mmult(BDB,B[J]);
>
> You could probably make the code substantially faster already if you moved 
> the 
> declaration of the B, tmpM, and BDB matrices out of the loop over the 
> quadrature points, and simply set them to zero inside the loop. Creating 
> these 
> variables requires allocating memory every time you do one iteration over 
> the 
> quadrature points, which is expensive.
>
> But in the end, compare what work you are doing in your loop with what the 
> equivalent code in step-8 is doing, and you shouldn't be surprised by your 
> observation. step-8 only does 2 (and if i==j, then 3) tensor contractions 
> for 
> each index i,j,q. You are allocating and releasing memory, and doing two 
> products between 6x3 and 3x3 matrices, plus a lot of index work. It does 
> not 
> surprise me that this is slow.
>
> Best
> W.
>
> -- 
> 
> Wolfgang Bangerth email: bang...@colostate.edu
> www: http://www.math.colostate.edu/~bangerth/
>
>

-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/00c2608c-7489-4538-8bd7-09c989d844c4n%40googlegroups.com.

Re: [deal.II] Efficient Matrix-based Calculation

[Wolfgang Bangerth]

Masoud,

Here is a part of the tutorial code for example 8 to calculate stiffness 
matrix of an element:

Screenshot 2021-10-15 at 12.48.39.png

And here is my equivalent code for the same purpose with a matrix-based 
viewpoint:
Screenshot 2021-10-15 at 23.12.13.png

it is based on the following formulation for a 3D elastic problem:
thumbnail_image.png
By "slow", I mean in comparison with the same code in tutorial (example 8). I 
measured the lines of the code, and I am sure the difference is for the 
mentioned part of the code. To be more accurate: the operations on matrices, 
as follow in the second image:

B[I].Tmmult(tmpM,D);
tmpM.mmult(BDB,B[J]);


You could probably make the code substantially faster already if you moved the 
declaration of the B, tmpM, and BDB matrices out of the loop over the 
quadrature points, and simply set them to zero inside the loop. Creating these 
variables requires allocating memory every time you do one iteration over the 
quadrature points, which is expensive.


But in the end, compare what work you are doing in your loop with what the 
equivalent code in step-8 is doing, and you shouldn't be surprised by your 
observation. step-8 only does 2 (and if i==j, then 3) tensor contractions for 
each index i,j,q. You are allocating and releasing memory, and doing two 
products between 6x3 and 3x3 matrices, plus a lot of index work. It does not 
surprise me that this is slow.


Best
 W.

--

Wolfgang Bangerth  email: bange...@colostate.edu
   www: http://www.math.colostate.edu/~bangerth/

--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups "deal.II User Group" group.

To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/b674a16b-bcbe-7bcf-5ccb-cb659154e1f7%40colostate.edu.

Re: [deal.II] Efficient Matrix-based Calculation

[Wolfgang Bangerth]

On 10/12/21 10:33 AM, Masoud Ahmadi wrote:


I'm trying to use a Matrix-based solution for a general 3D FE elastic 
problem (just like example-8). In that example the authors used a linear 
scaler solution for calculating the stiffness matrix of each element; 
but, due to some reasons I want to use a matrix-based solution so that 
the stiffness matrix of each element can be calculated as follows:

K_cell = B^T . D . B,
where B matrix follows from the differential operator for strain 
calculation and D is the usual elasticity matrix for 3D isotropic 
problems. The size of the B and D matrices are 6X3 and 6X6 respectively. 
I defined these matrices as "vector" of "FullMatrix" and try to 
do operate on them by FullMatrix operators as follows:

B[i].Tmmult(tmpMat,D);
tmpMat.mmult(BDB,B[j]);
and it works correctly; however, it works rather slowly.
My question is, is "FullMatrix" the most efficient and fast way in 
dealii to evaluate such like matrix calculations? Are there any 
alternatives?


Masoud,
it's difficult to say without seeing the actual code. As a general rule, 
it is impossible to tell for even good programmers where the "slow" 
parts of a code are unless one actually uses a profiling tool to measure 
each line of the code. So I would suggest that you figure out which 
line(s) is actually slow in your code, and then asking how that one part 
can be improved.


My best guess is that it is not actually the FullMatrix itself, but how 
you build the matrix. But as I said, it's not possible to say so without 
actual benchmarking.


Best
 W.

--

Wolfgang Bangerth  email: bange...@colostate.edu
   www: http://www.math.colostate.edu/~bangerth/

--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups "deal.II User Group" group.

To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/3230d7f8-06b8-1e55-2ce2-a8ed05455402%40colostate.edu.

[deal.II] Efficient Matrix-based Calculation

[Masoud Ahmadi]

Dear all,

I'm trying to use a Matrix-based solution for a general 3D FE elastic 
problem (just like example-8). In that example the authors used a linear 
scaler solution for calculating the stiffness matrix of each element; but, 
due to some reasons I want to use a matrix-based solution so that the 
stiffness matrix of each element can be calculated as follows:
K_cell = B^T . D . B,
where B matrix follows from the differential operator for strain calculation 
and D is the usual elasticity matrix for 3D isotropic problems. The size of 
the B and D matrices are 6X3 and 6X6 respectively. I defined these matrices 
as "vector" of "FullMatrix" and try to do operate on them by 
FullMatrix operators as follows:
B[i].Tmmult(tmpMat,D);
tmpMat.mmult(BDB,B[j]);
and it works correctly; however, it works rather slowly.
My question is, is "FullMatrix" the most efficient and fast way in dealii 
to evaluate such like matrix calculations? Are there any alternatives?

Best regards,
Masoud




-- 
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see 
https://groups.google.com/d/forum/dealii?hl=en
--- 
You received this message because you are subscribed to the Google Groups 
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email 
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/711e55f4-df0c-4673-8c29-b9cba6cc2185n%40googlegroups.com.