Convergence Study
Goal
Provide script-driven convergence analysis for verifying that numerical solutions converge at the expected rate as the mesh or timestep is refined.
Requirements
- Python 3.10+
- No third-party packages — scripts use only the Python standard library (
math).
| Input |
Description |
Example |
| Grid spacings |
Sequence of mesh sizes (coarse to fine) |
0.4,0.2,0.1,0.05 |
| Timestep sizes |
Sequence of dt values |
0.04,0.02,0.01 |
| Solution values |
QoI at each refinement level |
1.16,1.04,1.01,1.0025 |
| Expected order |
Formal order of the numerical scheme |
2.0 |
| Safety factor |
GCI safety factor (1.25 default) |
1.25 |
Script Outputs (JSON Fields)
| Script |
Key Outputs |
scripts/h_refinement.py |
results.observed_orders, results.mean_order, results.richardson_extrapolated_value, results.convergence_assessment |
scripts/dt_refinement.py |
Same as h_refinement but for temporal convergence |
scripts/richardson_extrapolation.py |
results.extrapolated_value, results.error_estimate, results.observed_order |
scripts/gci_calculator.py |
results.observed_order, results.gci_fine, results.gci_coarse, results.asymptotic_ratio, results.in_asymptotic_range, results.extrapolated_value, results.notes |
Workflow
- Run grid/timestep refinement study with at least 3 levels
- Compute observed convergence order with
h_refinement.py or dt_refinement.py
- Compare observed order to expected order of the scheme
- Estimate discretization error via Richardson extrapolation
- Report GCI for formal solution verification using
gci_calculator.py
- Document convergence results and any anomalies
Decision Guidance
Do you have 3+ refinement levels?
+-- YES --> Run h_refinement.py or dt_refinement.py
| +-- Observed order matches expected? --> Solution verified
| +-- Order too low? --> Check: pre-asymptotic, coding error, insufficient resolution
| +-- Order too high? --> Check: superconvergence or cancellation effects
+-- NO (only 2 levels) --> Use richardson_extrapolation.py with assumed order
(less reliable without order verification)
CLI Examples
# Spatial convergence with 4 grid levels
python3 scripts/h_refinement.py --spacings 0.4,0.2,0.1,0.05 --values 1.16,1.04,1.01,1.0025 --expected-order 2.0 --json
# Temporal convergence with 3 timestep levels
python3 scripts/dt_refinement.py --timesteps 0.04,0.02,0.01 --values 2.12,2.03,2.0075 --expected-order 2.0 --json
# Richardson extrapolation with assumed 2nd-order
python3 scripts/richardson_extrapolation.py --spacings 0.02,0.01 --values 1.0032,1.0008 --order 2.0 --json
# GCI for 3-mesh verification
python3 scripts/gci_calculator.py --spacings 0.04,0.02,0.01 --values 1.0128,1.0032,1.0008 --json
Error Handling
| Error |
Cause |
Resolution |
spacings and values must have the same length |
Mismatched input arrays |
Provide equal-length lists |
At least 2 refinement levels required |
Too few data points |
Add more refinement levels |
Exactly 3 refinement levels required |
GCI needs 3 levels |
Provide fine/medium/coarse |
Oscillatory convergence detected |
Non-monotone convergence |
Check mesh quality or scheme |
Interpretation Guidance
| Scenario |
Meaning |
Action |
| Observed order matches expected |
Strongest evidence of asymptotic range |
Report GCI, extrapolate |
| Observed order < expected |
Pre-asymptotic or coding bug |
Refine further or debug |
| Negative observed order |
Solution diverging |
Check implementation |
| GCI asymptotic ratio near 1.0 |
See caveat below |
Confirm with order comparison |
| GCI asymptotic ratio far from 1.0 |
Not in asymptotic range |
Refine further |
Asymptotic-ratio caveat (constant refinement ratios). When the refinement
ratios are equal (r21 == r32, the common case), the asymptotic ratio
AR = GCI_coarse / (r^p * GCI_fine) reduces algebraically to f1/f2. It then
only measures the relative gap between the two finest QoI values — not whether
the data follow the assumed power law — so an AR near 1.0 can give false
reassurance even when the observed order is far from expected. The
gci_calculator.py JSON emits a notes entry flagging this. For a real
asymptotic-range determination:
1. Compare the observed order p to the scheme's theoretical/expected order — a
match is the meaningful evidence of being in the asymptotic range.
2. For stronger verification, use 4+ systematically refined grids and check that
the observed order is consistent across successive grid triplets
(h_refinement.py reports one order per triplet plus mean_order).
Verification checklist
- [ ] Used >= 3 systematically refined grids/timesteps so
h_refinement.py / dt_refinement.py can report a results.mean_order; recorded the per-triplet results.observed_orders (a single-pair Richardson run does not verify order).
- [ ] Recorded
results.mean_order and confirmed results.convergence_assessment reads PASS (observed order within 10% of the scheme's expected order); a FAIL or unknown means the result is not yet verified.
- [ ] Confirmed
results.in_asymptotic_range is true and that no notes entry reports pre-asymptotic (>50% order variation), negative/non-positive order, or zero error differences before quoting any GCI or extrapolated value.
- [ ] Checked refinement ratios
r21/r32 from gci_calculator.py are >= 1.3 (round-off noise floor) and that no Oscillatory convergence detected error was raised.
- [ ] Recorded
results.gci_fine (with the safety factor used: 1.25 for >= 3 grids with verified order, 3.0 for 2 grids with assumed order) as the reported discretization uncertainty, plus results.extrapolated_value as the best estimate.
- [ ] For constant refinement ratios, did NOT treat
asymptotic_ratio near 1.0 as proof of asymptotic range (it degenerates to f1/f2); confirmed the asymptotic range via observed-order-vs-expected and read the gci_calculator.py notes caveat.
Common pitfalls & rationalizations
| Tempting shortcut |
Why it's wrong / what to do |
| "Two grids agree closely, so it's converged" |
Two levels cannot estimate observed order. Run h_refinement.py/dt_refinement.py with >= 3 levels; a 2-grid Richardson run uses an assumed order and needs safety factor 3.0, not 1.25. |
| "The asymptotic ratio is ~1.0, so we're in the asymptotic range" |
With constant refinement ratios AR = GCI_coarse/(r^p*GCI_fine) reduces to f1/f2 and only measures the gap between the two finest QoI values. Verify the asymptotic range by comparing observed order to the expected order (and 4+ grid consistency). |
| "GCI_fine is tiny, so the solution is grid-independent" |
A near-zero GCI can also mean the QoI is insensitive to refinement or the differences are in round-off noise (ratios < 1.3). Confirm in_asymptotic_range is true and refinement ratios are >= 1.3 first. |
| "Observed order is higher than expected, even better" |
Order well above the formal order usually signals superconvergence or error cancellation, not extra accuracy. Treat it as a flag (convergence_assessment is FAIL when >10% off) and verify with more grid levels. |
| "The script printed an extrapolated value, so use it" |
When the observed order is non-positive the solution is diverging and h_refinement.py returns richardson_extrapolated_value = null with a diverging note. Do not quote an extrapolated value or GCI in that case. |
| "Implicit/stable solver, so any timestep is fine for the study" |
Temporal stability is not temporal accuracy. dt_refinement.py still needs >= 3 systematically reduced timesteps to recover the scheme's order; a too-coarse dt sequence stays pre-asymptotic. |
Security
- All numeric parameters (
spacings, timesteps, values, expected-order, order) are validated as finite positive numbers
- Comma-separated value lists are length-matched (spacings and values must have equal length) and capped at 10,000 entries
- GCI calculator enforces exactly 3 refinement levels; Richardson extrapolation requires at least 2
- Safety factor is validated as a finite number not less than 1.0 (Roache uses Fs in {1.25, 3.0})
File Access
- Scripts read no external files; all inputs are provided via CLI arguments
- Scripts write only to stdout (JSON output); no files are created unless the agent explicitly uses the Write tool
- Bash: Used to execute the four Python analysis scripts (
h_refinement.py, dt_refinement.py, richardson_extrapolation.py, gci_calculator.py) with explicit argument lists
- Read: Used to inspect script source and reference documentation
Safety Measures
- No
eval(), exec(), or dynamic code generation
- All subprocess calls use explicit argument lists (no
shell=True)
- Scripts use only Python standard library (
math module); no pickle loading or deserialization of untrusted data
- Minimal tool surface (Bash and Read only) limits the agent's ability to modify the filesystem
References
references/convergence_theory.md - Formal convergence order, log-log analysis, asymptotic range
references/gci_guidelines.md - Roache's GCI method, ASME V&V 20, safety factors