[petsc-users] troubleshooting AMG, coupled navier stokes, large eigenvalues on coarsest level
Mark Lohry
mlohry at gmail.com
Tue Nov 7 14:00:39 CST 2017
I've now got gamg running on matrix-free newton-krylov with the jacobian
provided by coloring finite differences (thanks again for the help). 3D
Poisson with 4th order DG or higher (35^2 blocks), gamg with default
settings is giving textbook convergence, which is great. Of course coupled
compressible navier-stokes is harder, and convergence is bad-to-nonexistent.
1) Doc says "Equations must be ordered in “vertex-major” ordering"; in my
discretization, each "node" has 5 coupled degrees of freedom (density, 3 x
momentum, energy). I'm ordering my unknowns:
rho_i, rhou_i, rhov_i, rhow_i, Et_i, rho_i+1, rhou_i+1, ... e.g. row-major
matrix order if you wrote the unknowns [{rho}, {rhou}, ... ].
and then setting block size to 5. Is that correct? I've also tried using
the actual sparsity of the matrix which has larger dense blocks (e.g.
[35x5]^2), but neither seemed to help.
2) With default settings, and with -pc_gamg_square_graph,
pc_gamg_sym_graph, agg_nsmooths 0 mentioned in the manual, the eigenvalue
estimates explode on the coarsest level, which I don't see with poisson:
Down solver (pre-smoother) on level 1 -------------------------------
KSP Object: (mg_levels_1_) 32 MPI processes
type: chebyshev
eigenvalue estimates used: min = 0.18994, max = 2.08935
eigenvalues estimate via gmres min 0.00933256, max 1.8994
Down solver (pre-smoother) on level 2 -------------------------------
KSP Object: (mg_levels_2_) 32 MPI processes
type: chebyshev
eigenvalue estimates used: min = 0.165969, max = 1.82566
eigenvalues estimate via gmres min 0.0290728, max 1.65969
Down solver (pre-smoother) on level 3 -------------------------------
KSP Object: (mg_levels_3_) 32 MPI processes
type: chebyshev
eigenvalue estimates used: min = 0.146479, max = 1.61126
eigenvalues estimate via gmres min 0.204673, max 1.46479
Down solver (pre-smoother) on level 4 -------------------------------
KSP Object: (mg_levels_4_) 32 MPI processes
type: chebyshev
eigenvalue estimates used: min = 6.81977e+09, max = 7.50175e+10
eigenvalues estimate via gmres min -2.76436e+12, max 6.81977e+10
What's happening here? (Full -ksp_view below)
3) I'm not very familiar with chebyshev smoothers, but they're only for SPD
systems (?). Is this an inappropriate smoother for this problem?
4) With gmres, the preconditioned residual is ~10 orders larger than the
true residual; and the preconditioned residual drops while the true
residual rises. I'm assuming this means something very wrong?
5) -pc_type hyper -pc_hypre_type boomeramg also works perfectly for the
poisson case, but hits NaN on the first cycle for NS.
KSP Object: 32 MPI processes
type: fgmres
restart=30, using Classical (unmodified) Gram-Schmidt Orthogonalization
with no iterative refinement
happy breakdown tolerance 1e-30
maximum iterations=50, initial guess is zero
tolerances: relative=1e-10, absolute=1e-07, divergence=10.
right preconditioning
using UNPRECONDITIONED norm type for convergence test
PC Object: 32 MPI processes
type: gamg
type is MULTIPLICATIVE, levels=5 cycles=v
Cycles per PCApply=1
Using externally compute Galerkin coarse grid matrices
GAMG specific options
Threshold for dropping small values in graph on each level = 0.1
0.1 0.1
Threshold scaling factor for each level not specified = 1.
AGG specific options
Symmetric graph true
Number of levels to square graph 10
Number smoothing steps 0
Coarse grid solver -- level -------------------------------
KSP Object: (mg_coarse_) 32 MPI processes
type: preonly
maximum iterations=10000, initial guess is zero
tolerances: relative=1e-05, absolute=1e-50, divergence=10000.
left preconditioning
using NONE norm type for convergence test
PC Object: (mg_coarse_) 32 MPI processes
type: bjacobi
number of blocks = 32
Local solve is same for all blocks, in the following KSP and PC
objects:
KSP Object: (mg_coarse_sub_) 1 MPI processes
type: preonly
maximum iterations=1, initial guess is zero
tolerances: relative=1e-05, absolute=1e-50, divergence=10000.
left preconditioning
using NONE norm type for convergence test
PC Object: (mg_coarse_sub_) 1 MPI processes
type: lu
out-of-place factorization
tolerance for zero pivot 2.22045e-14
using diagonal shift on blocks to prevent zero pivot [INBLOCKS]
matrix ordering: nd
factor fill ratio given 5., needed 0.
Factored matrix follows:
Mat Object: 1 MPI processes
type: seqaij
rows=0, cols=0, bs=5
package used to perform factorization: petsc
total: nonzeros=1, allocated nonzeros=1
total number of mallocs used during MatSetValues calls =0
not using I-node routines
linear system matrix = precond matrix:
Mat Object: 1 MPI processes
type: seqaij
rows=0, cols=0, bs=5
total: nonzeros=0, allocated nonzeros=0
total number of mallocs used during MatSetValues calls =0
not using I-node routines
linear system matrix = precond matrix:
Mat Object: 32 MPI processes
type: mpiaij
rows=5, cols=5, bs=5
total: nonzeros=25, allocated nonzeros=25
total number of mallocs used during MatSetValues calls =0
not using I-node (on process 0) routines
Down solver (pre-smoother) on level 1 -------------------------------
KSP Object: (mg_levels_1_) 32 MPI processes
type: chebyshev
eigenvalue estimates used: min = 0.18994, max = 2.08935
eigenvalues estimate via gmres min 0.00933256, max 1.8994
eigenvalues estimated using gmres with translations [0. 0.1; 0.
1.1]
KSP Object: (mg_levels_1_esteig_) 32 MPI processes
type: gmres
restart=30, using Classical (unmodified) Gram-Schmidt
Orthogonalization with no iterative refinement
happy breakdown tolerance 1e-30
maximum iterations=10, initial guess is zero
tolerances: relative=1e-12, absolute=1e-50, divergence=10000.
left preconditioning
using PRECONDITIONED norm type for convergence test
estimating eigenvalues using noisy right hand side
maximum iterations=2, nonzero initial guess
tolerances: relative=1e-05, absolute=1e-50, divergence=10000.
left preconditioning
using NONE norm type for convergence test
PC Object: (mg_levels_1_) 32 MPI processes
type: sor
type = local_symmetric, iterations = 1, local iterations = 1, omega
= 1.
linear system matrix = precond matrix:
Mat Object: 32 MPI processes
type: mpiaij
rows=70, cols=70, bs=5
total: nonzeros=1550, allocated nonzeros=1550
total number of mallocs used during MatSetValues calls =0
using nonscalable MatPtAP() implementation
not using I-node (on process 0) routines
Up solver (post-smoother) same as down solver (pre-smoother)
Down solver (pre-smoother) on level 2 -------------------------------
KSP Object: (mg_levels_2_) 32 MPI processes
type: chebyshev
eigenvalue estimates used: min = 0.165969, max = 1.82566
eigenvalues estimate via gmres min 0.0290728, max 1.65969
eigenvalues estimated using gmres with translations [0. 0.1; 0.
1.1]
KSP Object: (mg_levels_2_esteig_) 32 MPI processes
type: gmres
restart=30, using Classical (unmodified) Gram-Schmidt
Orthogonalization with no iterative refinement
happy breakdown tolerance 1e-30
maximum iterations=10, initial guess is zero
tolerances: relative=1e-12, absolute=1e-50, divergence=10000.
left preconditioning
using PRECONDITIONED norm type for convergence test
estimating eigenvalues using noisy right hand side
maximum iterations=2, nonzero initial guess
tolerances: relative=1e-05, absolute=1e-50, divergence=10000.
left preconditioning
using NONE norm type for convergence test
PC Object: (mg_levels_2_) 32 MPI processes
type: sor
type = local_symmetric, iterations = 1, local iterations = 1, omega
= 1.
linear system matrix = precond matrix:
Mat Object: 32 MPI processes
type: mpiaij
rows=270, cols=270, bs=5
total: nonzeros=7550, allocated nonzeros=7550
total number of mallocs used during MatSetValues calls =0
using nonscalable MatPtAP() implementation
not using I-node (on process 0) routines
Up solver (post-smoother) same as down solver (pre-smoother)
Down solver (pre-smoother) on level 3 -------------------------------
KSP Object: (mg_levels_3_) 32 MPI processes
type: chebyshev
eigenvalue estimates used: min = 0.146479, max = 1.61126
eigenvalues estimate via gmres min 0.204673, max 1.46479
eigenvalues estimated using gmres with translations [0. 0.1; 0.
1.1]
KSP Object: (mg_levels_3_esteig_) 32 MPI processes
type: gmres
restart=30, using Classical (unmodified) Gram-Schmidt
Orthogonalization with no iterative refinement
happy breakdown tolerance 1e-30
maximum iterations=10, initial guess is zero
tolerances: relative=1e-12, absolute=1e-50, divergence=10000.
left preconditioning
using PRECONDITIONED norm type for convergence test
estimating eigenvalues using noisy right hand side
maximum iterations=2, nonzero initial guess
tolerances: relative=1e-05, absolute=1e-50, divergence=10000.
left preconditioning
using NONE norm type for convergence test
PC Object: (mg_levels_3_) 32 MPI processes
type: sor
type = local_symmetric, iterations = 1, local iterations = 1, omega
= 1.
linear system matrix = precond matrix:
Mat Object: 32 MPI processes
type: mpiaij
rows=1610, cols=1610, bs=5
total: nonzeros=55550, allocated nonzeros=55550
total number of mallocs used during MatSetValues calls =0
using nonscalable MatPtAP() implementation
using I-node (on process 0) routines: found 6 nodes, limit used
is 5
Up solver (post-smoother) same as down solver (pre-smoother)
Down solver (pre-smoother) on level 4 -------------------------------
KSP Object: (mg_levels_4_) 32 MPI processes
type: chebyshev
eigenvalue estimates used: min = 6.81977e+09, max = 7.50175e+10
eigenvalues estimate via gmres min -2.76436e+12, max 6.81977e+10
eigenvalues estimated using gmres with translations [0. 0.1; 0.
1.1]
KSP Object: (mg_levels_4_esteig_) 32 MPI processes
type: gmres
restart=30, using Classical (unmodified) Gram-Schmidt
Orthogonalization with no iterative refinement
happy breakdown tolerance 1e-30
maximum iterations=10, initial guess is zero
tolerances: relative=1e-12, absolute=1e-50, divergence=10000.
left preconditioning
using PRECONDITIONED norm type for convergence test
estimating eigenvalues using noisy right hand side
maximum iterations=2, nonzero initial guess
tolerances: relative=1e-05, absolute=1e-50, divergence=10000.
left preconditioning
using NONE norm type for convergence test
PC Object: (mg_levels_4_) 32 MPI processes
type: sor
type = local_symmetric, iterations = 1, local iterations = 1, omega
= 1.
linear system matrix followed by preconditioner matrix:
Mat Object: 32 MPI processes
type: mffd
rows=153600, cols=153600
Matrix-free approximation:
err=1.49012e-08 (relative error in function evaluation)
Using wp compute h routine
Does not compute normU
Mat Object: 32 MPI processes
type: mpiaij
rows=153600, cols=153600, bs=5
total: nonzeros=65280000, allocated nonzeros=65280000
total number of mallocs used during MatSetValues calls =0
using I-node (on process 0) routines: found 960 nodes, limit used
is 5
Up solver (post-smoother) same as down solver (pre-smoother)
linear system matrix followed by preconditioner matrix:
Mat Object: 32 MPI processes
type: mffd
rows=153600, cols=153600
Matrix-free approximation:
err=1.49012e-08 (relative error in function evaluation)
Using wp compute h routine
Does not compute normU
Mat Object: 32 MPI processes
type: mpiaij
rows=153600, cols=153600, bs=5
total: nonzeros=65280000, allocated nonzeros=65280000
total number of mallocs used during MatSetValues calls =0
using I-node (on process 0) routines: found 960 nodes, limit used is 5
1 SNES Function norm 5.917486103148e+05
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20171107/7a38b769/attachment-0001.html>
More information about the petsc-users
mailing list