[petsc-users] Incorrect Eigenvalues when Setting KSP and PC types

[Stefano Zampini]

Hong,

I personally believe that commit https://bitbucket.org/petsc/petsc/commits/
966c94c9cf4fa13d455726ec36800a577e00b171 should be reverted.
I agree on the fact that when the user sets an option (the hermitian one in
this case) and that feature is not supported we should throw an error (
https://bitbucket.org/petsc/petsc/commits/8f21f52c465b775a76cda90fe9c51d0c742078c7)
, but I don't like the fact that the user is forced to set on option to use
a feature that he already requested (as in https://bitbucket.org/petsc/
petsc/commits/966c94c9cf4fa13d455726ec36800a577e00b171).

Barry, what do you think?

2017-09-29 5:28 GMT+03:00 Hong <hzhang at mcs.anl.gov>:

> Greg:
>
>> Thanks so much for the detailed response. I am glad to know the reason
>> behind it--hopefully we eventually figure out why the solvers have this
>> behavior! Hong, I really appreciate you working on a patch to throw an
>> error in this case. It definitely bit me and some people using my code...
>> Hopefully it doesn't happen to anyone else! :)
>>
> I added an error flag for using MUMPS Cholesky factorization on Hermitian
> matrix. See branch hzhang/mumps-HermitianCholesky-errflag
> https://bitbucket.org/petsc/petsc/commits/8f21f52c465b775a76cda90fe9c51d
> 0c742078c7
>
> Jose,
> PETSc does not support Cholesky for Hermitian matrix.
>
>>
>>>> The linear solver table probably needs to be updated. I have tried
>>>> several Cholesky solvers. With mkl_pardiso I get an explicit error message
>>>> that it does not support Cholesky with complex scalars. The rest (PETSc,
>>>> MUMPS, CHOLMOD) give a wrong answer (without error message). The problem is
>>>> not related to your matrix, nor to shift-and-invert in SLEPc. I tried with
>>>> ex1.c under PETSC_DIR/src/ksp/ksp/examples/tutorials. The example
>>>> works in complex scalars, but the matrix is real. As soon as you add
>>>> complex entries Cholesky fails, for instance adding this:
>>>>   ierr = MatSetValue(A,0,1,1.0+PETSC_i,INSERT_VALUES);CHKERRQ(ierr);
>>>>   ierr = MatSetValue(A,1,0,1.0-PETSC_i,INSERT_VALUES);CHKERRQ(ierr);
>>>>
>>>
> In this case, you must call
> MatSetOption(A,MAT_HERMITIAN,PETSC_TRUE);
>
> Then, petsc will throw an error for '-pc_type cholesky'.
>
>>
>>>> I don't know if it is a bug or that the algorithm cannot support
>>>> complex Hermitian matrices. Maybe Hong can confirm any of these. In the
>>>> latter case, I agree that all packages should give an error message, as
>>>> mkl_pardiso does.
>>>>
>>>
> I also add an error flag for Cholesky/ICC if user does not set
> MatSetOption(A,MAT_SYMMETRIC,PETSC_TRUE) for aij matrix.
> See https://bitbucket.org/petsc/petsc/commits/
> 966c94c9cf4fa13d455726ec36800a577e00b171
>
> Let me know if you have any comments about this fix.
>
> Hong
>
>>
>>>> > El 25 sept 2017, a las 7:17, Greg Meyer <gregory.meyer at gmail.com>
>>>> escribió:
>>>> >
>>>> > Hi all,
>>>> >
>>>> > Hong--looking at your link, there may be no special algorithm for
>>>> Hermitian matrices in MUMPS, but that doesn't mean it can't solve them like
>>>> it would any matrix. Furthermore it appears that Cholesky of complex
>>>> matrices is supported from this link: https://www.mcs.anl.gov/petsc/
>>>> documentation/linearsolvertable.html
>>>> >
>>>> > So do you or anyone have any idea why I get incorrect eigenvalues?
>>>> >
>>>> > Thanks,
>>>> > Greg
>>>> >
>>>> > On Thu, Sep 21, 2017 at 5:51 PM Greg Meyer <gregory.meyer at gmail.com>
>>>> wrote:
>>>> > Ok, thanks. It seems that PETSc clearly should throw an error in this
>>>> case instead of just giving incorrect answers? I am surprised that it does
>>>> not throw an error...
>>>> >
>>>> > On Thu, Sep 21, 2017 at 5:24 PM Hong <hzhang at mcs.anl.gov> wrote:
>>>> > Greg :
>>>> > Yes, they are Hermitian.
>>>> >
>>>> > PETSc does not support  Cholesky factorization for Hermitian.
>>>> > It seems mumps does not support Hermitian either
>>>> > https://lists.mcs.anl.gov/mailman/htdig/petsc-users/2015-
>>>> November/027541.html
>>>> >
>>>> > Hong
>>>> >
>>>> >
>>>> > On Thu, Sep 21, 2017 at 3:43 PM Hong <hzhang at mcs.anl.gov> wrote:
>>>> > Greg:
>>>> >
>>>> > OK, the difference is whether LU or Cholesky factorization is used.
>>>> But I would hope that neither one should give incorrect eigenvalues, and
>>>> when I run with the latter it does!
>>>> > Are your matrices symmetric/Hermitian?
>>>> > Hong
>>>> >
>>>> > On Thu, Sep 21, 2017 at 2:05 PM Hong <hzhang at mcs.anl.gov> wrote:
>>>> > Gregory :
>>>> > Use '-eps_view' for both runs to check the algorithms being used.
>>>> > Hong
>>>> >
>>>> > Hi all,
>>>> >
>>>> > I'm using shift-invert with EPS to solve for eigenvalues. I find that
>>>> if I do only
>>>> >
>>>> > ...
>>>> >   ierr = EPSGetST(eps,&st);CHKERRQ(ierr);
>>>> >   ierr = STSetType(st,STSINVERT);CHKERRQ(ierr);
>>>> > ...
>>>> >
>>>> > in my code I get correct eigenvalues. But if I do
>>>> >
>>>> > ...
>>>> >   ierr = EPSGetST(eps,&st);CHKERRQ(ierr);
>>>> >   ierr = STSetType(st,STSINVERT);CHKERRQ(ierr);
>>>> >   ierr = STGetKSP(st,&ksp);CHKERRQ(ierr);
>>>> >   ierr = KSPGetPC(ksp,&pc);CHKERRQ(ierr);
>>>> >   ierr = KSPSetType(ksp,KSPPREONLY);CHKERRQ(ierr);
>>>> >   ierr = PCSetType(pc,PCCHOLESKY);CHKERRQ(ierr);
>>>> > ...
>>>> >
>>>> > the eigenvalues found by EPS are completely wrong! Somehow I thought
>>>> I was supposed to do the latter, from the examples etc, but I guess that
>>>> was not correct? I attach the full piece of test code and a test matrix.
>>>> >
>>>> > Best,
>>>> > Greg
>>>> >
>>>>
>>>>
>>>
>


-- 
Stefano
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20170929/4e5044cd/attachment.html>