Multigrid Questions
Barry Smith
bsmith at mcs.anl.gov
Wed Aug 26 15:36:48 CDT 2009
Please run with
-ksp_monitor_true_residual -snes_monitor -
mg_levels_ksp_monitor_true_residua -mg_coarse_ksp_type gmres -
mg_coarse_ksp_monitor_true_residual -snes_view
and send the output for one time-step.
Barry
On Aug 26, 2009, at 2:34 PM, Hal Finkel wrote:
> Hello,
>
> I've constructed a evolution scheme for a time-dependent 3D PDE
> system using the DMMG layer. The core of the code is essentially a
> loop around DMMGSolve in the spirit of snes/examples/tutorials/
> ex29.c. If I only run with one level, then the code runs and
> produces sensible results. Passing the option -mg_levels_0_ksp_type
> bcgsl increases the solver speed by about 30%, and -
> mg_levels_0_pc_type fieldsplit gives me another 5% improvement. Does
> that say something about the structure of my problem?
>
> If I run with two levels, the coarse solve seems to work, but the
> fine grid's KSP never seems to complete. If I run with -
> mg_levels_1_ksp_monitor_true_residual, I can see the residual does
> not seem to be converging at all:
>
> 0 KSP preconditioned resid norm 1.979915874387e+02 true resid
> norm 1.889855551272e+04 ||Ae||/||Ax|| 1.889855551272e+04
> 0 KSP preconditioned resid norm 3.833279709488e+06 true resid
> norm 3.665290817713e+08 ||Ae||/||Ax|| 3.665290817713e+08
> 0 KSP preconditioned resid norm 5.080465654647e+01 true resid
> norm 8.284057252719e+03 ||Ae||/||Ax|| 8.284057252719e+03
> 1 KSP preconditioned resid norm 6.057649668915e+00 true resid
> norm 5.136009463814e+02 ||Ae||/||Ax|| 5.136009463814e+02
> 0 KSP preconditioned resid norm 5.511785811223e+04 true resid
> norm 1.055162492334e+07 ||Ae||/||Ax|| 1.055162492334e+07
> 0 KSP preconditioned resid norm 2.164672340919e+01 true resid
> norm 2.999406329054e+04 ||Ae||/||Ax|| 2.999406329054e+04
> 0 KSP preconditioned resid norm 5.275983407134e+05 true resid
> norm 9.196351032040e+08 ||Ae||/||Ax|| 9.196351032040e+08
> 0 KSP preconditioned resid norm 3.338344860439e+01 true resid
> norm 5.133553722923e+03 ||Ae||/||Ax|| 5.133553722923e+03
> 1 KSP preconditioned resid norm 4.306407436023e+00 true resid
> norm 3.361745947223e+02 ||Ae||/||Ax|| 3.361745947223e+02
> ...
>
> If I run with -pc_mg_monitor I see, after the course solve, only
> things like:
> Residual norms for mg_levels_1_ solve.
> 0 KSP Residual norm 1.979915874387e+02
> Residual norms for mg_levels_1_ solve.
> 0 KSP Residual norm 3.833279709488e+06
> Residual norms for mg_levels_1_ solve.
> 0 KSP Residual norm 5.080465654647e+01
> 1 KSP Residual norm 6.057649668915e+00
> Residual norms for mg_levels_1_ solve.
> 0 KSP Residual norm 5.511785811223e+04
> Residual norms for mg_levels_1_ solve.
> 0 KSP Residual norm 2.164672340919e+01
> Residual norms for mg_levels_1_ solve.
> 0 KSP Residual norm 5.275983407134e+05
> Residual norms for mg_levels_1_ solve.
> 0 KSP Residual norm 3.338344860439e+01
> 1 KSP Residual norm 4.306407436023e+00
> ...
>
> running with -info shows:
> [0] KSPDefaultConverged(): user has provided nonzero initial guess,
> computing 2-norm of preconditioned RHS
> [0] KSPDefaultConverged(): Linear solver is diverging. Initial right
> hand size norm 0.0066162, current residual norm 23282.7 at iteration 0
> [0] KSPDefaultConverged(): user has provided nonzero initial guess,
> computing 2-norm of preconditioned RHS
> [0] KSPDefaultConverged(): Linear solver is diverging. Initial right
> hand size norm 0.00287711, current residual norm 38.2205 at
> iteration 0
> [0] KSPDefaultConverged(): user has provided nonzero initial guess,
> computing 2-norm of preconditioned RHS
> [0] KSPDefaultConverged(): Linear solver is diverging. Initial right
> hand size norm 0.00287711, current residual norm 809191 at iteration 0
> [0] KSPDefaultConverged(): user has provided nonzero initial guess,
> computing 2-norm of preconditioned RHS
> [0] KSPDefaultConverged(): Linear solver is diverging. Initial right
> hand size norm 0.00970186, current residual norm 186.983 at
> iteration 0
> [0] KSPDefaultConverged(): user has provided nonzero initial guess,
> computing 2-norm of preconditioned RHS
> [0] KSPDefaultConverged(): Linear solver is diverging. Initial right
> hand size norm 0.00970186, current residual norm 3.62292e+06 at
> iteration 0
>
> What is a good way to debug this?
>
> No matter how I run, I've arranged the code so that the fine grid is
> *always* the same size; it has to be because the initial conditions
> are loaded from a file.
>
> In general, however, I don't understand exactly what the mg
> preconditioner is doing. Is there a writeup somewhere which is more
> verbose than the manual? In the one-level case, how are the results
> of the solve used to precondition the Jacobian?
>
> In an 1-processor, 1-level case, running with -snes_view gives:
> SNES Object:
> type: ls
> line search variant: SNESLineSearchCubic
> alpha=0.0001, maxstep=1e+08, minlambda=1e-12
> maximum iterations=50, maximum function evaluations=10000
> tolerances: relative=1e-08, absolute=1e-50, solution=1e-08
> total number of linear solver iterations=36
> total number of function evaluations=3
> KSP Object:
> type: fgmres
> 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
> right preconditioning
> PC Object:
> type: mg
> MG: type is FULL, levels=1 cycles=v
> Coarse gride solver -- level 0 presmooths=1 postsmooths=1 -----
> KSP Object:(mg_levels_0_)
> type: bcgsl
> BCGSL: Ell = 2
> BCGSL: Delta = 0
> maximum iterations=1, initial guess is zero
> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
> left preconditioning
> PC Object:(mg_levels_0_)
> type: fieldsplit
> FieldSplit with MULTIPLICATIVE composition: total splits =
> 2, blocksize = 2
> Solver info for each split is in the following KSP objects:
> Split number 0 Fields 0
> KSP Object:(mg_levels_0_fieldsplit_0_)
> type: preonly
> maximum iterations=10000, initial guess is zero
> tolerances: relative=1e-05, absolute=1e-50,
> divergence=10000
> left preconditioning
> PC Object:(mg_levels_0_fieldsplit_0_)
> type: ilu
> ILU: 0 levels of fill
> ILU: factor fill ratio allocated 1
> ILU: tolerance for zero pivot 1e-12
> ILU: using diagonal shift to prevent zero pivot
> ILU: using diagonal shift on blocks to prevent zero pivot
> out-of-place factorization
> matrix ordering: natural
> ILU: factor fill ratio needed 1
> Factored matrix follows
> Matrix Object:
> type=seqaij, rows=19683, cols=19683
> package used to perform factorization: petsc
> total: nonzeros=531441, allocated nonzeros=531441
> not using I-node routines
> linear system matrix = precond matrix:
> Matrix Object:
> type=seqaij, rows=19683, cols=19683
> total: nonzeros=531441, allocated nonzeros=531441
> not using I-node routines
> Split number 1 Fields 1
> KSP Object:(mg_levels_0_fieldsplit_1_)
> type: preonly
> maximum iterations=10000, initial guess is zero
> tolerances: relative=1e-05, absolute=1e-50,
> divergence=10000
> left preconditioning
> PC Object:(mg_levels_0_fieldsplit_1_)
> type: ilu
> ILU: 0 levels of fill
> ILU: factor fill ratio allocated 1
> ILU: tolerance for zero pivot 1e-12
> ILU: using diagonal shift to prevent zero pivot
> ILU: using diagonal shift on blocks to prevent zero pivot
> out-of-place factorization
> matrix ordering: natural
> ILU: factor fill ratio needed 1
> Factored matrix follows
> Matrix Object:
> type=seqaij, rows=19683, cols=19683
> package used to perform factorization: petsc
> total: nonzeros=531441, allocated nonzeros=531441
> not using I-node routines
> linear system matrix = precond matrix:
> Matrix Object:
> type=seqaij, rows=19683, cols=19683
> total: nonzeros=531441, allocated nonzeros=531441
> not using I-node routines
> linear system matrix = precond matrix:
> Matrix Object:
> type=seqaij, rows=39366, cols=39366
> total: nonzeros=2125764, allocated nonzeros=2125764
> using I-node routines: found 19683 nodes, limit used is 5
> linear system matrix = precond matrix:
> Matrix Object:
> type=seqaij, rows=39366, cols=39366
> total: nonzeros=2125764, allocated nonzeros=2125764
> using I-node routines: found 19683 nodes, limit used is 5
>
> Thank you in advance,
> Hal
More information about the petsc-users
mailing list