I tried with -esp_type lapack, but it still tries to do LU
I attached the code without perturbation of B,
probably something is missing in my way of setting up the problem.
I also attached the eps_view with perturbation
Thanks
Eugenio
Nitsche_genEigen -esp_type lapack
[0]PETSC ERROR: --------------------- Error Message
--------------------------------------------------------------
[0]PETSC ERROR: Zero pivot in LU factorization:
https://www.mcs.anl.gov/petsc/documentation/faq.html#zeropivot
[0]PETSC ERROR: Zero pivot row 1 value 0. tolerance 2.22045e-14
[0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html for
trouble shooting.
[0]PETSC ERROR: Petsc Development GIT revision: v3.11.3-2263-gce77f2ed1a GIT
Date: 2019-09-26 13:31:14 -0500
[0]PETSC ERROR: Nitsche_genEigen on a arch-linux2-c-opt named linux-8biu by
eaulisa Wed Jul 8 10:09:31 2020
[0]PETSC ERROR: Configure options --with-debugging=0 --with-x=1 COPTFLAGS="-O3
-march=native -mtune=native" CXXOPTFLAGS="-O3 -march=native -mtune=native"
FOPTFLAGS="-O3 -march=native -mtune=native" --download-openmpi=1
--download-fblaslapack=1 --download-hdf5=1 --download-metis=1
--download-parmetis=1 --with-shared-libraries=1 --download-blacs=1
--download-scalapack=1 --download-mumps=1 --download-suitesparse
[0]PETSC ERROR: #1 MatPivotCheck_none() line 731 in
/home/eaulisa/software/petsc/include/petsc/private/matimpl.h
[0]PETSC ERROR: #2 MatPivotCheck() line 748 in
/home/eaulisa/software/petsc/include/petsc/private/matimpl.h
[0]PETSC ERROR: #3 MatLUFactorNumeric_SeqAIJ_Inode() line 1676 in
/home/eaulisa/software/petsc/src/mat/impls/aij/seq/inode.c
[0]PETSC ERROR: #4 MatLUFactorNumeric() line 3056 in
/home/eaulisa/software/petsc/src/mat/interface/matrix.c
[0]PETSC ERROR: #5 PCSetUp_LU() line 126 in
/home/eaulisa/software/petsc/src/ksp/pc/impls/factor/lu/lu.c
[0]PETSC ERROR: #6 PCSetUp() line 894 in
/home/eaulisa/software/petsc/src/ksp/pc/interface/precon.c
[0]PETSC ERROR: #7 KSPSetUp() line 377 in
/home/eaulisa/software/petsc/src/ksp/ksp/interface/itfunc.c
[0]PETSC ERROR: #8 STSetUp_Shift() line 119 in
/home/eaulisa/software/slepc/src/sys/classes/st/impls/shift/shift.c
[0]PETSC ERROR: #9 STSetUp() line 271 in
/home/eaulisa/software/slepc/src/sys/classes/st/interface/stsolve.c
[0]PETSC ERROR: #10 EPSSetUp() line 273 in
/home/eaulisa/software/slepc/src/eps/interface/epssetup.c
[0]PETSC ERROR: #11 EPSSolve() line 136 in
/home/eaulisa/software/slepc/src/eps/interface/epssolve.c
[0]PETSC ERROR: --------------------- Error Message
--------------------------------------------------------------
[0]PETSC ERROR: Argument out of range
[0]PETSC ERROR: Argument 2 out of range
[0]PETSC ERROR: See https://www.mcs.anl.gov/petsc/documentation/faq.html for
trouble shooting.
[0]PETSC ERROR: Petsc Development GIT revision: v3.11.3-2263-gce77f2ed1a GIT
Date: 2019-09-26 13:31:14 -0500
[0]PETSC ERROR: Nitsche_genEigen on a arch-linux2-c-opt named linux-8biu by
eaulisa Wed Jul 8 10:09:31 2020
[0]PETSC ERROR: Configure options --with-debugging=0 --with-x=1 COPTFLAGS="-O3
-march=native -mtune=native" CXXOPTFLAGS="-O3 -march=native -mtune=native"
FOPTFLAGS="-O3 -march=native -mtune=native" --download-openmpi=1
--download-fblaslapack=1 --download-hdf5=1 --download-metis=1
--download-parmetis=1 --with-shared-libraries=1 --download-blacs=1
--download-scalapack=1 --download-mumps=1 --download-suitesparse
[0]PETSC ERROR: #12 EPSGetEigenpair() line 399 in
/home/eaulisa/software/slepc/src/eps/interface/epssolve.c
2.96439e-323
Eugenio Aulisa
Department of Mathematics and Statistics,
Texas Tech University
Lubbock TX, 79409-1042
room: 226
http://www.math.ttu.edu/~eaulisa/
phone: (806) 834-6684
fax: (806) 742-1112
________________________________
From: Jose E. Roman <jro...@dsic.upv.es>
Sent: Wednesday, July 8, 2020 10:07 AM
To: Aulisa, Eugenio <eugenio.aul...@ttu.edu>
Cc: petsc-users@mcs.anl.gov <petsc-users@mcs.anl.gov>; Kara, Erdi
<erdi.k...@ttu.edu>
Subject: Re: [petsc-users] Help with generalized eigenvalue problem
Run with -eps_type lapack
Eigen is not interfaced from SLEPc.
Jose
> El 8 jul 2020, a las 17:06, Aulisa, Eugenio <eugenio.aul...@ttu.edu> escribió:
>
> yes the size varies from 6 to 100,
>
> I will try EPSLAPACK or EIGEN to see how it works
>
> Thanks again
> Eugenio
>
>
>
>
>
> Eugenio Aulisa
>
> Department of Mathematics and Statistics,
> Texas Tech University
> Lubbock TX, 79409-1042
> room: 226
> http://www.math.ttu.edu/~eaulisa/
> phone: (806) 834-6684
> fax: (806) 742-1112
>
>
>
>
> From: Jose E. Roman <jro...@dsic.upv.es>
> Sent: Wednesday, July 8, 2020 10:01 AM
> To: Aulisa, Eugenio <eugenio.aul...@ttu.edu>
> Cc: petsc-users@mcs.anl.gov <petsc-users@mcs.anl.gov>; Kara, Erdi
> <erdi.k...@ttu.edu>
> Subject: Re: [petsc-users] Help with generalized eigenvalue problem
>
> In a standard symmetric eigenproblem, a small perturbation of the matrix
> results in a small perturbation of the eigenvalues. In generalized
> eigenproblems, especially with singular matrices, I don't think perturbation
> theory guarantees this in all cases.
>
> What do you mean by small problems? Do you have millions of problems of size
> 6? In that case you should use EPSLAPACK. It will not factor the matrix (if
> solved as a non-symmetric problem).
>
> Jose
>
>
> > El 8 jul 2020, a las 15:55, Aulisa, Eugenio <eugenio.aul...@ttu.edu>
> > escribió:
> >
> > Dear Jose,
> >
> > Thank you for your answer.
> >
> > I tried MLU but it did not help.
> >
> > I also tried to eliminate the null space of B, but I do not know it apriori,
> > and at least to my knowledge of PETSC-SLEPC there is not an explicit way to
> > find it
> > but to find all the 0 eigenvectors of B, which again it requires the
> > inverse of B.
> >
> > I am solving similar small problems millions of times, and this works fine
> > most of the time,
> > but once in a while it fails with the same LU zero pivot exception.
> >
> > One thing I tried yesterday, and it seams to work all the times,
> > is perturbing the diagonal of B of a small factor, += 1.0e-10 ||B||,
> > but I do not know how theoretically it sounds with all these indeterminate
> > eigenvalues,
> > maybe you have a clue on this.
> >
> > Thanks,
> > Eugenio
> >
> >
> >
> > Eugenio Aulisa
> >
> > Department of Mathematics and Statistics,
> > Texas Tech University
> > Lubbock TX, 79409-1042
> > room: 226
> > http://www.math.ttu.edu/~eaulisa/
> > phone: (806) 834-6684
> > fax: (806) 742-1112
> >
> >
> >
> >
> > From: Jose E. Roman <jro...@dsic.upv.es>
> > Sent: Tuesday, July 7, 2020 11:01 AM
> > To: Aulisa, Eugenio <eugenio.aul...@ttu.edu>
> > Cc: petsc-users@mcs.anl.gov <petsc-users@mcs.anl.gov>; Kara, Erdi
> > <erdi.k...@ttu.edu>
> > Subject: Re: [petsc-users] Help with generalized eigenvalue problem
> >
> > I don't know what options you are using, shift-and-invert or not, which,
> > target. Anyway, if A and B have a common null space, then A-sigma*B will
> > always be singular, so the KSP solver will have to deal with a singular
> > system.
> >
> > If you pass the nullspace vectors to the EPS solver via
> > https://nam04.safelinks.protection.outlook.com/?url=https%3A%2F%2Fslepc.upv.es%2Fdocumentation%2Fcurrent%2Fdocs%2Fmanualpages%2FEPS%2FEPSSetDeflationSpace.html&data=02%7C01%7CEugenio.Aulisa%40ttu.edu%7C48eb040474e44a2af2f208d82350b5a7%7C178a51bf8b2049ffb65556245d5c173c%7C0%7C0%7C637298176878634547&sdata=lL%2B5GTRzFgeviU%2BedhKieyXDO%2Fgyko2CLYGU%2BxtilB0%3D&reserved=0
> > then these vectors will also be attached to the KSP solver.
> >
> > Also, since you have MUMPS, try adding
> > -st_pc_factor_mat_solver_type mumps
> > (see section 3.4.1 of the SLEPc users manual). MUMPS should be quite robust
> > with respect to singular linear systems.
> >
> > Jose
> >
> >
> > > El 7 jul 2020, a las 17:40, Aulisa, Eugenio <eugenio.aul...@ttu.edu>
> > > escribió:
> > >
> > > Can somebody help me figure it out the solver options for the attached
> > > generalized eigenvalue problem?
> > >
> > > A x = lambda B x
> > >
> > > where A and B have the same null space, with 5 zero eigenvalues.
> > >
> > > The only non indefinite eigenvalue should be 41.1892
> > >
> > > The problem is Hermite semi-positive so automatically slepc should purify
> > > it, but to be sure I added
> > >
> > > EPSSetPurify(eps, PETSC_TRUE);
> > >
> > > I guess with my non-solver-options slepc is still trying to do the
> > > inverse of the full
> > > B matrix with an LU decomposition and gets zero pivot.
> > >
> > > This is what I get as a PETSC error Message
> > >
> > > [0]PETSC ERROR: --------------------- Error Message
> > > --------------------------------------------------------------
> > > [0]PETSC ERROR: Zero pivot in LU factorization:
> > > https://nam04.safelinks.protection.outlook.com/?url=https%3A%2F%2Fwww.mcs.anl.gov%2Fpetsc%2Fdocumentation%2Ffaq.html%23zeropivot&data=02%7C01%7CEugenio.Aulisa%40ttu.edu%7C48eb040474e44a2af2f208d82350b5a7%7C178a51bf8b2049ffb65556245d5c173c%7C0%7C0%7C637298176878634547&sdata=MfAJt1HSEK5wvhRzx0BM9OEEGxWjnlI%2FwTlzhN4BDeE%3D&reserved=0
> > > [0]PETSC ERROR: Bad LU factorization
> > > [0]PETSC ERROR: See
> > > https://nam04.safelinks.protection.outlook.com/?url=https%3A%2F%2Fwww.mcs.anl.gov%2Fpetsc%2Fdocumentation%2Ffaq.html&data=02%7C01%7CEugenio.Aulisa%40ttu.edu%7C48eb040474e44a2af2f208d82350b5a7%7C178a51bf8b2049ffb65556245d5c173c%7C0%7C0%7C637298176878634547&sdata=NTc7lCzIo1Vj2sYmAehiE5b8ToXDTvN0wzDyh7JC30w%3D&reserved=0
> > > for trouble shooting.
> > > [0]PETSC ERROR: Petsc Development GIT revision:
> > > v3.11.3-2263-gce77f2ed1a GIT Date: 2019-09-26 13:31:14 -0500
> > > [0]PETSC ERROR: Nitsche_ex4a on a arch-linux2-c-opt named linux-8biu
> > > by eaulisa Tue Jul 7 10:17:55 2020
> > > [0]PETSC ERROR: Configure options --with-debugging=0 --with-x=1
> > > COPTFLAGS="-O3 -march=native -mtune=native" CXXOPTFLAGS="-O3
> > > -march=native -mtune=native" FOPTFLAGS="-O3 -march=native -mtune=native"
> > > --download-openmpi=1 --download-fblaslapack=1 --download-hdf5=1
> > > --download-metis=1 --download-parmetis=1 --with-shared-libraries=1
> > > --download-blacs=1 --download-scalapack=1 --download-mumps=1
> > > --download-suitesparse
> > > [0]PETSC ERROR: #1 MatLUFactor_SeqDense() line 633 in
> > > /home/eaulisa/software/petsc/src/mat/impls/dense/seq/dense.c
> > > [0]PETSC ERROR: #2 MatLUFactorNumeric_SeqDense() line 432 in
> > > /home/eaulisa/software/petsc/src/mat/impls/dense/seq/dense.c
> > > [0]PETSC ERROR: #3 MatLUFactorNumeric() line 3056 in
> > > /home/eaulisa/software/petsc/src/mat/interface/matrix.c
> > > [0]PETSC ERROR: #4 PCSetUp_LU() line 126 in
> > > /home/eaulisa/software/petsc/src/ksp/pc/impls/factor/lu/lu.c
> > > [0]PETSC ERROR: #5 PCSetUp() line 894 in
> > > /home/eaulisa/software/petsc/src/ksp/pc/interface/precon.c
> > > [0]PETSC ERROR: #6 KSPSetUp() line 377 in
> > > /home/eaulisa/software/petsc/src/ksp/ksp/interface/itfunc.c
> > > [0]PETSC ERROR: #7 STSetUp_Shift() line 119 in
> > > /home/eaulisa/software/slepc/src/sys/classes/st/impls/shift/shift.c
> > > [0]PETSC ERROR: #8 STSetUp() line 271 in
> > > /home/eaulisa/software/slepc/src/sys/classes/st/interface/stsolve.c
> > > [0]PETSC ERROR: #9 EPSSetUp() line 273 in
> > > /home/eaulisa/software/slepc/src/eps/interface/epssetup.c
> > > [0]PETSC ERROR: #10 EPSSolve() line 136 in
> > > /home/eaulisa/software/slepc/src/eps/interface/epssolve.c
> > >
> > >
> > >
> > > Thanks
> > > Eugenio
> > >
> > >
> > >
> > >
> > > Eugenio Aulisa
> > >
> > > Department of Mathematics and Statistics,
> > > Texas Tech University
> > > Lubbock TX, 79409-1042
> > > room: 226
> > > http://www.math.ttu.edu/~eaulisa/
> > > phone: (806) 834-6684
> > > fax: (806) 742-1112
> > >
> > >
> > >
> > > <genEigen.cpp>
> >
>
#include "slepceps.h"
#include <iostream>
int main(int argc, char** args) {
SlepcInitialize(&argc, &args, PETSC_NULL, PETSC_NULL);
double a[6][6] = {{0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
{0.0, 1.9247225877689356, -1.9247225877689356, -0.96236129388446778, 0.0, 0.96236129388446778 },
{0.0, -1.9247225877689356, 1.9247225877689356, 0.96236129388446778, 0.0, -0.96236129388446778},
{0.0, -0.96236129388446778, 0.96236129388446778, 0.48118064694223389, 0.0, -0.48118064694223389},
{0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
{0.0, 0.96236129388446778, -0.96236129388446778, -0.48118064694223389, 0.0, 0.48118064694223389}
};
double b[6][6] = {{0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
{0.0, 0.046728852764308028, -0.046728852764308028, -0.023364426382154014, 0.0, 0.023364426382154014},
{0.0, -0.046728852764308028, 0.046728852764308028, 0.023364426382154014, 0.0, -0.023364426382154014 },
{0.0, -0.023364426382154014, 0.023364426382154014, 0.011682213191077007, 0.0, -0.011682213191077007 },
{0.0, 0.0, 0.0, 0.0, 0.0, 0.0},
{0.0, 0.023364426382154014, -0.023364426382154014, -0.011682213191077007, 0.0, 0.011682213191077007 }
};
Mat A, B;
EPS eps;
int idx[6] = {0, 1, 2, 3, 4, 5};
MatCreateAIJ(PETSC_COMM_SELF,6,6,6,6,6,PETSC_NULL,0,PETSC_NULL,&A);
MatCreateAIJ(PETSC_COMM_SELF,6,6,6,6,6,PETSC_NULL,0,PETSC_NULL,&B);
MatSetValues(A, 6, &idx[0], 6, &idx[0], &a[0][0], INSERT_VALUES);
MatSetValues(B, 6, &idx[0], 6, &idx[0], &b[0][0], INSERT_VALUES);
MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
MatAssemblyBegin(B, MAT_FINAL_ASSEMBLY);
MatAssemblyEnd(B, MAT_FINAL_ASSEMBLY);
double nrm;
MatNorm(B, NORM_INFINITY, &nrm);
//MatShift(B, 1.0e-10 * nrm);
EPSCreate(PETSC_COMM_SELF, &eps);
EPSSetOperators(eps, A, B);
EPSSetWhichEigenpairs(eps, EPS_LARGEST_MAGNITUDE);
EPSSetPurify(eps, PETSC_TRUE);
EPSSetFromOptions(eps);
EPSSolve(eps);
double real;
EPSGetEigenpair(eps, 0, &real, PETSC_NULL, PETSC_NULL, PETSC_NULL);
std::cout << real << " " << std::endl;
EPSDestroy(&eps);
MatDestroy(&A);
MatDestroy(&B);
return 1;
}
EPS Object: 1 MPI processes
type: krylovschur
50% of basis vectors kept after restart
using the locking variant
problem type: generalized non-symmetric eigenvalue problem
selected portion of the spectrum: largest eigenvalues in magnitude
number of eigenvalues (nev): 1
number of column vectors (ncv): 6
maximum dimension of projected problem (mpd): 6
maximum number of iterations: 100
tolerance: 1e-08
convergence test: relative to the eigenvalue
BV Object: 1 MPI processes
type: svec
7 columns of global length 6
vector orthogonalization method: classical Gram-Schmidt
orthogonalization refinement: if needed (eta: 0.7071)
block orthogonalization method: GS
doing matmult as a single matrix-matrix product
DS Object: 1 MPI processes
type: nhep
ST Object: 1 MPI processes
type: shift
shift: 0.
number of matrices: 2
all matrices have different nonzero pattern
KSP Object: (st_) 1 MPI processes
type: preonly
maximum iterations=10000, initial guess is zero
tolerances: relative=1e-08, absolute=1e-50, divergence=10000.
left preconditioning
using NONE norm type for convergence test
PC Object: (st_) 1 MPI processes
type: lu
out-of-place factorization
tolerance for zero pivot 2.22045e-14
matrix ordering: nd
factor fill ratio given 5., needed 1.
Factored matrix follows:
Mat Object: 1 MPI processes
type: seqaij
rows=6, cols=6
package used to perform factorization: petsc
total: nonzeros=36, allocated nonzeros=36
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 2 nodes, limit used is 5
linear system matrix = precond matrix:
Mat Object: 1 MPI processes
type: seqaij
rows=6, cols=6
total: nonzeros=36, allocated nonzeros=36
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 2 nodes, limit used is 5
41.1892
