I'm solving a "linear elasticity with contact" problem using SNES. The matrix is symmetric, and it corresponds to the standard linear elasticity matrix, plus some extra nonlinear "penalty spring" terms due to the contact conditions.
The solve works fine if I use a direct solver (e.g. MUMPS). I've also found that CG usually works fine too, but in one specific case it fails to converge. (I thought that penalty springs + iterative solvers may not be a good combination, but it does seem to work fine in most cases.) Below I've pasted the (last part of the) output of "-ksp_type cg -ksp_norm_type natural -snes_view -ksp_monitor" for the case that fails. If anyone has some suggestions about other solver/preconditioner options that would be worth trying, that would be most appreciated. Thanks, David ------------------------------------------------------------------------------- 102 KSP Residual norm 8.910585826247e-05 103 KSP Residual norm 8.395178827803e-05 104 KSP Residual norm 7.758478472168e-05 105 KSP Residual norm 7.099534889984e-05 106 KSP Residual norm 6.615351528474e-05 107 KSP Residual norm 6.303646443688e-05 NL step 5, |residual|_2 = 1.996487e+01 0 KSP Residual norm 5.844456247091e+00 1 KSP Residual norm 1.308756943973e+00 2 KSP Residual norm 7.545133187463e-01 3 KSP Residual norm 1.346068942607e+00 SNES Object: 1 MPI processes type: newtonls maximum iterations=50, maximum function evaluations=10000 tolerances: relative=1e-08, absolute=1e-50, solution=1e-08 total number of linear solver iterations=535 total number of function evaluations=11 norm schedule ALWAYS SNESLineSearch Object: 1 MPI processes type: bt interpolation: cubic alpha=1.000000e-04 maxstep=1.000000e+08, minlambda=1.000000e-12 tolerances: relative=1.000000e-08, absolute=1.000000e-15, lambda=1.000000e-08 maximum iterations=40 KSP Object: 1 MPI processes type: cg maximum iterations=10000, initial guess is zero tolerances: relative=1e-05, absolute=1e-50, divergence=10000 left preconditioning using NATURAL 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 matrix ordering: natural factor fill ratio given 1, needed 1 Factored matrix follows: Mat Object: 1 MPI processes type: seqaij rows=967608, cols=967608 package used to perform factorization: petsc total: nonzeros=7.55539e+07, allocated nonzeros=7.55539e+07 total number of mallocs used during MatSetValues calls =0 using I-node routines: found 322536 nodes, limit used is 5 linear system matrix = precond matrix: Mat Object: () 1 MPI processes type: seqaij rows=967608, cols=967608 total: nonzeros=7.55539e+07, allocated nonzeros=7.55539e+07 total number of mallocs used during MatSetValues calls =0 using I-node routines: found 322536 nodes, limit used is 5 Number of nonlinear iterations: 5 Nonlinear solver convergence/divergence reason: DIVERGED_LINEAR_SOLVE