A 2-D diffusion eigenproblem with a one-parameter coefficient

This is the reference example for the whole method: a scalar diffusion eigenproblem on the unit square whose coefficient depends on a single scalar parameter \mu. The spatial problem is two-dimensional; the parameter it is swept over is one-dimensional, \mu \in [0.4, 1]. Those two dimensions are unrelated and worth keeping apart — everything below is “2-D in space, 1-D in parameter”.

The page is built in the order the method itself works. First the shared groundwork: the model, the matching operation on a single subinterval, and the discretisation. Then each of the two operating modes gets a full, self-contained treatment, run to its own convergence: Window mode — the faithful reproduction of the published 1-D grid — followed by FirstK mode — a fixed-band approximation. The page closes with an honest side-by-side comparison of the two.

Related pages

Why this example

It is the smallest problem that shows the full phenomenon the method exists to handle: eigenvalue branches that approach and exchange character as the parameter varies, and an adaptive grid that must find those events on its own and refine around them. It is small enough to solve densely for ground truth, so every claim the adaptive loop makes can be checked against a known answer — which is why the accuracy and robustness studies (ROM vs dense, tolerance perturbation) use it as their reference. If the method works here, the machinery is sound; larger problems (a 2-D parameter box, 6-D elasticity) add scale, not new ideas.

The model

Find (\lambda, u) \in \mathbb{R} \times H_0^1(\Omega) on the unit square \Omega = [0,1]^2 such that

-\nabla \cdot \bigl(c(\mu)\,\nabla u\bigr) \;=\; \lambda\,u \quad \text{in } \Omega, \qquad u = 0 \text{ on } \partial\Omega .

In weak form this is the generalised eigenproblem: find (\lambda, u) \in \mathbb{R} \times H_0^1(\Omega) with

a(u, v;\,\mu) \;=\; \lambda\, b(u, v) \qquad \forall\, v \in H_0^1(\Omega), a(u, v;\,\mu) = \int_\Omega \bigl(c(\mu)\,\nabla u\bigr)\cdot \nabla v \, \mathrm{d}x, \qquad b(u, v) = \int_\Omega u\,v \, \mathrm{d}x .

The coefficient is a global, space-constant 2\times 2 matrix:

c(\mu) \;=\; \begin{bmatrix} \mu^{-2} & 1 \\[2pt] 1 & 0.7^{-2} \end{bmatrix}, \qquad \mu \in [0.4,\, 1] .

Only the (1,1) entry varies with \mu; the other three are fixed. Because the weak form is affine in the entries of c, the stiffness matrix decomposes into four fixed tensor-component matrices K_{ij} (one per entry of c), assembled once:

K(\mu) \;=\; \sum_{i,j \in \{1,2\}} c_{ij}(\mu)\, K_{ij} \;=\; \mu^{-2}\,K_{11} \;+\; \bigl(K_{12} + K_{21} + 0.7^{-2}\,K_{22}\bigr).

Each parameter point then costs four cached component contractions plus one eigensolve — no re-meshing, no per-\mu quadrature. This affine structure is what the surrogate construction relies on.

A single subinterval

Before any sweep, consider the method on a single subinterval [\mu_1, \mu_2] = [0.4, 0.7] — the atomic operation the adaptive loop repeats on every edge of its grid, and the reference’s worked §3.3.1 example. Understand one edge and the whole loop is just “do this everywhere, bisect where it fails”.

Inside the spectral window I_\lambda = [0, 270] there are four eigenvalues at \mu_1 = 0.4 and nine at \mu_2 = 0.7 (the reference reports exactly this band count). Two steps:

  1. A-priori matching. Build the cost matrix D_{j\ell} pairing each \mu_1-band j with each \mu_2-band \ell — a weighted sum of eigenvalue distance (weight 1) and eigenvector misalignment (weight 200, measured by the M-inner-product principal angle). The optimal (Hungarian) assignment matches the four \mu_1-bands to \mu_2-bands \sigma^\star = (1, 2, 5, 8) (1-indexed).
  2. A-posteriori certification. Form the projection matrix \Pi_{j\ell} = v_j(\mu_1)^\top M(\bar\mu)\, v_\ell(\mu_2) on the matched core and apply the relative threshold t_\pi = 0.21 (the reference’s stated value for this worked subinterval). The first matched pair is uniquely supported → certified; the second has a competing off-pair entry → refine; the higher pairs sit in a shared near-degenerate cluster → cluster.

Subinterval matching on [0.4, 0.7]. Left: cost matrix D (4 bands at \mu=0.4 × 9 at \mu=0.7); the optimal assignment \sigma^\star = (1,2,5,8) is outlined in red. Right: |\Pi| on the matched core with each pair’s t_\pi = 0.21 verdict (green certified, red refine, orange cluster).

A refine verdict on any edge tells the loop to bisect it and re-test the two halves. Both sweeps below are exactly this operation on an adaptively growing grid — they differ only in how many bands each grid point carries.

Discretisation

Both runs share the same spatial discretisation and parameter domain; only the spectral bookkeeping (and so the tolerances tied to it) differs between modes.

Quantity Value Role
Mesh \texttt{mesh\_n} = 64, P1 triangulation spatial discretisation (rebuild choice; fine enough that the crossings below are mesh-independent)
Parameter box \mu \in [0.4,\, 1.0], d_{\text{param}} = 1 parameter domain
Max refinement levels 8 (Window terminates at L = 3, FirstK at L = 1) iteration cap

The two operating modes split here:

  • Window mode keeps whichever eigenvalues fall inside a fixed spectral window I_\lambda = [0, 270], so per-node rank is ragged (4 bands at \mu = 0.4, 9 at \mu = 0.7). This is how the reference 1-D example was actually computed.
  • FirstK mode keeps a fixed number n_{\text{eigs}} = 6 of the smallest bands at every node. This is the rebuild’s fixed-band approximation.

Each mode is treated in full below, on its own honest tolerances, before they are compared.

Window mode — the faithful reproduction

This is the run that reproduces the reference 1-D grid, because the §3.3.2 example was itself computed in window mode. Snapshots are ragged: each parameter point keeps whatever eigenpairs fall inside the fixed window

I_\lambda = [0,\, 270]

(the reference’s stated window — the band counts it reports, four eigenvalues at \mu = 0.4 and nine at \mu = 0.7, are reproduced exactly). The verdicts use the one a-posteriori (Hungarian) path with anchored-absolute pruning at (t_\pi = 0.21,\, t_\lambda = 0.085): t_\pi = 0.21 is the §3.3.1 subinterval value above, and t_\lambda = 0.085 is the MATLAB clustering tolerance — the pair is validated by reproducing the published L=0\ldots L=3 grid exactly (below). A surplus band that merely enters or leaves the window between adjacent nodes is certified by a row-\ell_2-norm test (the reference’s Algorithm 6 “certify, don’t refine” gate) rather than forcing refinement; those enter/exit events are recorded for introspection but do not drive the grid.

Reference spectrum

The ground truth: a dense \mu-sweep, solved directly on a window [0, 400] — each \mu-point returns exactly the eigenpairs inside that band (a ragged window solve, nothing above 400 computed, nothing clipped) — and drawn in character-tracked branch order: each curve follows one eigenfunction (matched by mass-overlap across the sweep), so where two branches exchange energetic order the figure draws a genuine crossing instead of sorting it away. The parameter is one-dimensional, so every crossing here is avoided — the gap stays positive, there is no true degeneracy (that needs two free parameters). But an avoided crossing still exchanges eigenfunction character: the two branches approach, swap which one is higher, and separate, so their curves cross. This example carries four such swaps — bands 3↔︎4 (near \mu \approx 0.45 and again \approx 0.86), 5↔︎6 (\approx 0.72), and 7↔︎8 (\approx 0.81) — visible as X crossings in the low bands. The grey strip is the algorithm’s kept window I_\lambda = [0, 270]; the sweep window runs a little higher so the bands that descend into the kept window as \mu grows are visible entering from above. The adaptive loop has no access to this; it must rediscover the crossings from local edge tests alone.

Reference spectrum, solved directly on the window [0, 400] (no clipping); the grey strip is the algorithm’s kept window I_\lambda=[0,270]. Branches are character-tracked, so the four avoided crossings show as X crossings where two curves exchange order. The band-7/8 avoided crossing near \mu \approx 0.81 — where two branches approach to a minimum gap of \approx 13.7 and separate, exchanging eigenfunction character without touching (an avoided crossing: the gap never reaches zero) — is the event refinement concentrates on. Bands descend into the window from above as \mu grows.

Adaptive grid evolution

Drag the slider (or use ←/→ when focused) to step through L = 0..3. Nodes are coloured by insertion level; the loop bisects only edges flagged for refinement, so nodes pile up around the crossing region. The grid grows 3 \to 5 \to 7 \to 8 nodes.

Adaptive grid per refinement level (Window)
Adaptive grid (Window) at L=0

Per-level surrogate

Solid coloured curves are the surrogate, reconstructed from the greedy grid’s own stored matchings (not re-matched) by the same fixed-slot reshuffle used for the reference spectrum; dashed grey is the dense reference. Branch k of the surrogate keeps one colour across every level and follows branch k of the reference through a crossing. Scroll through L = 0..3 to watch the surrogate converge as the grid refines. Where a band drops out of [0, 270] between adjacent nodes the window surrogate shows an honest gap in a fixed slot — the broken curve is one branch keeping its colour across a window exit/entry, never a connector to an unrelated branch.

Per-level surrogate spectrum (Window)
Piecewise-linear surrogate (Window) at L=0

A-posteriori projection matrices

Each panel shows eigenvalues scattered in the (\mu, \lambda) plane at the grid nodes active at that level. Solid coloured lines are the a-priori Hungarian-matched pairs — the algorithm’s candidate for “same branch, step left.” Black dashed lines are alternative matchings: off-diagonal entries |\Pi_{j\ell}| that survive the pruning threshold but are not the primary match. Their presence at an edge is the direct REFINE signal — it says “I cannot confidently assign a unique partner to this band.”

A-posteriori summary (Window) at L=0
A-posteriori summary (Window) at L=0

At L = 0, dashed lines cluster near the avoided crossing (\mu \approx 0.8), flagging those subintervals for refinement. By L = 3, all dashed lines are gone — every subinterval is certified and the run terminates.

The avoided crossing at [0.775,\,0.85]

This is what the whole apparatus exists to catch. Between \mu = 0.775 and \mu = 0.85, bands 7 and 8 (1-indexed) approach to a minimum gap of \approx 13.7 and exchange character — an avoided crossing (the gap stays positive; the branches never truly meet). Two independent signals flag the edge:

  1. Projection mass. The 2\times2 sub-block of \Pi over bands 7–8 is no longer diagonal — each band at the left endpoint has substantial M-overlap with both bands at the right endpoint, so neither matched pair is uniquely supported and the relative threshold keeps the off-pair entries.
  2. Eigenfunction overlap. The eigenfunctions themselves rotate into one another across the edge — visible directly below.

Eigenfunctions at the crossing — bands 7 (top row) and 8 (bottom row) at \mu = 0.775 (left) and \mu = 0.85 (right), solved at exactly the overlap pair. The band-7 and band-8 shapes swap their dominant lobe structure between the two \mu-values — the spatial signature of the avoided crossing.

Error table

The run, level by level — the reference Table 1D Error. This is the faithful reproduction: the reference 3 \to 5 \to 7 \to 8 node progression, terminating at L = 3 with every subinterval certified.

Level n_{\text{points}} n_{\text{subintervals}} n_{\text{uncertified}}
0 3 2 2
1 5 4 2
2 7 4 1
3 8 2 0

n_uncertified counts the subintervals flagged for refinement (REFINE) at each level; it falls to 0 at L = 3, where the run terminates with ran_out=False. There is no n_wrong_match column here: that FirstK diagnostic compares each node’s propagated band labels against a fixed-band truth sweep, and window mode’s ragged ranks have no fixed band index to compare against. The CSV behind this table is window/error_table.csv.

Verdict map

Verdict map (Window) — each tested edge as a horizontal segment at y = \text{level}, coloured by majority verdict. Band enter/exit events are recorded for introspection but, under the Algorithm 6 row-norm gate, do not by themselves drive refinement.

Summary. Window mode is the faithful reproduction of the reference §3.3.2 result: the reference 8-node Table 1D Error grid, computed on the reference window I_\lambda = [0, 270], surfacing the band-7/8 avoided crossing its ragged ranks are wide enough to see.

FirstK mode — the fixed-band approximation

The same problem, with different spectral bookkeeping: keep a fixed n_{\text{eigs}} = 6 smallest bands at every node, regardless of where the spectrum sits. This is the rebuild’s generic mode — it needs no problem-specific window — but on this example it is an approximation, not the reference method, and it is honest to say where that bites.

The tolerances are the reference’s stated 1-D values, which also coincide with the rebuild’s library defaults:

Quantity Value Role
Eigenvalues per point n_{\text{eigs}} = 6 fixed band count (rebuild choice)
Pruning threshold t_\pi = 0.21 relative cutoff for off-pair projection mass (reference value)
Eigenvalue cluster threshold t_\lambda = 0.001 groups near-degenerate bands (reference value)

Reference spectrum

The same dense ground-truth sweep as above, here without window shading — FirstK carries six bands, so the figure shows the two low-band avoided crossings that fall inside them (3↔︎4 near \mu \approx 0.45 and \approx 0.86, 5↔︎6 near \approx 0.72), each drawn as a character-tracked X crossing.

Reference spectrum — character-tracked eigenvalue branches from a dense sweep. FirstK with n_{\text{eigs}} = 6 carries the lowest six bands, so the 3↔︎4 and 5↔︎6 avoided crossings show as X crossings; the band-7/8 avoided crossing near \mu \approx 0.81 sits above its fixed cutoff, so FirstK sees only the boundary disturbance it casts onto band 6, never that crossing itself.

Adaptive grid evolution

Drag the slider to step through L = 0..1. The grid grows 3 \to 4 nodes and terminates at L = 1 — it refines less than the window run, for the reason explained in the comparison below.

Adaptive grid per refinement level (FirstK)
Adaptive grid (FirstK) at L=0

Per-level surrogate

Solid coloured curves are the surrogate, reconstructed by the same fixed-slot reshuffle as the window case above; dashed grey is the dense reference. Unlike the window surrogate, FirstK’s is continuous (fixed rank everywhere, no window gaps) but stops at band 6.

Per-level surrogate spectrum (FirstK)
Piecewise-linear surrogate (FirstK) at L=0

A-posteriori projection matrices

Same layout as the Window section above: solid coloured lines are the primary Hungarian match; black dashed lines are alternative matches that survive pruning and trigger refinement.

A-posteriori summary (FirstK) at L=0
A-posteriori summary (FirstK) at L=0

All matrices are square 6 \times 6 throughout — fixed rank on both sides — in contrast with the window run’s ragged, sometimes rectangular \Pi.

Error table

The run, level by level:

Level n_{\text{points}} n_{\text{subintervals}} n_{\text{wrong\_match}} n_{\text{uncertified}}
0 3 2 2 1
1 4 2 2 0

The run terminates at L = 1 with all edges certified (n_{\text{uncertified}} = 0) — on just 4 nodes, coarser than the window run’s 8. It flags one of the two L=0 subintervals, bisects it once, and both halves then certify.

Note

On n_{\text{wrong\_match}}. This column compares the pipeline’s propagated band labels against an independent continuous truth sweep — a strict consistency check available only under FirstK’s fixed bands. It can be positive even when the grid topology is correct, because it counts any single band whose propagated identity differs from the continuous-tracking identity; see the knob-sweep page.

What FirstK is blind to

FirstK’s fixed cutoff at band 6 sits right at the band-7/8 near-degeneracy. It never computes bands 7 or 8, so it cannot see the avoided crossing the window run resolves directly — nor the tighter near-degeneracy among the highest bands near \mu \approx 0.94. What reaches it is only the boundary disturbance that the uncomputed band 7 casts onto its last requested band 6. Under the reference-faithful additive pruning rule, that disturbance does not lift any off-pair projection entry above the matched pair’s by more than t_\pi = 0.21, so the ambiguous last band is certified rather than refined: FirstK flags a single subinterval, bisects once, and terminates at L = 1 on 4 nodes — coarser than the window run’s 8. Lacking the bands that carry the crossing, and with the additive prune declining to manufacture refinement from the boundary ambiguity, the fixed-band approximation simply under-resolves the region the window run refines around. This is the honest cost of approximating a ragged-window problem with a fixed band count.

Comparison: FirstK vs Window

The two modes solve the same PDE with different spectral bookkeeping, and they answer different questions:

  • Window asks “track every band inside [0, 270]” — this benchmark’s actual question, ragged ranks and all.
  • FirstK asks “track the six smallest bands everywhere” — a uniform approximation that needs no window but cannot reach the upper-band physics.

They terminate at different levels and on different grids:

Window FirstK
Bands per node ragged (4–9) fixed (6)
Spectral region I_\lambda = [0, 270] 6 smallest bands
Tolerances t_\pi = 0.21,\ t_\lambda = 0.085 t_\pi = 0.21,\ t_\lambda = 0.001
Terminal level 3 1
Grid nodes 8 4
Reproduces Table 1D Error yes no

The window run reaches 8 nodes and reproduces the reference 1-D grid. FirstK, with the same t_\pi = 0.21 and the reference-faithful additive pruning rule, terminates two levels earlier on just 4 nodes: its fixed 6-band cutoff sits at the band-7/8 near-degeneracy, but the additive prune certifies the resulting boundary ambiguity on band 6 rather than chasing it — so FirstK under-resolves exactly the region the window run refines around.

The two modes do not land on the same grid. Earlier versions of this page showed a different FirstK grid — first a manufactured 8-node match (a tuned t_\pi = 0.35 reverse-engineered to hit the window’s grid, a value found nowhere in the reference), then a 13-node over-refinement — because the pruning rule and preset tolerances had not yet been made reference-faithful. With the faithful additive prune and the stated t_\pi = 0.21, FirstK settles at 4 nodes. The window-mode run is the faithful reproduction; FirstK here is a fixed-band approximation that under-resolves the crossing it cannot see.

Adaptive grid (Window) — 8 nodes, the faithful grid.

Adaptive grid (FirstK) — 4 nodes, the fixed-band approximation terminating early.

The surrogate’s accuracy as a reduced-order model — versus a dense solve at arbitrary \mu — is measured for both modes on the ROM vs dense page.

Side-by-side with the original MATLAB

The panels below place the Python rebuild next to a freshly generated figure from the dissertation’s own MATLAB code (the frozen CodeSeptember dump), so the two implementations can be eyeballed directly. The MATLAB figures are produced headlessly by scripts/matlab_figures.py under MATLAB R2026a and committed as artifacts (CI has no MATLAB, so they are not regenerated on every build).

  • Source: GreedyAlgorithm/CodeSeptember/Main.m, choice1='1', InitialLevel=2, run through scripts/matlab_figures.py (interactive prompts, VideoWriter, and Windows-only paths are neutralised; the frozen dump is never edited — the grid-helpers lev2knots_doubling/tensor_grid are addpath’d from the Clean dump, which only fills gaps CodeSeptember lacks).
  • 1-parameter coefficient (this page): the frozen diffusion solver hard-wires the §3.3.3 2-parameter coefficient [µ₁⁻²; 0.8·µ₂⁻¹; µ₂⁻²] and cannot express the §3.3.2 eq-3.10 coefficient. The MATLAB figures here therefore use a documented, unfrozen adaptationmatlab/matched_eig1D_3_3_2.m, an editable copy of the frozen matched_eig2Dexample.m with only the coefficient replaced by [µ⁻²; 1; 0.7⁻²] (eq-3.10 exactly). Everything else (matching, a-posteriori gate) is byte-for-byte the frozen original.
  • Tolerances: the frozen choice1='1' values t_lambda=0.085, t_pi=0.57, window I_λ=[0,200]not the rebuild’s window preset, so the two grids are expected to differ in node count.

Rebuild — reference spectrum (Window mode), bands descending into the kept window as \mu grows.

MATLAB — matched eigenvalue curves vs \mu from the original code (eq-3.10 coefficient). Same qualitative structure: branches descend as \mu grows, higher bands enter the window near the band-7/8 avoided crossing at \mu\approx0.8.

Rebuild — adaptive grid (Window), 8 nodes, refinement concentrated around the avoided crossing.

MATLAB — adaptive grid from the original code, 6 nodes at the frozen choice1='1' tolerances (t_\lambda=0.085, looser than the window preset). Refinement concentrates in the same [0.7, 0.85] crossing region.

The two adaptive grids are not expected to be node-for-node identical — they run at different tolerances, and the MATLAB grid uses the node-centric forward_points refinement the rebuild only reproduces under its forward_points strategy. What the pairing shows is that both implementations place their refinement in the same crossing region and recover the same descending band structure. Rebuild-only figures on this page (per-level sliders, a-posteriori summaries, verdict map) have no MATLAB counterpart — the original code never draws them.