Lapack and Pardiso results are not the same.

[michaelmoc]

Hello!
I have a problem with Pardiso set up. Previously to solve my linear system i have used two lapack procedures dgetrf for factorization and dgetrs for solving system. Now i want to use Pardiso i have compressed my matrix and run solver over it, but results from Pardiso are different than this one obtained from (i have used default settings of Pardiso) does any one know maybe how should i set up Pardiso to obtain similar results ?

I want to add, that i have zeroes on diagonal and i am using following configuration:
iparm[0] = 0;
iparm[1] = 1;
iparm[2] = 1;
iparm[3] = 0;
iparm[4] = 0;
iparm[5] = 0;
iparm[6] = 0;
iparm[7] = 20;
iparm[8] = 0;
iparm[9] = 8;
iparm[10] = 1;
iparm[11] = 0;
iparm[12] = 1;
iparm[13] = 0;
iparm[14] = 0;
iparm[15] = 0;
iparm[16] = 0;
iparm[17] = -1;
iparm[18] = -1;
iparm[19] = 0;
iparm[21] = 0;
iparm[26] =1;
I am running PARDISO in three phases 11,22,33.

Best regards, Michael.

[Todd_R_Intel]

Hi Michael,

It's a little unclear to me the nature of the problem that you are trying to solve, but I see a few inconsistencies in your settings I believe.

• You set iparm[0] to 0 (iparm(1) in the reference manual or use this new PARDISO parameter table) to use default values, but the proceed to set the rest.
• You set iparm[1] to 1, but the default value is 2. Which value would you like?
• You set iparm[7] to 20, but it has a default value of 0.
• You set iparm[21] (which is undefined in the manual), but do not set iparm[20] to the default value of 1

Some of these might mean PARDISO is not doing what you wanted expected or wanted it to do.

Do you have a simple test case you can provide?

Todd

[michaelmoc]

Hi Todd
I pasted wrong initialization :/. I have removed all assignment to iparm and i set iparm[0](iparm(1)) to 0.
Then i want to solve a sparse matrix (which have some 0 values on diagonal) of size 80x80 (its just small sample). The matrix origins are from Finite Elements Method. Matrix type is real unsymetric, mtype=11 (probably is symetric but i am compressing it to unsymetric form). So now i have standard settings of PARDISO and following results (LPK - are results from Lapack and PRD -from PARDISO (sorry for that:) ):

LPK[0]: 0.786590 PAR[0]: -0.411848
LPK[1]: -0.021839 PAR[1]: 1.322760
LPK[2]: -0.056615 PAR[2]: -2.966042
LPK[3]: -0.391792 PAR[3]: 0.257903
LPK[4]: 0.009292 PAR[4]: -1.766065
LPK[5]: 0.038043 PAR[5]: 0.676759
LPK[6]: -1.685003 PAR[6]: -3.789905
LPK[7]: -0.890598 PAR[7]: -2.279749
LPK[8]: 0.410337 PAR[8]: 1.222382
LPK[9]: 1.202881 PAR[9]: 5.315492
LPK[10]: -1.635822 PAR[10]: 9.369398
LPK[11]: 1.577946 PAR[11]: 2.168596
LPK[12]: 1.026373 PAR[12]: -6.362875
LPK[13]: 0.846750 PAR[13]: 6.724189
LPK[14]: 0.448387 PAR[14]: 3.193257
LPK[15]: -0.207959 PAR[15]: -1.647199
LPK[16]: -0.614929 PAR[16]: -6.553944
LPK[17]: 0.792811 PAR[17]: -4.023097
LPK[18]: -0.793406 PAR[18]: -4.861817
LPK[19]: -0.498558 PAR[19]: 3.358251
LPK[20]: -0.096015 PAR[20]: 0.030691
LPK[21]: 0.001481 PAR[21]: -0.882419
LPK[22]: 0.009519 PAR[22]: -0.046186
LPK[23]: 0.213747 PAR[23]: 3.493633
LPK[24]: 0.111353 PAR[24]: 1.729788
LPK[25]: -0.051658 PAR[25]: -0.930219
LPK[26]: -0.157391 PAR[26]: -3.200128
LPK[27]: 0.196809 PAR[27]: -0.732047
LPK[28]: -0.202826 PAR[28]: -2.594379
LPK[29]: -0.125007 PAR[29]: 0.850935
LPK[30]: 0.073628 PAR[30]: -0.100059
LPK[31]: -0.000673 PAR[31]: 0.056713
LPK[32]: -0.024724 PAR[32]: -0.245212
LPK[33]: -0.160592 PAR[33]: -0.139697
LPK[34]: -0.089182 PAR[34]: 0.041665
LPK[35]: 0.043562 PAR[35]: -0.064295
LPK[36]: 0.130618 PAR[36]: 0.219946
LPK[37]: -0.104556 PAR[37]: 1.017802
LPK[38]: 0.147956 PAR[38]: 0.049902
LPK[39]: 0.065511 PAR[39]: -0.705388
LPK[40]: 0.786590 PAR[40]: -0.411848
LPK[41]: -0.056615 PAR[41]: -2.966042
LPK[42]: -0.021839 PAR[42]: 1.322760
LPK[43]: -0.391792 PAR[43]: 0.257903
LPK[44]: 0.038043 PAR[44]: 0.676759
LPK[45]: 0.009292 PAR[45]: -1.766065
LPK[46]: -1.685003 PAR[46]: -3.789905
LPK[47]: -1.635822 PAR[47]: 9.369398
LPK[48]: 1.026373 PAR[48]: -6.362875
LPK[49]: 1.577946 PAR[49]: 2.168596
LPK[50]: -0.890598 PAR[50]: -2.279749
LPK[51]: 1.202881 PAR[51]: 5.315492
LPK[52]: 0.410337 PAR[52]: 1.222382
LPK[53]: 0.846750 PAR[53]: 6.724189
LPK[54]: 0.792811 PAR[54]: -4.023097
LPK[55]: -0.498558 PAR[55]: 3.358251
LPK[56]: -0.793406 PAR[56]: -4.861817
LPK[57]: 0.448387 PAR[57]: 3.193257
LPK[58]: -0.614929 PAR[58]: -6.553944
LPK[59]: -0.207959 PAR[59]: -1.647199
LPK[60]: -0.096015 PAR[60]: 0.030691
LPK[61]: 0.009519 PAR[61]: -0.046186
LPK[62]: 0.001481 PAR[62]: -0.882419
LPK[63]: 0.213747 PAR[63]: 3.493633
LPK[64]: 0.196809 PAR[64]: -0.732047
LPK[65]: -0.125007 PAR[65]: 0.850935
LPK[66]: -0.202826 PAR[66]: -2.594379
LPK[67]: 0.111353 PAR[67]: 1.729788
LPK[68]: -0.157391 PAR[68]: -3.200128
LPK[69]: -0.051658 PAR[69]: -0.930219
LPK[70]: 0.073628 PAR[70]: -0.100059
LPK[71]: -0.024724 PAR[71]: -0.245212
LPK[72]: -0.000673 PAR[72]: 0.056713
LPK[73]: -0.160592 PAR[73]: -0.139697
LPK[74]: -0.104556 PAR[74]: 1.017802
LPK[75]: 0.065511 PAR[75]: -0.705388
LPK[76]: 0.147956 PAR[76]: 0.049902
LPK[77]: -0.089182 PAR[77]: 0.041665
LPK[78]: 0.130618 PAR[78]: 0.219946
LPK[79]: 0.043562 PAR[79]: -0.064295

As you can see the diff is quite big, i have tried many options of PARDISO, but is always resulting in identical solution.

Best regards, Michael.

[Gennady_F_Intel]

Michael,

Could you please check the relativeresidualof the solution or please give us the original matrix you are using.

I will check the the problem on our side.

--Gennady

[michaelmoc]

Hello Gennady!
It tooked me a little bit, because of easters :).
I have taken two snapshots of this matrixes, one with dense matrix and second with sparse matrix for PARDISO. I will also try with smaler examples. Thank for your help.

Best regards, Michael.

[Gennady_F_Intel]

Michael,how many RHS in this case?

--Gennady

[michaelmoc]

It is one RHS. In both files it is written under RHS tag.

Best regards, Michael.

[michaelmoc]

I was able to run PARDISO and Lapack on small example from MKL manual.
PARDISO was running with default values, results was the same as results from Lapack.

[michaelmoc]

Hello Gennady,
Where you able to look into my problem? If you want i can write a simple application in which i will hardcode mentioned matrixes and i will resolve them using both lapack and PARDISO to compare results.

Best regards, Michael.

[Gennady_F_Intel]

Michael, apologies - I was not able to get to this problem. If you give the hardcoded example - I will look into it asap.

--Gennady

[michaelmoc]

Hello Gennady,
I have written a simple example which resolve 80x80 matrix (which is hardocded). Compression is performed via mkl_ddnscsr method.Than PARDISO and Lapack procedures are run, yet still i am not able to obtain even similar results. I really appreciate if you could help me out with this.

Best regards, Michael.

[Gennady_F_Intel]

Michael,

ooking at your code, I cannot see how you initializing the iparm[64] arrays.

I can see only

iparm[0] = 1;

iparm[2] = 1;

how about others?

tips: for the first step, lease try to setiparm[0]= 0, and

enable matrix checker iparm[26] = 1,

then check the error status after each calculation phase.

--Gennady

[michaelmoc]

Hello Gennady,
It was my mistake (i was experimenting with code and i must give you later version). Change iparm[0] to 0 and run it. I have done many combinations of iparm configuration and yet i was not able to obtain similar results. If you need additional code or smaller code example i will provide it (It is quite urgent for me) :). Also i have noticed that iparm(7) after execution of 33 (with default parameters) is 0 (no iterative refinement performed?). Question - iparm(14) - reports number of non successfull pivots performed or number of all performed pivots ?
Thank You in advance for your help.

Best regards, Michael.

[Gennady_F_Intel]

Michael,

did you enable matrix checker?

--Gennady

[michaelmoc]

Hi,
Yes i have enabled checker (I dont think that matrix is wrongly compressed, because i have used MKL for compression - you can see it in code).

Regards, Michael.

[michaelmoc]

Hello Gennady,
It seems that i have resolved my issue (it is not problem with PARDISO). I will describe the problem later when i will perform some tests on my code. Sorry fro problem.

Best regards, Michael.

[Gennady_F_Intel]

Hello Michael,

not a problem! Anyway, thanks for the update. Please let us know the reason of this problem.

--Gennady

[michaelmoc]

Hello Gennady
I my application i can't obtain directly whole matrix so i am assembling it in few steps and than i am passing my CSR matrix to PARDISO (and as you know, results was not as supposed to be). Wrong results was, because i compressed columnwise stored matrix as rowwise matrix. Due to that i obtained a matrix which was transposed to the original :/. It was the source of problem. Later i run Lapack with transposition option enabled and it work out that the results are the same.

Thank Your for your help ans sorry for problem :).

Best regards, Michael.

[slavaua]

Thought I would bring this thread up instead of creating new one.

I have 3 programs to solve linear equations system - Pardiso, Lapack and one written with OpenMP by hand.

2nd and 3rd always give identical results, but Pardiso results are same with 2nd and 3rd programs only when system is not big, i.e. 100 equations or so.

After 200 there is huge difference in results, example:

Pardiso:

[bash]x[190] = 0.148965665581836 
x[191] = 0.070085873857829 
x[192] = -0.081967888801244 
x[193] = -0.004636397132979 
x[194] = 0.144152764617072 
x[195] = 0.068992992279192 
x[196] = -0.075972409973426 
x[197] = -0.002248602377340 
x[198] = 0.139676431782044 
x[199] = 0.067989446179749 
x[200] = -0.043945941466850[/bash]

Lapack:

[bash]x[190] = 0.000054015242395 
x[191] = -0.000008779165926 
x[192] = 0.000032009486079 
x[193] = 0.000074904450179 
x[194] = 0.000016875289852 
x[195] = -0.000043032478423 
x[196] = -0.000004206083729 
x[197] = 0.000036651462645 
x[198] = -0.000018741413415 
x[199] = -0.000075950902342 
x[200] = -0.000038955173214[/bash]

Hand Made:

[bash]x[190] = 0.000054015242395 
x[191] = -0.000008779165926 
x[192] = 0.000032009486079 
x[193] = 0.000074904450179 
x[194] = 0.000016875289852 
x[195] = -0.000043032478423 
x[196] = -0.000004206083729 
x[197] = 0.000036651462645 
x[198] = -0.000018741413415 
x[199] = -0.000075950902342 
x[200] = -0.000038955173214[/bash]

My Pardiso settings are:

[bash]        for (i = 0; i < 64; i++) {
                iparm = 0;
        }
        iparm[0] = 1; /* No solver default */
        iparm[1] = 2; /* Fill-in reordering from METIS */
        /* Numbers of processors, value of OMP_NUM_THREADS */
        iparm[2] = 2;
        iparm[3] = 0; /* No iterative-direct algorithm */
        iparm[4] = 0; /* No user fill-in reducing permutation */
        iparm[5] = 0; /* Write solution into x */
        iparm[6] = 0; /* Not in use */
        iparm[7] = 2; /* Max numbers of iterative refinement steps */
        iparm[8] = 0; /* Not in use */
        iparm[9] = 13; /* Perturb the pivot elements with 1E-13 */
        iparm[10] = 1; /* Use nonsymmetric permutation and scaling MPS */
        iparm[11] = 0; /* Not in use */
        iparm[12] = 0; /* Maximum weighted matching algorithm is switched-on (default for non-symmetric) */
        iparm[13] = 0; /* Output: Number of perturbed pivots */
        iparm[14] = 0; /* Not in use */
        iparm[15] = 0; /* Not in use */
        iparm[16] = 0; /* Not in use */
        iparm[17] = -1; /* Output: Number of nonzeros in the factor LU */
        iparm[18] = -1; /* Output: Mflops for LU factorization */
        iparm[19] = 0; /* Output: Numbers of CG Iterations */
        maxfct = 1; /* Maximum number of numerical factorizations. */
        mnum = 1; /* Which factorization to use. */
        msglvl = 0; /* Print statistical information in file */
        error = 0; /* Initialize error flag */
        for (i = 0; i < 64; i++) {
                pt = 0;
        }
[/bash]

Which is pretty similar to the ones in example.

Once again, results are same in all 3 programs when system has 100 equations or under.

[Gennady_F_Intel]

Slava, what is the matrix type?

did you check the condition number of this matrix?

--Gennady