[petsc-users] KSPSolve convergence and norms
Matthew Knepley
knepley at gmail.com
Fri Oct 25 11:59:45 CDT 2013
On Fri, Oct 25, 2013 at 11:33 AM, Torquil Macdonald Sørensen <
torquil at gmail.com> wrote:
> Hi!
>
> I'm using PETSc 3.4.3, with only default runtime settings, to solve Ax=b,
> for complex numbers (I have used --with-scalar-type=complex).
>
> After reading sections 4.3.1, 4.3.2, 4.3.3 in the manual, I decided to
> investigate the convergence in my program, since I'm having some problems
> when the system is larger. So, I'm running it with the options
>
> -ksp_monitor_true_residual -ksp_converged_reason -ksp_view
>
> The output shows that the relative tolerance required for convergence is
> 1e-5, but none of the reported norms are that small. Despite this,
> convergence is reported (CONVERGED_RTOL).
>
> I though that CONVERGED_RTOL should indicate that ||r||_2/||b||_2 < 1e-5,
> where r is the preconditioned residual?
>
> Here is the output of my program (which includes some code to compute the
> l_2 norm of b, just before KSPSolve is run):
>
> *************************************************************
> Norm of b: 0.179829
> 0 KSP preconditioned resid norm 2.508789303280e+04 true resid norm
> 1.798290175843e-01 ||r(i)||/||b|| 1.000000000000e+00
> 1 KSP preconditioned resid norm 1.856991132478e+00 true resid norm
> 3.186462510386e-01 ||r(i)||/||b|| 1.771940120227e+00
> 2 KSP preconditioned resid norm 1.704678606994e-01 true resid norm
> 4.004071321891e-02 ||r(i)||/||b|| 2.226599119363e-01
>
So the relative preconditioned residual is:
1.7e-01 / 2.5e+04 < 1.0e-5
while the relative true residual is 2.2e-1. This happens because you are
using ILU(0). It can be really
crappy, and in this case the preconditioner is very ill-conditioned. I
can't imagine a problem where I
would recommend using ILU, but its our default because its the only
black-box PC that exists.
Thanks
Matt
> Linear solve converged due to CONVERGED_RTOL iterations 2
> KSP Object: 1 MPI processes
> type: gmres
> GMRES: restart=30, using Classical (unmodified) Gram-Schmidt
> Orthogonalization with no iterative refinement
> GMRES: happy breakdown tolerance 1e-30
> maximum iterations=10000, initial guess is zero
> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
> left preconditioning
> using PRECONDITIONED norm type for convergence test
> PC Object: 1 MPI processes
> type: ilu
> ILU: out-of-place factorization
> 0 levels of fill
> tolerance for zero pivot 2.22045e-14
> using diagonal shift on blocks to prevent zero pivot [INBLOCKS]
> matrix ordering: natural
> factor fill ratio given 1, needed 1
> Factored matrix follows:
> Matrix Object: 1 MPI processes
> type: seqaij
> rows=18, cols=18
> package used to perform factorization: petsc
> total: nonzeros=144, allocated nonzeros=144
> total number of mallocs used during MatSetValues calls =0
> using I-node routines: found 9 nodes, limit used is 5
> linear system matrix = precond matrix:
> Matrix Object: 1 MPI processes
> type: seqaij
> rows=18, cols=18
> total: nonzeros=144, allocated nonzeros=360
> total number of mallocs used during MatSetValues calls =18
> using I-node routines: found 9 nodes, limit used is 5
>
> ************************************************************************************************
>
> Here is the code I used to get the l_2 norm of b:
>
> PetscReal val;
> VecNorm(b, NORM_2, &val);
> PetscPrintf(PETSC_COMM_WORLD, "Norm of b: %G\n", val);
>
> Best regards
> Torquil Sørensen
>
--
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/20131025/edcfb9c8/attachment.html>
More information about the petsc-users
mailing list