# Solving a large non-hermitian generalised eigenvalue problem from a linear stability analysis using SLEPc

I have a generalised matrix problem: $A x = \lambda B x$ from a spectral method on a linear stability analysis problem. My matrix B is diagonal and positive semi-definite. A is non-hermitian and complex.

My problem is essentially that when using the SLEPc generalised eigenvalue solver I get the error "zero pivot in LU factorisation". The rest of the below is details about the problem and things I have tried so far. Thanks for the help!

## Details of the problem

The matrix will be at its largest about 48000 by 48000, and I want to find the eigenvalues. The eigenvalues I am interested in are ones with the largest real part near 0+0i. Ideally, I want to be able to find them even if they are internal (i.e when there are other eigenvalues with larger positive real part in the spectrum). However, I would be happy if I could get it to work for problems where all eigenvalues have real parts < 0 apart from the eigenvalue of interest.

At the moment I have used the scipy linalg.eig and sparse.eigs functions. As far as I know, these use LAPACK and ARPACK respectively to do the heavy lifting. I have decided to see if I can achieve better performance through using the SLEPc library. If this is a bad decision, let me know!

I want to move onto using PETSc with the SLEPc eigenvalue solvers. I have been trying out SLEPc using the examples provided as part of the tutorial. Exercise 7 (http://www.grycap.upv.es/slepc/handson/handson3.html) reads matricies A and B from a file and outputs the solutions. I got this to work fine using the matrices provided. However, if I substitute a smaller sized test version of my problem (6000x6000), I get a variety of errors depending on the command line arguments I supply.

The main problem I have is the error: "zero pivot in LU factorisation!" when I use the default settings.

I think this might be related to the fact that B contains rows of zeros, although my understanding of linear algebra is somewhat basic. Is this true?

I have tried setting the options suggested on the petsc website, -pc_factor_shift_type NONZERO etc but all I get is an additional warning that these options were not used

I assumed that this was a problem with the preconditioner, so I tried setting -eps_target to 0.1 and both with and without specifying -st_type sinvert and shift. Still I get the same error.

Then I tried -st_pc_type jacobi and st_pc_type bjacobi. jacobi runs, but does not produce any eigenvalues. Block jacobi does an LU factorisation and gives me the same error again.

The default method is krylov-schur, so I have experimented with the -eps_type gd and -eps_type jd options. Unfortunately these seem to produce nonsense eigenvalues, which do not appear on the spectrum at all when I solve using LAPACK in scipy.

I know my matrix problem is not singular, because I can solve it using scipy.

Do you know of any books/guides I might need to read besides the PETSC and SLEPC manuals to understand the behaviour of all these different solvers?

The output from the case with no command line options is given below.

Thanks a lot for taking the time to read my first post!

Kind Regards, Toby

## Terminal output from SLEPc

tobymac:SLEPC toby$mpiexec ./ex7 -f1 LHS-N7-M40-Re0.0-b0.1-Wi5.0-amp0.02.petsc -f2 RHS-N7-M40-Re0.0-b0.1-Wi5.0-amp0.02.petsc -eps_view Generalized eigenproblem stored in file. [0]PETSC ERROR: --------------------- Error Message ------------------------------------ [0]PETSC ERROR: Detected zero pivot in LU factorization: see http://www.mcs.anl.gov/petsc/documentation/faq.html#ZeroPivot! [0]PETSC ERROR: Empty row in matrix: row in original ordering 2395 in permuted ordering 3600! [0]PETSC ERROR: ------------------------------------------------------------------------ [0]PETSC ERROR: Petsc Release Version 3.3.0, Patch 5, Sat Dec 1 15:10:41 CST 2012 [0]PETSC ERROR: See docs/changes/index.html for recent updates. [0]PETSC ERROR: See docs/faq.html for hints about trouble shooting. [0]PETSC ERROR: See docs/index.html for manual pages. [0]PETSC ERROR: ------------------------------------------------------------------------ [0]PETSC ERROR: ./ex7 on a arch-darw named tobymac by toby Thu Jul 25 10:20:40 2013 [0]PETSC ERROR: Libraries linked from /opt/local/lib [0]PETSC ERROR: Configure run at Tue Jul 23 15:11:27 2013 [0]PETSC ERROR: Configure options --prefix=/opt/local --with-valgrind-dir=/opt/local --with-shared-libraries --with-scalar-type=complex --with-clanguage=C++ --with-superlu-dir=/opt/local --with-blacs-dir=/opt/local --with-scalapack-dir=/opt/local --with-mumps-dir=/opt/local --with-metis-dir=/opt/local --with-parmetis-dir=/opt/local --COPTFLAGS=-O2 --CXXOPTFLAGS=-O2 --FOPTFLAGS=-O2 --LDFLAGS=-L/opt/local/lib --CFLAGS="-O2 -mtune=native" --CXXFLAGS="-O2 -mtune=native" [0]PETSC ERROR: ------------------------------------------------------------------------ [0]PETSC ERROR: MatLUFactorSymbolic_SeqAIJ() line 334 in /opt/local/var/macports/build/_Users_toby_MyPorts_scienceports_math_petsc/petsc/work/petsc-3.3-p5/src/mat/impls/aij/seq/aijfact.c [0]PETSC ERROR: MatLUFactorSymbolic() line 2750 in /opt/local/var/macports/build/_Users_toby_MyPorts_scienceports_math_petsc/petsc/work/petsc-3.3-p5/src/mat/interface/matrix.c [0]PETSC ERROR: PCSetUp_LU() line 135 in /opt/local/var/macports/build/_Users_toby_MyPorts_scienceports_math_petsc/petsc/work/petsc-3.3-p5/src/ksp/pc/impls/factor/lu/lu.c Number of iterations of the method: 0 Number of linear iterations of the method: 0 Number of requested eigenvalues: 1 Stopping condition: tol=1e-08, maxit=750 Number of converged approximate eigenpairs: 0 [0]PETSC ERROR: PCSetUp() line 832 in /opt/local/var/macports/build/_Users_toby_MyPorts_scienceports_math_petsc/petsc/work/petsc-3.3-p5/src/ksp/pc/interface/precon.c [0]PETSC ERROR: KSPSetUp() line 278 in /opt/local/var/macports/build/_Users_toby_MyPorts_scienceports_math_petsc/petsc/work/petsc-3.3-p5/src/ksp/ksp/interface/itfunc.c [0]PETSC ERROR: PCSetUp_Redundant() line 176 in /opt/local/var/macports/build/_Users_toby_MyPorts_scienceports_math_petsc/petsc/work/petsc-3.3-p5/src/ksp/pc/impls/redundant/redundant.c [0]PETSC ERROR: PCSetUp() line 832 in /opt/local/var/macports/build/_Users_toby_MyPorts_scienceports_math_petsc/petsc/work/petsc-3.3-p5/src/ksp/pc/interface/precon.c [0]PETSC ERROR: KSPSetUp() line 278 in /opt/local/var/macports/build/_Users_toby_MyPorts_scienceports_math_petsc/petsc/work/petsc-3.3-p5/src/ksp/ksp/interface/itfunc.c [0]PETSC ERROR: STSetUp_Shift() line 94 in src/st/impls/shift/shift.c [0]PETSC ERROR: STSetUp() line 280 in src/st/interface/stsolve.c [0]PETSC ERROR: EPSSetUp() line 204 in src/eps/interface/setup.c [0]PETSC ERROR: EPSSolve() line 109 in src/eps/interface/solve.c tobymac:SLEPC toby$

• Think about what a row of all zeros means for a little bit. Write out the equation decide if its helpful. :) – meawoppl Feb 23 '14 at 3:27
• I have a similar problem. But did u try ./ex7 -st_type sinvert .. it works. But I am not getting the correct eigenvector, but m getting eigenvalue. Did u manage to solve ur problem ? – han17 Mar 23 '14 at 12:21

In problems of this form, you typically need that $B$ is non-singular, i.e., that it is invertible. You don't say where $A,B$ come from, so I don't know anything about your $B$ but you need to make sure that it is not singular.
• I believe that scipy is doing $A^{-1}Bx = \lambda^{-1}x$, which ought to be ok? – Toby Searle Jul 25 '13 at 15:55
• Assume $B$ is singular, then there will be vectors $x$ in its kernel so that $A^{-1}Bx=0$, i.e., they are eigenvectors corresponding to $\lambda^{-1}=0$ which clearly presents a problem for the definition of the generalized eigenvalues of the original problem. Consequently, one typically wants $B$ to be non-singular. – Wolfgang Bangerth Jul 25 '13 at 22:14