Pardiso Incorrect Implementation of Symmetric Indefinite Systems

[blamm]

First, the documentation for Pardiso is paltry, at best. In fact, it's quite incomplete, and not just the MKL documentation of Pardiso, but the original documentation put together by the authors is also quite bad.

But, since I'm using the Intel MKL version of Pardiso, I'll just stick to that, so whenever I refer to "documentation", read "Intel MKL documentation" please.

First, the doc reads that unsymmetric systems are not yet supported. I fed pardiso a symmetric matrix but set mtype to unsymmetric, and it produced no input errors and correct solution. So, unless I'm missing something, MKL Pardiso DOES support unsymmetric system matrix.

Now for the not so good news, in my opinion anyway.

MKL pardiso apparently does not support symmetric indefinite systems. I fed a symmetric matrix to pardiso, set mtype to symmetric indefinite, and it failed. It gave no indication of input errors (another matter I'll get to in a moment). Pardisoallegedy supports symmetric indefinite systems, so it should also support solution of symmetric definite systems using an (the) unsymmetric algorithm. Going from encompassing set to subsets, we have symmetric indefinite => symmetric semi-definite => symmetric positive definite, according as the eigenvalues are less and/or equal and/or greater than zero, equal to and/or greater than zero, or greater than zero, respectively.

So, is it MKL's pardiso, or pardiso in general, that doesn't correctly support symmetric indefinite systems? It's one or the other, certainly not neither.

It would be somewhat naive to expect a user to know every time pardiso is called whether or not the system is indefinite or positive definite. I attemp to use pardiso as the linear solver for structural analysis finite element equations, and sometimes during the analysis the system could be symmetric indefinite ( rigid body modes manifesting ) and sometimes symmetric positive definite (no rigid body modes manifesting). I haven't yet tried to"fool" pardiso that a symmetric indefinite system is just an unsymmetric system, maybe that could work, but then I have to use much more memory since the entire non-zero system has to be made available. I suppose I could code to determine the entire system's eigenvalues to determine if the system was indefinite or positive definite, but then it would be a small step to go ahead and determine the eigenvectors as well, and then just perform a modal decomposition which would cast the equations in completely uncoupled form, and then I would have no need of a linear solver such as pardiso.

The point is symmetric indefinite systems are not correctly handled, at least in MKL pardiso, as far as I can determine.

Next, the fact that an error condition "inconsistent input" exists begs THAT THE POSSIBLE INCONSISTENT INPUT OUGHT TO BE PUBLISHED!

In all my experiments to prove pardiso does not correctly handle symmetric indefinite systems, I output whether or not pardiso reported an error in inconcsistent input, and never received a positive on that. So incorrect handling by pardiso of symmetric indefinite systems was not because of incorrect input on my part.

So, Intel needs to check this out, and also properly update their documentation of pardiso since MIKL pardiso does indeed handle unsymmetric systems ( at least pardiso handles symmetric systems that it's told are unsymmetric ), and also correct mis-handling of systems that pardiso is told are symmetric indefinite.

Just "feed" MKL pardiso a symmetric system and tell it the system is symmetric indefinite and it will fail. The fact that MKL pardiso handles a symmetric system correctly when it's told that it's actually an unsymmetric system I think proves MKL pardiso is ready for unsymmetric systems.

Someone expert in the implementation of pardiso needs to explicitly indicate to all users precisely what constitutes inconsistent input, because, obvioulsy, pardiso is coded to determine if input is inconsistent, SO WHY NOT LET THE USERS IN ON THAT INFORMATION? It would save a GREAT DEAL OF HEARTBURN, I'm sure, for a lot of us users.

Brian Earl Lamm

Message Edited by blamm on 08-08-2004 10:24 AM

Message Edited by blamm on 08-08-2004 10:26 AM

Message Edited by blamm on 08-08-2004 11:04 AM

[blamm]

Thanks so much for the overwhelming amount responses!

Good Day.

Brian Earl Lamm

[basel]

Brian

the solver PARDISO Version 1.2.2. works fine and it solves symmetric indefinite systems very fast and accurately. Please look at:

http://www.computational.unibas.ch/cs/scicomp/publications/oschenk-etna.pdf

your information on the algorithm and results for a lot of highly
indefinite systems

and

ftp://ftp.numerical.rl.ac.uk/pub/talks/sct.cerfacs04.pdf

for solver performance comparison with e.g. WSMP, UMFPACK, HSL Solvers, Boeing Solvers.

The Intel MKL PARDISO solver solves indefinite systems without pivoting also very fast, but if you really need advanced pivot than you might consider other options.

Olaf