[petsc-users] Debugging SNES when "everything" looks okay

[Ellen M. Price]

Hi PETSc users,

I'm having trouble applying SNES to a new problem I'm working on. I'll
try to be as complete as possible but can't post the full code because
it's ongoing research and is pretty long anyway.

The nonlinear problem arises from trying to solve a set of two coupled
ODEs using a Galerkin method. I am using Mathematica to generate the
objective function to solve and the Jacobian, so I *think* I can rule
out human error on that front.

There are four things I can easily change:

- number of DMDA grid points N (I've tried 100 and 1000)
- preconditioner (I've tried LU and SVD, LU appears to work better, and
SVD is too slow for N = 1000)
- linear solver (haven't played with this much, but sometimes smaller
tolerances work better)
- nonlinear solver (what I'm having trouble with)

Trying different solvers should, as far as I'm aware, give comparable
answers, but that's not the case here. NEWTONTR works best of the ones
I've tried, but I'm suspicious that the function value barely decreases
before SNES "converges" -- and none of the options I've tried changing
seem to affect this, as it just finds another reason to converge without
making any real progress. For example:

  0 SNES Function norm 7.197788479418e+00
    0 KSP Residual norm 1.721996766619e+01
    1 KSP Residual norm 5.186021549059e-14
  Linear solve converged due to CONVERGED_STEP_LENGTH iterations 1
  1 SNES Function norm 7.197777674987e+00
    0 KSP Residual norm 3.296112185897e+01
    1 KSP Residual norm 2.713415432045e-13
  Linear solve converged due to CONVERGED_STEP_LENGTH iterations 1
.....
 50 SNES Function norm 7.197376803542e+00
    0 KSP Residual norm 6.222518656302e+02
    1 KSP Residual norm 9.630659996504e-12
  Linear solve converged due to CONVERGED_STEP_LENGTH iterations 1
 51 SNES Function norm 7.197376803542e+00
Nonlinear solve converged due to CONVERGED_SNORM_RELATIVE iterations 50

I've tried everything I can think of and the FAQ suggestions, including
non-dimensionalizing the problem; I observe the same behavior either
way. The particularly strange thing I can't understand is why some of
the SNES methods fail outright, after just one iteration (NCG, NGMRES,
and others) with DIVERGED_DTOL. Unless I've misunderstood, it seems
like, for the most part, I should be able to substitute in one method
for another, possibly adjusting a few parameters along the way.
Furthermore, the default method, NEWTONLS, diverges with
DIVERGED_LINE_SEARCH, which I'm not sure how to debug.

So the only viable method is NEWTONTR, and that doesn't appear to
"really" converge. Any suggestions for further things to try are
appreciated. My current options are:

-pc_type lu -snes_monitor -snes_converged_reason -ksp_converged_reason
-snes_max_it 10000 -ksp_monitor -snes_type newtonls

Thanks in advance,
Ellen Price