TL;DR — Key Takeaways
- Divergence and instability fall into four categories: temporal (CFL/mesh), spatial (mesh quality), boundary/initial conditions, and solver settings.
- Residuals alone don’t prove convergence. Monitor points (forces, temperatures, velocities) must also be checked. Residuals can flatline while a localized variable blows up.
- Diagnose one variable at a time. Change only one setting, then rerun. If you change two things and the simulation stabilizes, you don’t know which fix actually worked.
- The boundary-swap trick (temporarily set all boundaries to fixed or default and rerun) is the fastest way to isolate whether divergence originates from boundary conditions or from time-step stability.
- Negative eigenvalues in FEA often signal buckling, not a solver error. The correct fix is an arc-length (Riks) method or implicit dynamics, not a tighter tolerance.
Why Simulations Diverge — And Why Diagnosis Matters
Every numerical simulation has the same terrifying signature: at some timestep or iteration, values explode toward NaN or infinity, the solver stops making progress, and you’re left staring at a terminal window full of Inf messages.
If you immediately restart with a smaller time step and the simulation runs fine, you got lucky. If you change something and the simulation still crashes, you’re now debugging with less information than before.
The problem isn’t that simulations are unreliable. They follow well-understood mathematical rules. The problem is that when a simulation diverges, the error message almost never tells you what went wrong — it just tells you that something went wrong.
Understanding the taxonomy of numerical failure is the first step toward diagnosing it systematically. Numerical failures break into four root-cause categories:
- Temporal instability — the time step is too large, the CFL condition is violated, or the solver steps outside its stability region.
- Spatial (mesh) pathology — highly skewed cells, extreme aspect ratios, or negative volumes create unphysical gradients.
- Boundary or initial condition errors — conflicting pressures, missing reference values, or inconsistent starting values create numerical contradictions.
- Solver configuration issues — wrong discretization schemes, aggressive under-relaxation factors, or inappropriate tolerances.
When you know which category your failure belongs to, diagnosis becomes a process of elimination. When you don’t, you end up guessing — and guessing is how you waste hours of compute time.
The Four Root-Cause Categories
1. Temporal Instability
Temporal instability is the most common cause of divergence in time-dependent simulations. The primary culprit is the CFL (Courant-Friedrichs-Lewy) condition, which states that for an explicit scheme, the time step must satisfy:
$$
C = \frac{u \Delta t}{\Delta x} \leq 1
$$
When $C > 1$, information propagates faster than the numerical scheme can track it, and errors compound exponentially. This is not a minor numerical nuisance — it’s a fundamental mathematical constraint.
What to look for:
- Simulation diverges immediately after the first few time steps
- The error message references “time step too small” or “stability limit exceeded”
- The same mesh and boundary conditions converge when you reduce the time step
What to fix:
- Reduce the time step until the CFL number drops below 0.5 as a safe starting point
- Use an adaptive time-stepping strategy that scales based on local Courant numbers
- Verify material property consistency — a high-velocity fluid in a dense medium may have a much tighter CFL limit than you expect
2. Spatial (Mesh) Pathology
Mesh quality is the second most common cause of divergence. Unlike temporal instability, which typically manifests early, mesh pathology can cause divergence at any stage — especially in nonlinear or multiphysics simulations where the mesh interacts with the solver iteratively.
Problem indicators:
- High cell skewness (> 0.9 on a 0–1 scale)
- Extreme aspect ratios (> 100:1)
- Negative cell volumes (the most severe error)
Mesh-generated divergence is particularly insidious because the mesh is usually correct for your physics — it’s only the numerical discretization that’s unhappy. A mesh that looks fine in a pre-processor may have cells with extreme aspect ratios in regions of rapid gradient change, causing unphysical diffusion terms.
What to fix:
- Inspect mesh quality metrics before the first run, not after
- Use a mesh quality dashboard (most CFD pre-processors include one) and verify skewness < 0.7, aspect ratio < 50 for most solvers
- Refine or restructure the mesh in high-gradient regions
- Never ignore negative cell volume warnings — they are mathematical impossibilities that the solver cannot resolve
3. Boundary and Initial Condition Errors
Boundary condition errors are the category where experienced researchers make the same mistake as beginners: the setup looks correct, but the physics contradicts itself.
Consider a common CFD scenario. You set a pressure inlet and a pressure outlet. Both have fixed pressure values. The solver needs a reference pressure to establish the absolute pressure field. Without a pressure reference, the solver tries to balance two incompatible constraints — it diverges.
Similarly, in structural simulations, a boundary condition that’s physically impossible — like constraining a node that also needs to displace freely — creates a contradiction the solver cannot resolve.
Common boundary condition mistakes:
- Conflicting pressure values at inlet and outlet (missing reference pressure)
- Backflow at an outlet (solver expects flow direction to match BC type)
- Missing reference values (temperature, pressure, electric potential)
- Inconsistent initial conditions (initial velocity does not match boundary velocity)
What to fix:
- Verify every boundary condition against the physics — does the boundary type match the expected flow direction?
- Add a pressure reference node when working with pressure boundaries
- Ensure initial conditions are consistent with boundary conditions and the expected physical state
4. Solver Configuration Issues
The solver settings category is where you can inadvertently destabilize a perfectly good simulation by using the wrong tool for the job. This category includes:
- Wrong discretization schemes — second-order upwind on a coarse mesh produces oscillations; first-order schemes are more stable but less accurate
- Aggressive under-relaxation — factors below 0.1 can cause divergence; typical values are 0.1–0.7 depending on the equation
- Inappropriate tolerances — too tight tolerances cause excessive iterations; too loose tolerances accept unphysical results
- Wrong solver type — using a pressure-based solver for incompressible flow vs a density-based solver for compressible flow
What to fix:
- Start with first-order discretization, then ramp up to second-order once the solution stabilizes
- Use default under-relaxation factors initially, then tune if needed
- Check solver type against your physics regime (compressible vs incompressible, turbulent vs laminar)
Systematic Diagnostic Workflow
When your simulation diverges, do not change multiple settings at once. This is the single most important debugging principle across CFD, FEA, and every numerical solver.
Follow this workflow:
Step 1: Make the simulation fail faster
If your simulation runs for 1000 iterations before diverging, you’re wasting time. Instead, start by using a larger time step or fewer iterations per step to force failure quickly. This gives you a tighter loop for diagnosis.
Step 2: Isolate the failure
Identify which variable diverges first — velocity, pressure, temperature, force? Use residual monitoring and monitor points (specific quantities you track at every iteration). The variable that diverges first tells you which equation and which part of the domain is problematic.
Critical insight: Residuals alone can hide divergence. As LEAP Australia’s CFD troubleshooting guide emphasizes, a solver may report converged residuals while monitor points (forces, velocities, temperatures) continue to drift. Always check monitor points alongside residuals.
Step 3: Change one thing, verify, then change another
After identifying the failure category, modify only one setting. Rerun the simulation. If it stabilizes, you’ve found a fix — but you don’t know if it was the right fix. Document what you changed and the result.
Then change only a second setting. Repeat. This controlled approach lets you isolate which change actually matters.
Step 4: Verify with the right diagnostics
Different simulation types require different diagnostic tools:
For CFD:
- Monitor residuals, monitor points, and velocity/pressure profiles
- Use residual contour visualization to see where residuals cluster spatially (enable expert mode:
solve set expert yes yes yes, then run one iteration:solve iterate 1) - Check continuity and momentum balance at domain boundaries
For FEA:
- Inspect Newton-Raphson convergence plots
- Look for negative eigenvalue warnings (which often indicate buckling rather than solver error)
- Check force-deflection curves for non-monotonic behavior
- Use arc-length methods (Riks) for snap-through or buckling problems
For general numerical code:
- Run a unit test with a manufactured solution (Method of Manufactured Solutions, MMS)
- Check for sign errors on diffusion or convection terms — the simplest cause of catastrophic divergence is often a misplaced minus sign
- Trace NaN or Inf back to the exact timestep where they first appear
Case Studies
Case 1: The Single Sign Error
A researcher debugging an advection-diffusion simulation found that the solution diverged completely — velocity values oscillated and grew exponentially. The mesh was fine, the CFL number was conservative, and boundary conditions were standard.
The root cause was a single sign error in the diffusion term. The discretization had +D instead of -D, effectively reversing the direction of diffusion. Instead of damping perturbations, the equation was amplifying them.
Lesson: Verify discretization signs first. Before checking mesh quality or solver settings, confirm that the mathematical form of your equations is correct. A sign error in a diffusion term is numerically equivalent to negative diffusion — the simulation diverges by design.
Case 2: Boundary Condition Swap (The FDTD Diagnostic Trick)
In an FDTD simulation, a user reported that changing boundary conditions fixed the divergence. The diagnostic workflow from Ansys FDTD documentation illustrates why:
Step 1: Temporarily set all boundary conditions to “Metal” (fixed) and rerun the simulation. If the simulation becomes stable, the divergence is PML-related.
Step 2: If the simulation still diverges with Metal boundaries, the divergence is due to $dt$ stability — the time step is too large for the solver to handle.
This boundary-swap trick is a practical diagnostic workflow that immediately narrows the search space. Instead of checking mesh, solver settings, material properties, and boundary conditions simultaneously, you isolate the boundary conditions as a category and fix the actual cause.
Case 3: Residual Contour Diagnostics
LEAP Australia’s expert TUI workflow demonstrates a powerful visualization technique. Instead of looking at scalar residual values, you generate residual contour plots and iso-surfaces that show exactly where in the domain high residuals cluster.
In one case, residuals flatlined at $10^{-3}$ while the simulation was actually diverging in a localized region. The residual contours showed extreme values near an inlet boundary where the mesh quality was poor (skewness > 0.85). The scalar residual masked the localized problem because it was averaged over the entire domain.
Lesson: Always visualize residuals spatially, not just numerically. Contour plots reveal where divergence originates.
Case 4: Non-Monotonic Stiffness in FEA
Enterfea’s nonlinear FEA case study covers a common failure mode. A structural simulation diverged at the equilibrium point where stiffness dropped sharply. A force-controlled solver cannot track the equilibrium point — it oscillates and fails.
The fix was switching from force control to displacement control or using an arc-length method (Riks). This is not a solver setting issue — it’s a fundamental mismatch between the solver strategy and the physics.
Lesson: For nonlinear problems with snap-through, buckling, or contact, use arc-length or displacement-control methods. Force-controlled solvers fail at equilibrium points where stiffness is non-monotonic.
Common Pitfalls — What Leads to Divergence
These mistakes are well-documented across CFD, FEA, and general numerical solvers. They all share one characteristic: the simulation looks correct until it diverges.
1. Changing Multiple Settings Simultaneously
When you change the time step, discretization scheme, and under-relaxation factor all at once, and the simulation stabilizes, you have no idea which change actually worked. Always change one setting at a time. This is the golden rule of numerical debugging.
2. Trusting Residuals Blindly
Residuals can flatline while a localized variable diverges. As the LEAP Australia checklist warns, you must also inspect monitor points — specific quantities like forces, velocities, and temperatures. Residuals measure the equation error, not the physical correctness of your solution.
3. Using Second-Order Schemes on Coarse Meshes
Second-order discretization produces oscillations on coarse meshes. Start with first-order discretization, stabilize the solution, then ramp up to second-order. This is not theoretical — it’s a practical workflow used by experienced CFD practitioners.
4. Ignoring Negative Cell Volumes
Negative cell volumes are mathematical impossibilities. They occur when a mesh cell is inverted or collapsed. If your pre-processor or mesh generator reports negative volumes, fix the mesh before attempting any simulation. No solver can resolve negative volumes.
5. Assuming Negative Eigenvalues Are Solvers Errors
In Abaqus and other FEA solvers, persistent negative eigenvalue messages often indicate buckling behavior, not numerical instability. The correct fix is an arc-length (Riks) method or implicit dynamics solver — not a tighter tolerance or different solver algorithm.
Diagnostic Tools and Techniques
Residual Monitoring
Most solvers report residual values for each equation (continuity, momentum, energy, turbulence). Monitor residuals over iterations. A well-converged residual typically drops by 3–4 orders of magnitude, from $10^{-1}$ to $10^{-4}$ or lower.
Important: If residuals drop initially but then plateau without reaching the target, the simulation may have converged numerically without reaching physical accuracy. Check monitor points for this case.
Monitor Points
Monitor points are specific quantities you track at every iteration. In CFD, these are typically forces (lift, drag), temperatures, velocities, or pressures at specific locations. In FEA, these are reaction forces, displacements, or contact pressures.
Monitor points reveal divergence that residuals hide. As the LEAP Australia guide emphasizes, a simulation with flatlined residuals can still have diverging monitor points. Always run residuals and monitor points together.
Solver Diagnostics
Different solvers provide different diagnostic tools:
- Newton-Raphson convergence plots (FEA) — show whether the solver is making progress toward equilibrium or oscillating
- Eigenvalue tracking (FEA) — negative eigenvalues indicate buckling or unstable contact
- Treedepth and R-hat diagnostics (MCMC/Bayesian solvers) — identify divergent transitions in Hamiltonian Monte Carlo
- BFMI (Basal Full Information) monitoring (Stan/MCMC) — detects regions of pathological curvature
Residual Contour Visualization
Enable expert mode in your solver and run a single iteration. Most solvers will generate residual contour plots showing where residuals cluster spatially. This reveals the domain region where divergence originates — an inlet boundary, a mesh transition zone, or a localized high-gradient region.
NaN and Inf Tracing
Identify the exact timestep or iteration where NaN or Inf values first appear. Then trace back to the previous step and examine the values that led to the NaN. This is particularly useful for identifying whether divergence is spatial (one region of the domain) or temporal (the solver stepping outside its stability region).
When to Escalate
Sometimes, a simulation diverges despite correct mesh, boundary conditions, and solver settings. In these cases, the problem may be:
- Physical impossibility — the boundary conditions or initial conditions describe a scenario that is physically impossible for the solver’s assumptions (e.g., incompressible flow with compressible materials)
- Solver limitation — the solver cannot handle the specific physics regime (e.g., multiphase flow with extreme density ratios)
- Code-level bug — a discretization error, sign mistake, or missing term in your implementation
When to escalate:
- You’ve verified mesh quality, boundary conditions, time step, and solver settings
- The simulation diverges even with a conservative time step and first-order discretization
- Monitor points diverge while residuals converge (or vice versa)
- You suspect a discretization sign error but cannot find it
In these cases, consider reaching out to solver support, consulting the literature for similar cases, or using the Method of Manufactured Solutions to verify your discretization code line by line.
Summary and Next Steps
Debugging numerical simulations is not about finding the right setting — it’s about following a systematic diagnostic workflow. The four root-cause categories (temporal, spatial, boundary/initial conditions, solver settings) give you a taxonomy to narrow the search space. The boundary-swap trick, residual contour visualization, and monitor point inspection give you practical tools that experienced practitioners use to isolate the actual cause.
Your next steps:
- When your simulation diverges, start by forcing failure faster (larger time step, fewer iterations)
- Isolate which variable diverges first using monitor points
- Change one setting at a time, verify, and document
- Use residual contour visualization to find where divergence originates spatially
- For nonlinear FEA, switch to arc-length methods when force control fails at equilibrium points
If you’re working with stiff problems specifically, the guide on stiff problems and stiff solvers covers when and why explicit methods fail. For general code debugging workflows beyond numerical diagnosis, see the Python debugging guide for scientific code.
Related Guides
- Python Debugging for Scientific Code — print statements, pdb, profiling, and testing workflows
- Stiff Problems and Stiff Solvers — when explicit methods fail and why stiff problems break solvers
- Reproducibility and Its Role in Debugging — versioned workflows for reproducible simulation debugging
- In-Situ Visualization — using in-situ visualization as an early warning system for numerical divergence
- Mesh Quality and Convergence Studies — mesh diagnostics that affect stability and convergence
Practical Checklist: What to Try First
Use this checklist when your simulation diverges. Work through each step before moving to the next.
- [ ] Verify mesh quality — skewness < 0.7, aspect ratio < 50, no negative volumes
- [ ] Check CFL number — reduce time step until $C < 0.5$
- [ ] Verify boundary conditions — inlet/outlet consistency, reference pressure, flow direction
- [ ] Set initial conditions — ensure they match boundary conditions and expected physical state
- [ ] Use first-order discretization — ramp up to second-order after stabilization
- [ ] Check under-relaxation factors — use defaults, then tune if needed
- [ ] Monitor residuals AND monitor points — flatlined residuals don’t guarantee convergence
- [ ] Generate residual contour plots — visualize where residuals cluster spatially
- [ ] For FEA: check eigenvalues — negative eigenvalues often indicate buckling, not solver error
- [ ] For FEA: switch to arc-length method (Riks) if force-controlled solver fails at equilibrium
- [ ] Verify discretization signs — a single sign error can cause catastrophic divergence
- [ ] Change only one setting at a time — always, never two simultaneously
Bottom line: Divergence is diagnosable. The four categories (temporal, spatial, boundary/initial conditions, solver settings) provide a framework. The techniques (boundary-swap, residual contours, monitor points, one-change-at-a-time) provide the tools. Use them systematically, and you’ll find the root cause — not just mask the symptom.
External Sources
The following external sources were consulted during research and are recommended for further reading:
- LEAP Australia — CFD Troubleshooting Checklist — Systematic CFD diagnostic workflow and monitor point setup
- Ansys FDTD — Troubleshooting Diverging Simulations — Boundary-swap diagnostic and dt-stability guidance
- ANSYS Fluent — Errors in Simulations — Comprehensive error categorization and solver settings
- Enterfea — Nonlinear Analysis Convergence — Arc-length methods and force vs displacement steering
- Abaqus Model Debugging — Newton-Raphson debugging and eigenvalue analysis
- SciComp StackExchange — Debugging Numerical Code — Discretization sign errors and Method of Manufactured Solutions
- NASA WIND — CFD Error Classification — CFD uncertainty and error taxonomy