Multigrid Questions
Hal Finkel
hal.finkel at yale.edu
Wed Aug 26 14:34:20 CDT 2009
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