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

sparse Ax=b with different right sides

Alexey1978
New Contributor I
552 Views

Hi,

 

A have an iteration process.

At each step I need to solve Ax=b with the same sparse big matrix A, but with different b.

What is more effective? 1. To solve Ax=b using pardiso at each step or 2. find invesrion A^{-1} and multiply at each step.

 

The size of A ~ 50 000x50 000.

 

thanks a lot,

Alex

0 Kudos
5 Replies
mecej4
Honored Contributor III
541 Views

If what you want to do is to solve linear equations, do not form the inverse! Forming the inverse is equivalent to solving A.X = B, where B = I, the identity matrix. In your case, that means solving with 50,000 right hand sides. In certain cases, it can also be less accurate to form and multiply by the inverse.

Call Pardiso with phase=12 to analyze and factor the matrix. Then, call Pardiso with phase = 33. If you know m values of b, you set nrhs = m in this call, after entering b1, b2, ... into columns-1, 2, ... of B. If you only know one value of b at a time, call with nrhs = 1 when you have that b ready and wish to solve. Remember that you can do this only if the matrix itself remains unchanged between calls.

0 Kudos
Alexey1978
New Contributor I
535 Views
Thanks a lot.
 
The problem is that phase=12  and phase = 33 will be realized in different functions. Where can I store the results of factorization (phase=12) and how to transfer it to further functions realizing phase=33.
 
Namely 
 
func1
{
phase=12  
PARDISO (pt, &maxfct, &mnum, &mtype, &phase,
             &n, a, ia, ja, &idum, &nrhs, iparm, &msglvl, &ddum, &ddum, &error);
}
 
and
 
func2
{
phase=33
PARDISO (pt, &maxfct, &mnum, &mtype, &phase,
             &n, a, ia, ja, &idum, &nrhs, iparm, &msglvl, b, x, &error);
}
 
As I understand, func1 I will called 1 time, when func2 will be called many times. Matrix a is constant, but right side b will be calculated in iteration process.
 
What arguments should I store and transfer from func1 into func2 in my case.
 
Thank you.
0 Kudos
mecej4
Honored Contributor III
528 Views

As long as you do not change arguments such as pt, n, a, ia, ja, between successive calls to Pardiso, the solver state will be preserved. Arguments with "dum" are not used in the calculation and are just place holders. Pardiso manages storing the factors of the matrix, the pivot array, etc., internally, and you do not need to participate in those tasks.

How do you preserve the arguments to Pardiso that do need to be saved? That depends on the programming language that you use and your programming style. In Fortran, for example, those could be variables that are placed in a module that is used in the various functions func1, func2, etc. In C, they could be global variables. Or, they could be added to the argument lists of func1, func2, etc. Or, you could put all the calls to Pardiso in one function. In between those calls, you could call your func1, func2, etc. to do their specific tasks.

0 Kudos
Alexey1978
New Contributor I
523 Views

Thaks a lot. The problem is solved.

0 Kudos
ShanmukhS_Intel
Moderator
510 Views

Hi,


Glad to know that your issue is resolved. If you need any additional information, please submit a new question as this thread will no longer be monitored.


Best Regards,

Shanmukh.SS


0 Kudos
Reply