[petsc-users] Zero pivot in LU factorisation
Toby
T.W.Searle at sms.ed.ac.uk
Thu Jul 25 08:44:27 CDT 2013
Thanks, I have corrected the shift amount option and I get the same
error as when using the -st_pc_factor_shift_amount 1
I have been using the PETSc IO module to convert my matrix from sparse
coordinate format into PETSc binary. I am guessing there is some problem
with this step. I use the function:
PetscBinaryIO.PetscBinaryIO().writeMatSciPy( fh, matrix)
to convert into Petsc binary format. Is there something wrong with this?
Do I need to manually put zeros into the diagonal elements of zero rows
of B?
When I solve the problem using scipy I use dense matrices. There are
definitely entries in the matrix A at row 2395. I used python to confirm
this before I convert it to petsc format:
column | value
LHS: 1795 10.9955742876j
1198 (76+0j)
1196 (72+0j)
595 0.02j
However, the matrix is all zeros at row 2395 in B. As far as I
understand, this doesn't make my problem singular. I am sure that the
problem is not rank deficient when I solve it in scipy.
Thanks again,
Toby
On 25/07/2013 12:51, Matthew Knepley wrote:
> On Thu, Jul 25, 2013 at 6:46 AM, Toby <T.W.Searle at sms.ed.ac.uk
> <mailto:T.W.Searle at sms.ed.ac.uk>> wrote:
>
>
> Ok, it looks like the -st_pc_factor_shift_type option is now being
> used, thanks. Unfortunately it hasn't fixed the problem.
>
> if I use the NONZERO type and a shift amount, the shift amount
> option is not recognised and I still get the zero pivot error. If
> I use the POSITIVE_DEFINITE type I also get the pivot error.
> Output is below
>
> Thanks,
> 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 -st_pc_type lu
> -st_pc_factor_shift_type NONZERO -st_pc_shift_type_amount 1
>
> Generalized eigenproblem stored in file.
>
> Reading COMPLEX matrices from binary files...
> [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!
>
>
> 1) You misspelled the option: -st_pc_factor_shift_amount 1
>
> 2) It was used, and fixed you problem on row 0, now you have a problem
> on row 2395, namely that it is missing. No
> LU factorization routine will succeed here. You must at least put
> a diagonal 0.0 on all rows.
>
> 3) Same goes for the run below.
>
> 4) A completely empty row means your matrix is indeed rank deficient.
> I have no idea what results are coming out of
> scipy, but they are not using factorization.
>
> Matt
>
> [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 12:41:57 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
> [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
> 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
>
> WARNING! There are options you set that were not used!
> WARNING! could be spelling mistake, etc!
> Option left: name:-st_pc_shift_type_amount value: 1
>
> 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 lu
> -st_pc_factor_shift_type POSITIVE_DEFINITE
>
> Generalized eigenproblem stored in file.
>
> Reading COMPLEX matrices from binary files...
> [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 12:40:20 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
> [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
> 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
>
>
>
>
> On 25/07/2013 12:34, Jose E. Roman wrote:
>
> El 25/07/2013, a las 13:30, Toby escribió:
>
> 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
>
> By default, ST is using PCREDUNDANT. In order to use the above
> options you must change it to e.g. PCLU.
> That is, add -st_pc_type lu
>
> Jose
>
>
> 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
>
>
>
>
> --
> 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/cf15c74d/attachment-0001.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/cf15c74d/attachment-0001.ksh>
More information about the petsc-users
mailing list