Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.
6977 Discussions

Intel MKL PARDISO : reordering once, factorizing multiple times

mlucio89
Beginner
1,366 Views

Hi all,

I have been using MKL PARDISO for two months now to solve linear systems within a Fortran-coded non-linear Finite Element software. At the moment I try to optimize the code and in particular the calls to MKL PARDISO.

In my problem, the coefficient matrix (tangent stiffness matrix) has a constant sparsity pattern, however the coefficients values change at each nonlinear iteration.

If I am correct, in this case it makes sense to reorder and symbolically factorize (phase 11) only once -at the beginning- and then factorize (phase 22) each time the coefficients change.

If I do so, I get no errors but the results become wrong after the first iteration.

Conversely, the results are correct if I perform the two phases 11+22 at each nonlinear iteration.

Is my hypothesis fundamentally wrong, to want to perform phase 11 only once? Or is there something else going on here that I do not see?

Thank you very much for your help!

ML

PS: I am NOT using the DSS interface

PPS: I have already tried to set iparm(5)=2 to save the permutation vector computed in phase 11

0 Kudos
1 Solution
Kirill_V_Intel
Employee
1,340 Views

Hello,

Your general line of thinking is correct but there are important details to understand.

1) If you use matching or scaling, then the result of phase 1 actually depends on your matrix values. So it is almost impossible (and PARDISO cannot do it now) re-use the phase 1's information in this case. So, you need to turn off these features.

2) Second, without special notification PARDISO cannot work like you tried, with calling phase 22 once, then "silently" changing the matrix values and then calling phase 22 again.

3) What PARDISO can do for you is three things.
First, if you want to keep older matrices and the number of different matrices is small, you can use mnum/maxfct parameters to force PARDISO store multiple factorizations in a single handle. I believe though it is not your case.
Second, there is so called "low rank" feature of PARDISO, see https://software.intel.com/content/www/us/en/develop/documentation/onemkl-developer-reference-c/top/sparse-solver-routines/onemkl-pardiso-parallel-direct-sparse-solver-interface.html (section "Low Rank Update"). You can specify which part of matrix values change (possibly, all of them) and then call phase 22 with this information. Then factorization will be re-computed only for the part which changed. But you need to know the location of the changes. If the changes are for all matrix values, it won't bring any benefit.

Third, if the matrix values are changing smoothly and not much, you can only do iterative refinement (essentially, use the previous matrix factors as a good approximation to the inverse). Then you only call phase 33 with iterative refinement, without re-doing factorization at all. For exmaple, you can re-compute factorization from scratch once you matrix changes are too large and not every nonlinear iteration. Just set the number of iteration steps to be reasonable (iparm(8)).

I hope there is something useful for you in my reply.

Best,
Kirill

 

 

View solution in original post

0 Kudos
3 Replies
Kirill_V_Intel
Employee
1,341 Views

Hello,

Your general line of thinking is correct but there are important details to understand.

1) If you use matching or scaling, then the result of phase 1 actually depends on your matrix values. So it is almost impossible (and PARDISO cannot do it now) re-use the phase 1's information in this case. So, you need to turn off these features.

2) Second, without special notification PARDISO cannot work like you tried, with calling phase 22 once, then "silently" changing the matrix values and then calling phase 22 again.

3) What PARDISO can do for you is three things.
First, if you want to keep older matrices and the number of different matrices is small, you can use mnum/maxfct parameters to force PARDISO store multiple factorizations in a single handle. I believe though it is not your case.
Second, there is so called "low rank" feature of PARDISO, see https://software.intel.com/content/www/us/en/develop/documentation/onemkl-developer-reference-c/top/sparse-solver-routines/onemkl-pardiso-parallel-direct-sparse-solver-interface.html (section "Low Rank Update"). You can specify which part of matrix values change (possibly, all of them) and then call phase 22 with this information. Then factorization will be re-computed only for the part which changed. But you need to know the location of the changes. If the changes are for all matrix values, it won't bring any benefit.

Third, if the matrix values are changing smoothly and not much, you can only do iterative refinement (essentially, use the previous matrix factors as a good approximation to the inverse). Then you only call phase 33 with iterative refinement, without re-doing factorization at all. For exmaple, you can re-compute factorization from scratch once you matrix changes are too large and not every nonlinear iteration. Just set the number of iteration steps to be reasonable (iparm(8)).

I hope there is something useful for you in my reply.

Best,
Kirill

 

 

0 Kudos
mlucio89
Beginner
1,312 Views

Hi @Kirill_V_Intel ,

Thank you very much for your answer. As you pointed out, I think my only option for the moment is probably to perform phase 1 and 2 together, each time.

However, I'll give iterative refinement a thought for the future.

Have a nice day,

ML

0 Kudos
PrasanthD_intel
Moderator
1,253 Views

Hi Marco,


Thanks for the confirmation.

As your issue has been resolved, we are closing this thread. If you require additional assistance from Intel, please start a new thread. Any further interaction in this thread will be considered community only


Regards

Prasanth


0 Kudos
Reply