Allen–Cahn Phase-Field Void Shrinking – Why Does Mesh Adaptivity Cause “Steady-State” Behavior?

[PengWei97]

Dear MOOSE team,

I am running a simple test case of curvature-driven void shrinkage using a pure Allen–Cahn phase-field model. The physics is extremely standard and the problem runs perfectly on a uniform mesh.

However, once mesh adaptivity is activated, the solution immediately becomes “steady”—the void stops evolving, the nonlinear residual becomes extremely small (~1e-12), and the time step keeps shrinking until dt < dtmin, eventually causing the simulation to fail.

I would like to ask for guidance on:

1. Why adding mesh adaptivity causes the void to stop shrinking?
2. How to correctly configure the nonlinear/linear solvers to work with adaptivity?

---

1️⃣ Governing equation & material parameters

The simulation uses a standard Allen–Cahn phase-field equation:

[
\frac{\partial \xi}{\partial t}
= -L\left( g'(\xi) - \kappa \nabla^2 \xi \right),
\qquad
g(\xi)=w ,\xi^2(1-\xi)^2
]

where

• xi = 1 → solid, xi = 0 → void
• Mobility L = 1e-6
• Gradient energy coefficient kappa = 4.5e-7
• Double-well height w = 3.5e6

Initial condition is a diffuse circular void (radius=15 μm, int_width=1 μm).

---

2️⃣ Baseline case (uniform 250×250 mesh): works perfectly

Below is the link to the full input (no adaptivity):
👉 inp2_pf_with_void0.i

Mesh block & Executioner excerpt

[Mesh]
  [./mesh]
    type   = GeneratedMeshGenerator
    dim    = 2
    nx     = 250
    ny     = 250
    xmin   = 0.0
    xmax   = 100e-6     # [m]
    ymin   = 0.0
    ymax   = 100e-6     # [m]
    elem_type = QUAD4
  [../]
[]
  

[Executioner]
  type        = Transient
  scheme      = bdf2
  solve_type  = PJFNK

  # Linear solver parameters
  petsc_options_iname  = '-pc_type -pc_hypre_type -ksp_gmres_restart -pc_hypre_boomeramg_strong_threshold'
  petsc_options_value  = 'hypre boomeramg 31 0.7'

  l_tol      = 1e-4
  l_max_its  = 30
  nl_max_its = 25
  nl_rel_tol = 1e-7

  end_time = 250.0
  dtmax    = 50.0
  dtmin    = 1e-8

  [./TimeStepper]
    type               = IterationAdaptiveDT
    dt                 = 1e-1         # initial timestep
    growth_factor      = 1.3
    cutback_factor     = 0.7
    optimal_iterations = 8
  [../]
[]

Log excerpt (uniform mesh)

Allen–Cahn behaves as expected. Residual decreases smoothly; Newton fully converges:

Time Step 37, time=21.2391
 0 Nonlinear |R| = 1.48e-13
 ...
 5 Nonlinear |R| = 1.91e-23
Solve Converged!

Void evolution result (works correctly)

(👉 Here I will insert the figures showing the void shrinking over time)

3️⃣ After enabling Mesh Adaptivity → void stops shrinking + solver “pseudo-diverges”

Adaptivity is added here (link):
👉 inp2_pf_with_void0_adap.i

Adaptivity block excerpt

[Adaptivity]
  initial_steps = 4
  initial_marker = err_xi

  max_h_level = 4
  marker = err_bnds

  [./Indicators]
    [./ind_xi]
      type = GradientJumpIndicator
      variable = xi
    [../]
    [./ind_bnds]
      type = GradientJumpIndicator
      variable = bnds
    [../]
  [../]

  [./Markers]
    [./err_xi]
      type = ErrorFractionMarker
      refine = 0.95
      coarsen = 0.3
      indicator = ind_xi
    [../]
    [./err_bnds]
      type = ErrorFractionMarker
      refine = 0.95
      coarsen = 0.3
      indicator = ind_bnds
    [../]
  [../]
[]

Everything else (Kernels, Materials, ICs, dt control) is unchanged.

---

❌ Problem: residual stagnates at 1e-12, Newton stalls, dt shrinks → failure

Typical log after first adapt step:

15 Linear |R| = 1.589079e-12
...
30 Linear |R| = 1.589035e-12
Linear solve did not converge due to DIVERGED_ITS (30 iterations)

6 Nonlinear |R| = 1.589102e-12
0 Linear |R| = 1.589102e-12
...
30 Linear |R| = 1.589035e-12
Linear solve did not converge due to DIVERGED_ITS

After several repeats:

Solve failed, cutting timestep
...
dt < dtmin
Simulation aborted

At this point xi field no longer evolves, and the void stays frozen.

---

4️⃣ Attempted fix from ChatGPT suggestions: use only absolute tolerances

I tested the following solver settings:

l_tol      = 0.0
l_abs_tol  = 1e-10
l_max_its  = 50

nl_rel_tol = 0.0
nl_abs_tol = 1e-10
nl_max_its = 20

This removes the relative tolerance criterion.

However, the result becomes:

Time Step 8, time = 2.38577, dt = 0.627485
 0 Nonlinear |R| = 2.341290e-12
 Solve Converged!

The nonlinear residual is again ~1e-12 at the FIRST Newton evaluation, so Newton takes 0 iterations, accepts the step, and the void still does not evolve at all.

This is extremely puzzling.

---

❓ My questions to the MOOSE team

Q1. Why does adding mesh adaptivity cause the Allen–Cahn void to immediately fall into an apparent “steady-state” where the residual is ~1e-12 and Newton cannot update the solution?

• On the uniform mesh, the void shrinks exactly as expected.
• After the first adaptive mesh update + solution projection, the system suddenly behaves as if it has reached a discrete steady state.
• The chemical driving force seems to vanish numerically even though the physical model has not changed.

What is the typical cause?

• Projection of xi onto refined mesh creating a locally flat interface?
• Indicator/marker design?
• Over-refinement / over-coarsening near the interface?
• Impact of BDF2 history terms after adaptivity?
• Something about ACInterface + adaptivity?

---

Q2. Solver configuration: how should linear/nonlinear tolerances be chosen for Allen–Cahn + adaptivity?

I suspect the issue may be solver-related because:

• Linear residual remains almost constant (1e-12 → 1e-12), GMRES cannot reduce it further.
• Absolute tolerances accept the solution immediately.
• Relative tolerances never converge (DIVERGED_ITS).
• Newton does 0 updates and time-stepping thinks this is “converged”.

What is the recommended convergence strategy for Allen–Cahn when mesh adaptivity is active?
Should I:

• tighten/loosen l_abs_tol?
• force a minimum number of Newton iterations?
• switch to implicit Euler instead of BDF2?
• disable coarsening near interfaces?
• avoid projecting both history solutions in BDF2?

---

📌 I would greatly appreciate guidance.

This behavior (void shrinks on uniform mesh, but becomes steady immediately after first adapt) seems nonphysical. I would like to understand whether:

• this is a known issue with adaptivity + Allen–Cahn,
• an artifact of my adaptivity setup,
• a solver configuration issue,
• or a misunderstanding of how adaptivity transforms the AC dynamics.

Thanks very much for your help!

Wei

[GiudGiud]

Hello

If the solver thinks it has converged already it won't do anything and essentially skip the step. This could explain the problem.

I would recommend you force a solve using nl_forced_its in the executioner block.

Also the starting residual is very small. Most physics are happily converged with 1e-12

You could try to increase that by either using a different unit system for your equations OR using large scaling factors (applied in the Variables block for each variable, or using automatic scaling in the Executioner)

[PengWei97]

Hi @GiudGiud,

Thank you very much for your detailed explanation and suggestions — they are really helpful.

I tried

I would recommend you force a solve using nl_forced_its in the executioner block.

nl_forced_its = 5   # The number of forced nonlinear iterations

in the Executioner block, and with this setting the simulation can run to completion. However, this is not exactly what I would like, because my ultimate goal is a tightly coupled multiphysics simulation (phase-field + diffusion + mechanics + electric field). I would prefer to understand and fix the underlying convergence issues rather than relying on forced nonlinear iterations.

“most physics are happily converged with 1e-12”.

It is extremely valuable modeling experience for me. Do you happen to know of any references or learning materials that systematically explain how to choose tolerances, interpret residuals, and tune nonlinear/linear solvers and preconditioning in MOOSE (or PETSc)?

In this project I’ve noticed that the “basic” work (e.g. deriving equations, implementing kernels and material classes) is relatively quick, while most of the time is spent on making the coupled solve robust and efficient (avoiding solver failure, improving convergence and performance). This seems to be largely an Executioner (and PETSc) issue.

If you think learning PETSc more deeply would help, could you suggest which specific aspects are most relevant for MOOSE users? For example, should I focus on topics like:

• nonlinear solvers options and strategies
• linear solvers and preconditioners
• Scaling of residuals and variables
• Any other PETSc features that are particularly important for multiphysics problems in MOOSE?

using a different unit system for your equations
this turned out to be very effective in my case — I’ve already verified that changing the unit system significantly improves the behavior of the nonlinear solve.

I still have a question about how to choose units when multiple fields are coupled. For example, I want to couple a phase-field, a diffusion field, and mechanics. In the example hex_grain_growth_2D_eldrforce.i, the phase-field part uses nm and ns, while LinearElasticStress uses MPa.

Is there a recommended strategy for designing the unit system in such a multiphysics context? For instance:

• More generally, how do you usually select characteristic length and time scales for phase-field and diffusion variables when they are coupled to elasticity?

using a different unit system for your equations OR using large scaling factors (applied in the Variables block for each variable, or using automatic scaling in the Executioner)
As for scaling: both help and can make the solve run, but in my tests they still seem less effective than choosing a better unit system.

When I couple the phase-field and diffusion equations and also enable mesh adaptivity, the simulation currently fails to converge; after reading your reply, I suspect that optimizing the unit system might be a more fundamental solution than just tweaking scaling. Do you have any best practices or recommendations for:

Thank you again for your time and for sharing your experience — it’s been extremely helpful for me in understanding where to focus my efforts in making these multiphysics simulations robust.

Wei

[GiudGiud]

Do you happen to know of any references or learning materials that systematically explain how to choose tolerances, interpret residuals, and tune nonlinear/linear solvers and preconditioning in MOOSE (or PETSc)?

This page might be helpful.
https://mooseframework.inl.gov/moose/source/systems/NonlinearSystemBase.html#scaling

Scaling of residuals and variables

This one should be the easiest to learn and fastest to yield results imo

linear solvers and preconditioners

The PETSc website has good resources. You can get scaling improvements from picking a better choice

nonlinear solvers options and strategies

We support newton and PJFNK. Anything else for the nonlinear solver would be possible but a lot of work to implement

More generally, how do you usually select characteristic length and time scales for phase-field and diffusion variables when they are coupled to elasticity?
I suspect that optimizing the unit system might be a more fundamental solution than just tweaking scaling. Do you have any best practices or recommendations for:

Picking units is actually quite similar to scaling the equations. I suspect it has roughly the same effect.

[PengWei97]

Thank you for your detailed answers. I will study further bases on the materials you provided. And I will ask more questions on Github if I have any new questions.