[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