Computation of the Schur complement in MKL

[Vasileios_K_]

Hello everyone,

I recently came across a project in which I have to compute the Schur complement of a complex symmetrix matrix. I know that starting from Intel® MKL 11.2 update 1,  MKL supports the computation of Schur complements. However, I have two questions that doesn't seem to have been answered somewhere. Note: I am writing a C code.

a) Does this functionality supports complex arithmetic? I know that the PARDISO implementation in MKL supports complex arithmetic, but does this apply to the Schur complement computation also?

b) It is not clear how to tune the parameters in order to obtain the Schur complement via PARDISO. What I am trying is the following :

PARDISO(pt, &maxfct, &mnum, &mtype, &phase, &sizeSchur, acsr, ia, ja, &idum, &one, iparm, &msglvl, &ddum, SS, &error);

where SS is a sizeSchur * sizeSchur complex vector and sizeSchur is the size of the Schur complement. The routine returns error = -1. Is there anyone who can provide an example of how to set the parameters for the above routine? I have seen the documentation provided in intel's website but it was not very helpful for me. I think the iparm vector is fine, I currently set iparm[35] = 1 (I have tried many combinations but none solves my issue).

I can provide a real example if the above are not clear. Thanks!

[Alexander_K_Intel2]

Hi,

there is a bit misunderstanding here. Just call simple pardiso without changing main parameters like size. You need to set iparm[35]=1 and perm array that must be filled in the following manner: set 1 values for the components that correspond to schur complement matrix and 0 otherwise. 

Thanks,

Alex

[Vasileios_K_]

Hi Alexander,

Thank you very much for your response. I changed the code according to your advice and indeed, now PARDISO returns error = 0. However, the Schur complement returned is just the zero vector. Thus, just as a follow-up question, I would like your advice on my inputs.

I set phase = 12 and set the number of right-hand sides to one. Moreover, I provide a solution vector of size s*s, where s is the size of the Schur complement. Regarding the permutation vector you suggested, I provide a vector "perm" of length n (the size of A) in which all values are set to zero except the last s ones, which I set to one.

All in all, my PARDISO call looks like

mtype = 3, phase = 12, iparm[35] = 1 (iparm has other values set too)

PARDISO(pt, &one, &one, &mtype, &phase, &n, a, ia, ja, perm, &one, iparm, &msgvl, &ddum, SS, &error);

What I actually want to achieve is exactly what takes place in the following link:

https://software.intel.com/en-us/articles/intel-mkl-support-to-new-functionality-schur-complement

If the above look correct, then I probably did some error somewhere before I reach this point.

Thank you so much for your help!

[Alexander_K_Intel2]

Hi,

Everything looks correct. Could you send reproducer of your problem to me to play with it on my side?

Thanks,

Alex

[Vasileios_K_]

Hi,

I attach a small working example for a 2D Laplacian matrix. The matrix is loaded, converted to CSR format, and shifted by a complex shift. Then, I compute and factorize the Schur complement.

[Fanny_D_]

Hi,

what is the status of the bug please?

With the MKL 2016 (11.3 Update 2) I tested the computation of the Schur complement on a complex matrix and I also obtained zero complex values as mentioned above.

Thanks in advance for your help.

[Alexander_K_Intel2]

HI,

I believe that it is passed on last version of MKL but let me check

Thanks,

Alex

[Alexander_K_Intel2]

Hi Fanny,

Ok, I've checked your example - the issue that you see is not caused by schur complement matrix but inpur parameters on second call. You set input matrix as zero based but call pardiso with one based (iparm[34] = 0;) Just set iparm[34] to 1 or set iparm[26] and msglvl to 1 to see matrix checker info. By the way, after change iparm[34] to 1 I've got error equal to 0 on second call

Thanks,

Alex

[Fanny_D_]

Hi Alex,

thanks for your answer.

But I have another example (please find the attached source files) where iparm[34] is well defined but which does not give good results on the Schur complement of a complex matrix.

The source file pardiso_schur_real.cpp returns the right real Schur complement:

-1.800000e+001 0.000000e+000

0.000000e+000 -3.200000e+001

The source file pardiso_schur_complex.cpp (which corresponds to the same matrix than the real case) returns a zero complex Schur complement.

I've tested with MKL 2016 (11.3 Update 2).

Could you please test with the last version and tell me if I made some mistakes?

Thanks.

[Fanny_D_]

Hi Alex,

Did you manage to test my program with the last version of MKL?

Thanks in advance.

[Shirman__Yuri]

I have a case which perfectly runs for float data (11) and returns results close to machine zero for complex data (13).

mtype options:
 

11 -> real and nonsymmetric

13 -> complex and nonsymmetric

I checked this with MKL-2018 and MKL-2019

Has anything changed in the last versions?