[petsc-users] Matrix missing diagonal entry

[Florian Lindner]

Hello,

when I try to KSPSolve my matrix petsc prints out that a diagonal entry is missing:

 KSPSetOperators(solver, _matCLU.matrix, _matCLU.matrix );
 KSPSolve(solver, vin.vector, vp.vector);

And petsc is perfectly right about that specific entry being zero. So this may be more of a mathematical problem then pure petsc.

The algorithm that I try to reimplement does the following:

  lu(_matrixCLU, _pivotsCLU);  // Compute LU decomposition
  // Left out: Permutate in according to pivotsCLU
  forwardSubstitution(_matrixCLU, in, y); // CLU^-1 * in = y  (lower triangle of CLU)
  backSubstitution(_matrixCLU, y, p);     // CLU^-1 * y = p (upper triangle of CLU)
  multiply(_matrixA, p, out );            // out = A * p (not relevant here)

I want to replace to lu, forward- and backwardSubstituation by using a petsc KSP solver. Before doing any of these operations _matrixCLU and _matCLU are equal:

=============== CLU Before LU =======
==== Dynamic Matrix of size 7 x 7 ====
 1.0000 0.3679 0.1353 0.3679 1.0000 0.0000 0.0000
 0.3679 1.0000 0.3679 0.1353 1.0000 1.0000 0.0000
 0.1353 0.3679 1.0000 0.3679 1.0000 1.0000 1.0000
 0.3679 0.1353 0.3679 1.0000 1.0000 0.0000 1.0000
 1.0000 1.0000 1.0000 1.0000 0.0000 0.0000 0.0000
 0.0000 1.0000 1.0000 0.0000 0.0000 0.0000 0.0000
 0.0000 0.0000 1.0000 1.0000 0.0000 0.0000 0.0000

Mat Object:CLU 1 MPI processes
  type: seqsbaij
row 0: (0, 1)  (1, 0.367879)  (2, 0.135335)  (3, 0.367879)  (4, 1)  (5, 0)  (6, 0) 
row 1: (1, 1)  (2, 0.367879)  (3, 0.135335)  (4, 1)  (5, 1)  (6, 0) 
row 2: (2, 1)  (3, 0.367879)  (4, 1)  (5, 1)  (6, 1) 
row 3: (3, 1)  (4, 1)  (5, 0)  (6, 1) 
row 4:
row 5:
row 6:

LU decomposition of matrixCLU works, no rank deficiency (no 0 on the main diagonal after LU and I tested for full rank of matrixCLU). Petsc complains with this matrix:

[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
[0]PETSC ERROR: Object is in wrong state
[0]PETSC ERROR: Matrix is missing diagonal entry 4
[0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.
[0]PETSC ERROR: Petsc Release Version 3.5.2, unknown 
[0]PETSC ERROR: ../build/debug/binprecice on a arch-linux2-c-debug named helium by lindnefn Mon Oct  6 16:26:49 2014
[0]PETSC ERROR: Configure options 
[0]PETSC ERROR: #1 MatICCFactorSymbolic_SeqSBAIJ() line 2430 in /data2/scratch/lindner/petsc/src/mat/impls/sbaij/seq/sbaijfact2.c
[0]PETSC ERROR: #2 MatICCFactorSymbolic() line 6266 in /data2/scratch/lindner/petsc/src/mat/interface/matrix.c
[0]PETSC ERROR: #3 PCSetup_ICC() line 20 in /data2/scratch/lindner/petsc/src/ksp/pc/impls/factor/icc/icc.c
[0]PETSC ERROR: #4 PCSetUp() line 902 in /data2/scratch/lindner/petsc/src/ksp/pc/interface/precon.c
[0]PETSC ERROR: #5 KSPSetUp() line 305 in /data2/scratch/lindner/petsc/src/ksp/ksp/interface/itfunc.c
[0]PETSC ERROR: #6 KSPSolve() line 417 in /data2/scratch/lindner/petsc/src/ksp/ksp/interface/itfunc.c
[0]PETSC ERROR: #7 map() line 696 in src/mapping/PetRadialBasisFctMapping.hpp
  numerical equality test failed: 0.00000000000000000000e+00 instead of 1.00000000000000000000e+00
  file: src/mapping/tests/PetRadialBasisFctMappingTest.cpp       line: 196
  statement: value==1.0
  difference:  1.00000000000000000000e+00


What's indeed perfectly right. But what is the best way to cope with that?

Thanks,
Florian