Poor nonlinear solver performance for VEP rheology + particles + continous compositional fields

[naliboff]

While updating the continental extension bookmark (see attached PRM), I found that switching between using continuous and discontinuous compositional fields (set Use discontinuous composition discretization = true or false) in combination with particles for track advected properties (viscoelastic stress, strain, and chemical compositions) has a large effect on nonlinear solver convergence behavior.

For set Use discontinuous composition discretization = true, the Newton solver will converge in anywhere from ~ 8-16 nonlinear iterations.

For set Use discontinuous composition discretization = false, the Newton Solver on many time steps will often reach the maximum number of allowed nonlinear iterations (200) and will clearly not converge to the selected tolerance with further iterations.

After discussion with @gassmoeller, this makes intuitive sense as interpolating the particle values to a discontinuous compositional fields makes sense for certain (cell average used here) or all interpolation schemes (cell average was used here).

I have not checked whether using a DG field for a model tracking properties on fields has the same effect, but will follow up with results checking that. I will also test whether this choice influences other viscoelastic benchmarks without plasticity.

Overall, somewhat of a mute point as #4370 enforces using DG compositional fields with viscoelasticity, but certainly useful information in the short term and potentially for future issues.

continental_extension_vep_particles.prm.txt

[naliboff]

After a bit more testing (see attached plot), it looks like the best convergence behavior for the modified Continental Extension cookbook (elasticity added, etc) with particles is obtained by using the cell average (CAV) particle interpolation scheme (+ DG composition field ~ DGC = True).

Using the Quadratic or Bilinear Least Squares + limiters (QLS + LIM, BLS + LIM) produces poorer convergence behavior.

Using a continuous composition field (default option) with a cell average particle interpolation scheme fails to produce convergence within 200 nonlinear iterations (max value) on multiple time steps.

I can do further testing to see if #5685 (@bangerth) or #5754 (@YiminJin) help.

At this stage, I think it is sufficient to recommend that users just use a discontinuous composition field with VEP models, as this will be required when #4370 is merged. I will make a note of this in the update continental extension cookbook documentation (working on PR now).

[gassmoeller]

Just to summarize the discussion at the hackathon: @anne-glerum and @MFraters could not reproduce the problems after switching off the operator splitting for compositional fields, but the behavior remains for @naliboff. @naliboff: Could you try the new particle interpolator in #5815 and check if the convergence is still problematic?

[anne-glerum]

An update based on further testing:

• Even with operator splitting off for the particles, @naliboff finds that on the same machine, with the same number of processes, the timestep of bad convergence differs depending on the version of deal.II and underlying libraries.
• A run by @MFraters with deal.II 9.6.0-pre and CAV on 1 process failed to converge in timestep 48. With BLSQ, this run converges fine for 53 timesteps. With twice the number of particles, 36 timesteps converge (then run was stopped).
• A run on 1 process with deal.II 9.5.0 and BLSQ by @anne-glerum fails at timestep 17.

[naliboff]

A more detailed summary of the tests mentioned by @anne-glerum, which all use the same PRM file as @MFraters and are broken down into the following cases:

Case 1

• Run on branch that combined fix_stresses_elasticity (#4370) + commits (cherry picked) from #5689

• ASPECT build using deal.II 9.5.0

• Run on 32 processors

• Result: Bad convergence first appears on time step 9, and occurs again periodically (time steps 21, 37, ...)

Case 2

• Same as case 1, except the branch was built using deal.II main (a9dad19480)

• Result: Bad convergence first appears on time step 9, but this time subsequently reoccurs on different time steps (11, 38)

Case 3

• Run on a branch from @MFraters, which added a few commits to the branch used in cases 1 and 2. The changes in this branch were mainly to provide additional output about the Newton solver, but may also include changes that affect the system behavior.

• Run on 1 processor

• Built using deal.II 9.5.0

• Result: Bad convergence first appears on time step 13, and then time steps 86, 139, 143, 158, 163, and on

Case 4

• Same as case 3, except the branch was built using deal.II main (a9dad19480)

• Result: Bad convergence first appears on time step 30, and then time steps 108 and 129 (stopped there)

The main takeaways here are:

1. Small variations between my system (ubuntu 22.04, gnu compilers, etc) and deal.II version (+ dependencies) and those used by @MFraters lead to differences in the model results.
2. Both switching between (1) deal.II 9.5.0 and main and (2) my branch and Menno's branch (added commits) produces changes in the results, at least in terms of which time steps exhibit poor convergence behavior.