← Research

Geometric multigrid fails on one-cell-thick shells, and two bugs it takes to find out

Anderson Brunsvold

Status: internal measurement, not peer reviewed. See §7. Tier: a negative result, with numbers — plus two construction findings.


1. Summary

We built a geometric multigrid V-cycle preconditioner for the momentum solve of a voxel-grid Stokes solver, verified it symmetric and positive definite by materialization, and measured it a regression at every resolution we can reach on real part geometry: 0.45× / 0.59× / 0.81× at 6k / 37k / 292k cells.

The cause is geometric and, we think, general to the application: a molded part is a shell, and at every voxel budget that is affordable its thinnest section is one cell. A one-cell wall cannot be coarsened — the coarse mask stops resembling the geometry after a single level — so the V-cycle pays in full for a correction that has no geometry left to correct.

The literature contains the complementary gap. The best available positive result for algebraic multigrid on stair-step voxel geometry (Adams et al., voxelized trabecular bone) resolves features with three to four cells across; a survey of the area concluded that no paper tests multigrid on a domain exactly one cell thick. This note supplies the one-cell measurement from the geometric side.

Two construction findings are recorded because both are silent, both were caught only by materializing the operator, and both would recur in any reimplementation.


2. The gap

A targeted literature survey of multigrid on thin and stair-step domains (conducted with AI assistance — see the colophon) returned a clear negative: no paper we could reach tests multigrid on a domain exactly one cell thick. The closest positive result is Adams et al. on voxelized trabecular bone, which resolves features with three to four cells across. We restate that finding in our own voice rather than quoting the survey, and note in §7 that the search budget was exhausted and several sources were paywalled.

The survey also drew a distinction that is essential and frequently conflated:

thin elementsthin domain
geometrycells of aspect ratio 10–10⁴cubic cells, domain 1–4 cells thick
troubleoperator anisotropy; relaxation smooths in one direction onlynothing intrinsic; only geometric coarsening breaks
literature verdictmultigrid genuinely strugglesworks well

Our case is the second. Uniform cubic voxels give an isotropic stencil; the operator is fine. What fails is the coarsening hierarchy. Most published "multigrid fails on thin problems" results are about the first column, or about shell bending operators whose conditioning catastrophe (κ ~ t⁻²) has no analogue in Stokes flow.


3. Method

Operator. -div(η grad u) on MAC face fields, masked to fluid faces, no-slip at faces touching non-fluid cells, harmonic interface viscosity.

V-cycle. Piecewise-constant prolongation (repeat each coarse value over its eight children) with its exact adjoint as restriction (sum over children); damped Jacobi smoothing (ω = 2/3, V(2,2)); a coarse cell is active when any of its eight children is; rediscretization at the coarse spacing, scaled (see §5.2). Used as a preconditioner inside the existing conjugate-gradient momentum solve, never as a solver — so a poor V-cycle can cost iterations but cannot weaken a convergence claim.

Geometries. (a) A solid n³ box, n = 8…64. (b) A shelled 210 × 50 × 24 mm tray with a 3.5 mm wall, voxelized at three budgets, timing one complete Stokes solve.

Provenance of the tray: a synthetic fixture authored for this engine's own test corpus, generated parametrically and dimensioned from the pipeline's own constraints (mold-base tie-bar coverage, clamp envelope, ejector-pin threshold) — realistic in character, free of any third-party origin.


4. Results

On a solid box the preconditioner works, CG iterations to 1e-8:

n8163264
plain CG3363119234
V-cycle10141926

9× at n = 64, with the advantage widening. Note this is not mesh independence — piecewise-constant transfer is aggregation multigrid and degrades slowly with h (~h^-0.45 against plain CG's h^-1). We do not claim otherwise.

On real part geometry it is a regression, one Stokes solve, plain vs preconditioned:

cellscells across wallplainV-cycleratio
6,10213.6 s8.1 s0.45×
36,908147.1 s79.9 s0.59×
291,7985624 s775 s0.81×

Confirmed at whole-fill level as well (2.3–2.9× slower on every row of a four-point cost curve).

The trend is the finding. At one cell across the wall the V-cycle is less than half the speed of no preconditioner at all. At five cells — the first budget where coarsening becomes possible, and one that costs 292k cells on a small part — it reaches 0.81×: approaching parity, never winning.

A second observation that redirects the effort. The same runs show the outer Uzawa iteration count growing on real geometry (116 → 180 → 214) where it is flat on the box (12 → 15). Preconditioning the inner solves cannot touch an outer loop that is itself scaling. This is treated separately in the companion note on Schur-complement degradation.


5. Two construction findings

Both are silent at uniform viscosity, which is the configuration most test suites use, and both were caught by materializing the operator into a dense matrix and inspecting it directly.

5.1 Coarsening the interface coefficients independently destroys symmetry

The six directional coefficient arrays of the stencil are not six independent fields. They are built so that the coefficient on a face is shared by the two rows it couples, c₊[i] == c₋[i−1], and that identity is precisely why the operator is symmetric. Averaging each direction array separately over a coarse cell's children breaks it — the mean of c₊ over one aggregate has no relation to the mean of c₋ over its neighbour — and the coarse operator comes out asymmetric, so the preconditioner violates CG's precondition.

Measured: 28% asymmetry on a 16³ grid at a 400× viscosity contrast, and exactly zero at uniform viscosity, where all six arrays are equal.

The fix is to coarsen the field the coefficients are derived from — the viscosity — and re-derive, so the shared-face identity is reconstructed at every level by the same function that guarantees it on the fine grid.

5.2 The aggregation scale is exactly 2^(d+1), and omitting it inverts the sign of the benefit

With piecewise-constant prolongation and its exact adjoint, the Galerkin coarse operator PᵀAP is exactly 16× a rediscretization at spacing 2h in three dimensions — four fine faces couple across each coarse face while the rediscretized coefficient carries only 1/(2h)². We measured the ratio at 16.000 on both the diagonal and the off-diagonal.

Unscaled, the coarse correction is 16× too large. All variants, re-measured at one identical n-set, tolerance (1e-8) and right-hand side so they are directly comparable:

n8163264
plain CG3363119234
V-cycle, scaled (as built)10141926
V-cycle, unscaled2560171508
V-cycle, linear transfer14254690

The unscaled V-cycle is not uniformly worse than no preconditioner — it beats plain CG at n = 8 and 16. It becomes worse from n = 32 and is more than twice as bad at n = 64, and crucially its iteration count grows faster than plain CG's (20× across the range against plain's 7×), so the defect worsens with refinement. That is the failure mode: a wrongly-scaled coarse correction degrades rather than merely underperforming.

A caution for anyone tempted by better interpolation: linear transfer (binomial-smoothed) beats plain CG but is consistently worse than the scaled piecewise-constant V-cycle — 14 / 25 / 46 / 90 against 10 / 14 / 19 / 26. The diagnostic explains why: with smoothed transfer the Galerkin operator is no longer proportional to a rediscretization at all (ratio 1.95 on the diagonal, −0.78 off it), so no single scale reconciles them and the coarse level solves a different problem. Making linear transfer work requires forming PᵀAP explicitly.

(An earlier draft of this note quoted "25 / 60 / 93" as the plain-CG baseline here. That figure was wrong — it came from a different script with a different n-set and was never plain CG. The table above supersedes it and every row was re-run together.)


6. What this implies

  1. For thin-shell domains at production voxel budgets, geometric multigrid is the wrong instrument, because the hierarchy dies at the first level.
  2. The literature's positive results for algebraic multigrid on stair-step voxel geometry (Adams et al., 3–4 cells across a trabecula, constant iteration count, 57M elements in 9 minutes against 5 weeks for the incumbent) are promising precisely because AMG coarsens by matrix connectivity rather than geometry — on a one-cell wall, aggregation can only coarsen in-plane, which is automatic semicoarsening and the exact capability the geometric hierarchy lacks. This remains untested at one cell and is the obvious next experiment.
  3. Whichever multigrid is used, it should be applied to the operator that actually limits the solve. On real geometry that is the pressure Schur complement, not the momentum block.

7. Threats to validity

  1. One implementation, one solver. The mechanism (a one-cell wall cannot be coarsened) follows from the geometry, but the magnitudes do not transfer.
  2. The tray timings use a single-front-cell pressure pin, which is not how a real fill pins the front. The plain-vs-V-cycle comparison is valid because both sides are identical in that respect, and the fill-level measurement agreed closely with the proxy (2.2× predicted vs 2.26× measured at one budget) — but absolute iteration counts from that setup should not be quoted.
  3. Aggregation multigrid, not the best available. A smoothed-aggregation or properly-Galerkin scheme might behave differently. We showed one alternative (linear transfer) is worse in this construction, not that no transfer works.
  4. Coarsening rule not swept. "A coarse cell is active if any child is"; the conservative alternative ("all children") was not measured.
  5. The literature claim is a negative, inherited from a survey whose search budget was exhausted. Absence in our search is not proof of absence.

8. What would make this publishable

  • The algebraic counterpart: pyamg or equivalent on the same operators and the same one-cell geometry, which would close the gap from both sides and either confirm or refute the in-plane-aggregation argument of §6.2.
  • A characterization of hierarchy collapse: how coarse-mask fidelity (some measure of geometric agreement) degrades with level as a function of cells-across-wall, which would predict the onset rather than observe it.
  • Measurement at more than one wall resolution — we have one (1) and five (5) cells and nothing between, so "approaching parity" rests on two points.
  • Independent reproduction, and a prior-art pass for the 2^(d+1) aggregation scale, which is standard aggregation-multigrid material and is here measured rather than cited.

9. Sources

  • Adams, Bayraktar, Keaveny & Papadopoulos, "Applications of Algebraic Multigrid to Large-Scale Finite Element Analysis of Whole Bone Micro-Mechanics" (SC'04) — smoothed-aggregation AMG on micro-CT voxel meshes with stair-step surfaces and 3–4 elements across a trabecula; constant iteration count; 57M elements in under 9 minutes.
  • Trottenberg, Oosterlee & Schüller, Multigrid (Academic Press, 2001) — §5.1–5.2 on the failure of pointwise relaxation under anisotropy, and the distinction between anisotropic operators and thin domains.
  • Tamstorf, Jones & McCormick, "Smoothed Aggregation Multigrid for Cloth Simulation" (SIGGRAPH Asia 2015) — the source of the widely-cited "multigrid fails on thin shells" claim; note their failure is strong inter-variable coupling in bending-dominated elasticity, and they show standard AMG works once components are decoupled. It does not describe a Stokes momentum block.
  • Clevenger & Heister, "Comparison Between Algebraic and Matrix-free Geometric Multigrid for a Stokes Problem with Variable Viscosity", arXiv:1907.06696 — for the assembly-cost accounting that makes matrix-free attractive when geometric multigrid works.

The implementation is retained in-tree but not wired into the solver, with this negative result in its module docstring, so the two construction findings of §5 are not lost.


Colophon — authorship and AI assistance

Author (responsible natural person): Anderson Brunsvold, FutureMold.

AI assistance disclosure. This note was produced with substantial AI assistance, and the extent is material to how it should be read:

  • Literature search and synthesis: AI-conducted. The surveys behind the "no published source" claims were run by AI agents. Their search budgets were exhausted and several primary sources were paywalled and read only via abstract or secondary description. The novelty claims rest on those searches and inherit their limits — see the threats-to-validity section.
  • Measurements: AI-executed, human-directed. The experiments were designed, scripted and run by an AI agent (Claude) working in the FutureMold repository, under direction from the author. The raw outputs are in the session record and the repository's reference documentation.
  • Drafting: AI-written, human-reviewed. The prose is AI-drafted. It was reviewed by the author, and that review caught substantive errors — including, in this series, a claim that contradicted its own data table and a headline figure that conflated two mechanisms. Corrections are marked in place rather than silently applied.
  • Not independently reproduced. No result here has been reproduced by a second implementation or a second party.

Standing caveat. Absence of a result in an AI-conducted literature search is weak evidence of absence in the field. Any novelty claim in this note should be re-checked against paywalled and non-indexed sources before it is relied upon.

Literature searches, experiment design, scripting, execution, and drafting were AI-performed (Claude, Anthropic) under the direction of the named author, who reviewed the work and takes responsibility for it. Each note carries a fuller provenance section, including the limits that AI-conducted searches place on its novelty claims.

Clearance: no-customer-data.