[petsc-users] Zero pivot in LU factorisation
Toby
T.W.Searle at sms.ed.ac.uk
Thu Jul 25 06:30:25 CDT 2013
Hi Matt,
Thanks for the speedy reply. I tried using the options you suggested:
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 -st_pc_factor_shift_type
NONZERO -st_pc_shift_type_amount 1
But I still get the warnings:
Option left: name:-st_pc_factor_shift_type value: NONZERO
Option left: name:-st_pc_shift_type_amount value: 1
I have tried it with just the -st_pc_type jacobi and -st_ksp_view. This
gives me some eigenvalues, but I don't believe them (they do not appear
on the spectrum which I solve using LAPACK where all eigenvalues have a
real part less than 0). The output was very large, but consists of
repetitions of this:
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 -st_pc_type jacobi -st_ksp_view
Generalized eigenproblem stored in file.
Reading COMPLEX matrices from binary files...
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: jacobi
linear system matrix = precond matrix:
Matrix Object: 1 MPI processes
type: seqaij
rows=6000, cols=6000
total: nonzeros=3600, allocated nonzeros=3600
total number of mallocs used during MatSetValues calls =0
using I-node routines: found 4080 nodes, limit used is 5
...
Number of iterations of the method: 189
Number of linear iterations of the method: 1520
Number of requested eigenvalues: 1
Stopping condition: tol=1e-08, maxit=750
Number of converged approximate eigenpairs: 2
k ||Ax-kBx||/||kx||
----------------- ------------------
1388.774454+0.001823 i 0.751726
1388.774441+0.001820 i 0.912665
Thanks again,
Toby
On 25/07/2013 11:26, Matthew Knepley wrote:
> On Thu, Jul 25, 2013 at 4:38 AM, Toby <T.W.Searle at sms.ed.ac.uk
> <mailto:T.W.Searle at sms.ed.ac.uk>> wrote:
>
> Dear all,
>
> I have a generalised matrix problem: A x = lambda B x. My matrix B
> is diagonal and positive semi-definite. My matrix A is a
> non-hermitian complex matrix.
>
> 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.
>
>
>
> 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 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
>
>
> 1) You probably need the correct prefix for this options, e.g.
> -st_pc_factor_shift_type NONZERO
>
> 2) We would like to see the output of -st_ksp_view, but you probably
> need -st_pc_type jacobi for it to finish
>
> Matt
>
> 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 this!
>
> Kind Regards,
> Toby
>
>
> 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$
>
>
> --
> The University of Edinburgh is a charitable body, registered in
> Scotland, with registration number SC005336.
>
>
>
>
> --
> What most experimenters take for granted before they begin their
> experiments is infinitely more interesting than any results to which
> their experiments lead.
> -- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20130725/fd116ee1/attachment.html>
-------------- next part --------------
An embedded and charset-unspecified text was scrubbed...
Name: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20130725/fd116ee1/attachment.ksh>
More information about the petsc-users
mailing list