# ACCESS-AIS3 inversion work log
Diagnosis and repair of the present-day fricti…on / rheology inversion.
Goal: a well-initialised Antarctic state with **grounded velocity RMSE < 50 m/yr** and
spatially-white residuals (the reference result achieved with shelf-only rheology
inversion plus a friction inversion, grounded B from temperature).
---
## 1. Headline progress
| configuration | grounded RMSE (m/yr) | notes |
|---|---|---|
| null model (`RMS(v_obs)`) | 310.6 | baseline to beat |
| original setup | 323.5 | **worse than null** — model was effectively static |
| joint friction + B, 70 iters | 171.2 | 70 % of grounded within ±50 |
| **friction-only, 200 iters** | **110.6** | **83 % within ±50; median residual −0.1** |
Residual structure at 110.6:
```
median -0.1 mean -4.9 std 110.5 |res| < 50 : 83 %
grounded, >100 km inland : RMSE 55.3 ratio 1.00 <- already at target
grounded, 0-5 km from GL : RMSE 185.9 ratio 1.18 <- grounding zone
Transantarctic Mountains : RMSE 269.6 ratio 2.00 <- worst region
shelves (out-of-sample) : RMSE 237.9 ratio 1.03
```
**The error is highly concentrated: 9.01 % of grounded vertices carry 90.6 % of the total
squared error.** The interior is already at the target; what remains is a small number of
identifiable regions.
---
## 2. Bugs found and fixed
### 2.1 Friction coefficient initialised in a dead zone
With Budd `p=q=3`, `tau_b = C^2 · N · u^(1/3)`, so **`u ∝ C^-6`**. The initial guess
`C = 100` put the model six orders into a region where `∂u/∂C ≈ 0`; grounded ice was
static (median velocity 0.0 m/yr against 71.9 observed).
Fixed: `C_init = 1.8`, bounds `[0.1, 10]`, `control_scaling_factors = 1.8`.
**Consequence:** every earlier tuning result is void. The "best config" from the original
regularisation campaign was optimised against `fastGrnd RMSE ≈ 498`, which was literally
`RMS(v_obs) = 499.1` — the metric never measured model skill because the model did not move.
### 2.2 `rheology_B` was a uniform placeholder
`cuffey(273.15 − 20) = 2.027e8` everywhere, collapsing a ~17× real spatial viscosity range
(8.5e7 → 1.43e9 over the RACMO range) to one value. 51 % of grounded vertices were wrong by
> 20 %. No friction inversion can correct this — friction acts on the bed, the error is in
the ice.
Fixed in `config/ais_0.1_param.py`: B derived from the RACMO surface-temperature field.
Caveat retained in the code: surface temperature is a proxy for depth-averaged, biased stiff.
### 2.3 Wrong cost-function ID for friction regularisation
The friction L-curve used **502** (`RheologyBbarAbsGradient`), which is not a control in a
friction inversion — the term was inert and the L-curve flat. Fixed to **501**
(`DragCoefficientAbsGradient`) in `config/ais_0.1.py`.
### 2.4 Element-based friction mask stripped drag from grounded ice
The zero-friction mask flagged *all* vertices of any element containing one floating
vertex, removing basal drag from a band of genuinely grounded ice at the grounding line.
92.7 % of the velocity blow-up lived in that band.
Fixed: vertex-based test `ocean_levelset < 0`.
### 2.5 Absolute-misfit weighting ignored the slow interior
`cf101 = 1000 / cf103 = 0.1` weights absolute error, so the slow interior contributed
almost nothing. Fixed: log-weighted `cf101 = 10 / cf103 = 100`.
---
## 3. Dead ends (with evidence)
### 3.1 Friction regularisation (501) is unusable
Sweeping 501 upward produced **five byte-identical non-starters** — `C moved 0.0 %`,
roughness exactly 0.0000.
Root cause: the friction field steps `1.8 → ~0` across the grounding line because floating
C is pinned. At the initial guess that single ring of elements (2.3 % of the mesh) supplied
**100.0000 %** of `|grad C|^2`; grounded ice contributed exactly **0**. Since floating C is
bounded to `[0,0]`, the penalty is irreducible — a fixed cost with a gradient ~10^7 times
the misfit gradient, which kills the m1qn3 line search for every `501 >= 1e2`.
Masking the grounding-line ring out of the functional removed the artefact (verified:
`|grad C|^2` at init `1.423e5 → 0`) and the runs then iterated 17–33 steps. But they still
stalled early with a **ruined fit** (RMSE 2900–3300, worse than null; C moved only ~25 %).
The smooth fields were smooth because they had barely moved, not because they were coherent.
**Conclusion:** the data term wants sharp features at ice-stream margins, the smoothing term
fights them, and the line search cannot reconcile the two. Abandoned; 501 left inert.
### 3.2 Joint friction + B inversion is worse than friction-only
| | abs misfit (101) | log misfit (103) | grounded RMSE |
|---|---|---|---|
| joint (stopped, dxmin @ 70) | 684.1 | 1.399e6 | 171.2 |
| friction-only (200 iters) | 615.9 | 6.852e5 | **110.6** |
B freedom let the optimiser **stop early at a worse velocity fit**. A previously-reported
"invariant fast-ice deficit of 0.70–0.74 across every configuration" was an artefact of
this early stopping — friction-only running to 200 iterations moved it to **0.87** on its own.
### 3.3 Dropping the thin-ice floor entirely
Removing the floor *and* excluding thin ice (H < 50 m) from the misfit gave
**RMSE 1917.4, maxvel 1.6e5** — a runaway. The lower cost function was an illusion created
by the smaller misfit set: excluding thin ice did not fix it, it merely stopped constraining
it. Thin ice must stay **in** the misfit.
---
## 4. The Transantarctic Mountains failure
Confirmed geographically: 30,739 grounded vertices, **RMSE 269.6, ratio 2.00** (2× too fast).
Residual binned by physical regime (baseline joint run):
| regime | n | RMSE | ratio |
|---|---|---|---|
| N < 1e5 Pa | 18,355 | 678.9 | **14.56** |
| H floored to 100 m | 71,578 | 368.7 | 2.19 |
| slope > 0.1 | 88,152 | 357.7 | 2.56 |
| H 0–150 m | 85,783 | 342.8 | 1.93 |
| bed 1000–4000 m | 66,695 | 265.3 | 1.58 |
| *real ice, H > 100 m* | 1,158,642 | 150.8 | 1.01 |
Every failing regime is **too fast**, and they all describe thin ice on steep, high bedrock.
### Root cause: a self-inflicted geometry error (in the experiment scripts, not production)
The experiment scripts floored thickness at 100 m and then rebuilt
`surface = bed + H`, which **raised the surface** wherever real ice was thinner and corrupted
the **surface slope**. Since `tau_d = rho·g·H·grad(s)`, that fabricates the driving stress —
worst exactly in the mountains.
Production (`config/ais_0.1_param.py:62`) does it correctly, flooring at 10 m and absorbing
the change into the **base**, leaving the observed surface intact.
### Hypothesis that was refuted
Effective-pressure starvation was proposed as the cause and **disproved**: above sea level
`N/overburden = 1.016`, i.e. already correct. Genuine N-starvation exists but only in the
small 18k-vertex group (`N/overburden = 0.007`) responsible for the 14.56× outliers.
---
## 5. Shelf rheology is inconsistent with the coupled model
`md.extract()` imposes Dirichlet BCs from `vel_obs` on **every** new boundary node
([`Model.py:558-560`](../../../pyISSM/src/pyissm/model/Model.py)). For the floating-only
extraction used by the shelf inversion, that clamped **both** the grounding line and the
**calving front** to observations. The coupled run instead uses a **Neumann** shelf front.
So B was fitted while artificially restrained; released, the shelf accelerates:
| distance from GL (km) | RMSE | median residual | ratio |
|---|---|---|---|
| 0–5 | 227.4 | −2.9 | 1.08 |
| 15–30 | 260.9 | +7.5 | 1.03 |
| 60–120 | 343.4 | +108.6 | 1.19 |
| 120–250 | 296.5 | **+226.9** | **1.26** |
The error **grows** with distance from the grounding line — i.e. it is worst at the calving
front, not inherited from grounded ice.
**A shelf-only domain cannot be made BC-consistent.** Freeing the calving front on floating
ice (which has no basal drag) leaves the system nearly singular; the solve fails outright
(`Recovery solver failed`, 0 iterations). The original workflow only worked because the
all-Dirichlet extraction was propping it up — the very artefact that biased B soft.
Therefore shelf B must be inverted **inside the coupled domain**, where the front BC and the
grounding-line inflow are both exactly what production uses.
---
## 6. Production config status
Committed:
- `ais_0.1_param.py` — temperature-derived `rheology_B`
- `ais_0.1.py` — 502 → 501 fix, `extract_friction_inversion_domain()` (Neumann fronts),
`friction_coupling` variable
**Still carrying known-wrong values:**
| parameter | production | validated |
|---|---|---|
| `friction.coefficient` init | 10.0 | **1.8** (bounds `[0.1, 10]`) |
| `friction_cf101 / cf103` | 1000 / 0.1 | **10 / 100** |
The `1000 / 0.1` weighting comes from the invalidated campaign (§2.1). Not yet promoted
because the configuration is still being revised.
Correct as-is: thickness floor 10 m into the base; `N_floor_frac = 0.07`; Budd `p=q=3`;
`coupling = 3`.
---
## 7. Open items
1. **N-starved vertices** — 18,355 with `N/overburden = 0.007`, running 14.56× too fast.
2. **Grounding zone** — grounded ice 0–5 km from the GL is 18 % too fast (RMSE 185.9 vs 55.3 inland).
3. **Steep-slope tail** — slope > 0.1, 2.56× too fast.
4. **C field is noisy** — roughness 0.164, 3 % at bounds. Not blocking the fit, but 501
cannot fix it; would need a different smoother (e.g. inverting `log C`).
5. **Promote validated parameters** to the production config.
6. **Forward / relaxation validation** — not yet run.
## 8. Practical notes
- `qcat -o <jobid>` streams a running job's stdout; ISSM writes `.outlog` only at exit.
- 192 cores (4 nodes, 760 GB) gives **6.5×** the throughput of 32 cores: 35.2 vs 5.45
iterations/hour, essentially linear scaling. Memory *rises* with rank count
(130 GB @ 32 → 450 GB @ 192) due to halo duplication.
- Friction-only inversions need ≥ 200 iterations; the 110.6 run hit the maxsteps ceiling
still improving.
- Cost-function values are **not comparable** between runs with different misfit masks —
always compare grounded RMSE on a common mask.