MODELLING & SIMULATION
Five Cases, With Their Limits
Coupled thermal, hydraulic and mechanical simulation using PFLOTRAN and OpenGeoSys, on source terms taken from published data or derived with open-source depletion. Each case states what was asked, how it was modelled, what came out, and what the result does not support. Where a comparison failed, it is listed as a failure.
01 — THERMAL ARRAYS
Thermal Interaction in Canister Arrays
Problem. A single canister checked in isolation will pass a temperature limit that the same canister fails once its neighbours are emplaced. The layout question — how close, in which direction — is an array question, and the answer changes with the waste form.
Method. A finite-height axisymmetric coupled thermal–hydraulic model of one package in granite with a bentonite buffer, solved in PFLOTRAN; array response assembled by superposition of that response over the emplacement lattice. Three waste forms were compared on one identical 6 m × 40 m model with only the decay-heat source changed between them. A separate fast screen swept tunnel spacing with canister spacing held at 6.0 m.
Finding. Peak inner-buffer temperature was 69.408 °C for PWR spent fuel, 23.518 °C for a GTHTR300-12 form and 18.947 °C for HTR-PM pebbles — same repository geometry, a 50 °C spread. Across tunnel spacing, the PWR case passed a 70 °C limit at 35 m (69.617 °C) and marginally failed at 30 m (70.113 °C), then rose monotonically to 106.409 °C as spacing closed to 4 m. Both gas-reactor forms stayed below 30 °C even at the tightest spacing screened.
Limit. Screening, verified and not validated — no measured data enters these runs. The 4 m spacing is the tightest geometry retained for sensitivity, not a proposed design, and no minimum pitch is nominated. The pebble decay-heat curve is sourced only to year 20 and extrapolated beyond it on a conceptual package, so the near-source temperature is not credible in absolute terms even where the buffer result is informative.
02 — ENGINEERED BARRIERS
Buffer Resaturation in a Finnish Deposition Hole
Problem. After closure the bentonite buffer takes up water from the host rock. How long that takes, and which uncertain input controls it, determines when the barrier reaches its design state — and which parameter is worth measuring better.
Method. A single KBS-3V deposition hole on an Olkiluoto basis, run in PFLOTRAN: five 200-year single-phase cases spanning two sensitivity families — host-rock permeability across a factor of 43, buffer permeability across a factor of 7.24 — plus one 20-year two-phase comparison run.
Finding. Peak temperature landed between 63.1 and 64.9 °C at 0.8–1.6 years in every combination tested: the near field is conduction-dominated and close to insensitive to permeability. Saturation timing is the opposite. Time to reach S ≥ 0.999 in the inner buffer moved 1.72× across the rock bracket but 4.13× across the buffer bracket, which is six times narrower. The buffer is the bottleneck; the host rock never was, staying 27–194× more permeable throughout.
Limit. No result here has been compared with measured data. The elapsed times are saturation-threshold times against a stated numerical criterion, not a repository resaturation duration — no sourced transition criterion exists for that. The two-phase comparison has been run once and is directional only, pending a numerics pass.
03 — COMPARISON WITH MEASURED DATA
FEBEX In-Situ Heater Test — Cross-Code Comparison
Problem. Analytical benchmarks confirm that a code solves its equations. Only a real experiment tests whether the model represents the system. The FEBEX in-situ heater test, with its published multi-code benchmark and instrument record, is the hardest such check available for a bentonite near field.
Method. An OpenGeoSys thermo-hydro-mechanical model of the test, compared against the published measured dataset and the benchmark code ensemble. Two corrections earned the result: the heater was driven by source-backed electrical power instead of a prescribed liner temperature — the figure previously used as a boundary condition turned out to be a maximum at a hot spot, not a uniform surface — and the domain was extended radially and axially until it stopped truncating the heat being injected.
Finding. Temperature: +2.37 °C bias, 2.65 °C RMSE and 4.49 °C maximum deviation across five outer-ring thermocouples, with 136 of 136 comparison points inside ±5 °C and 60 of 136 inside the ±2.2 °C band taken from the source. Moisture: a failure, and recorded as one. Near-heater drying is reproduced to within about 2 percentage points of the measured degree of saturation, but resaturation from the rock side runs roughly 20 points low at 4.93 years — the model wets the outer buffer too slowly, and the tested formulation does not close that gap.
Limit. The temperature agreement is frozen to the first 400 days and claims nothing beyond it: the thermal diffusion length reaches roughly 59 m by year 18, so boundaries that are adequate at 400 days would be truncated at the full test duration, and any longer-time claim needs a fresh domain-convergence study first. The moisture leg remains open — a disclosed model limitation, reported rather than tuned away.
04 — VERIFICATION
Numerical Verification and Source Audits
Problem. A coupled result inherits errors from three places that never appear in the output: the deck, the coupling scheme, and the publication the parameters were copied from. A plot that looks physical is not evidence that any of the three is sound.
Method. A published PFLOTRAN thermo-hydro-mechanical verification suite was reproduced file-for-file — 170 of 170 verified — with the solver commit and PETSc build pinned, then documented as a single coupling contract mapping every governing equation to its input card, parameter value, units and sign convention. Source terms are separately audited back to the primary publication before they are allowed into a deck.
Finding. The fixed-stress sequential coupling sequence is documented and qualified against its parent case. Two audit results changed what went into later models. A Biot coefficient of 1.0 across all six benchmark cases removes a term from the porosity update entirely — no case in the suite exercises that branch, so a study relying on it inherits an untested path. And a published 1700 W canister load reconciles with its own published 96.15 W·m⁻² surface flux only over the full closed cylinder, 17.681 m²; applying that power to the lateral wall alone runs about 10 % high on near-field peak temperature, and the source never states which surface it used.
Limit. Verification, not validation. These comparisons are against analytical solutions, a published benchmark suite and source documents — no field measurement is involved, and passing them says nothing on its own about whether a model represents a real site.
05 — INDEPENDENT RECONSTRUCTION
Reconstructing a Published Thermal Result
Problem. A published thermal analysis of a horizontal-emplacement used-fuel concept reports two temperature peaks — an early one near 86 °C and a second, in the same 75–85 °C band, at 900–1000 years. A late second peak changes what a design has to tolerate, so before reusing the result we asked which length scale in the repository could deliver heat that late.
Method. One axisymmetric finite-length conduction kernel in PFLOTRAN, then linear superposition over four nested geometries — the isolated container, a 1.5 m-pitch row, two stacked tiers aligned and staggered, and the full 20 m-pitch room lattice out to the 1400 m placement extent, 29,110 sources at the largest case. The heat curve was digitized from the publication and gated to 10-16 relative against the delivered volumetric source. The second-peak test was declared before the runs: ordered maximum, minimum, maximum, prominence at least 0.5 °C, later peak at least ten times the first peak’s time.
Finding. One peak, then monotone decline, in every case. Magnitudes converge on the published values — 80.8 °C at the room-lattice scale, 84.1 °C on the conservative bracket, against a reported 86 °C — so the reconstruction reaches the right temperatures and not the reported shape. Repository extent is not the missing ingredient: widening the placement area from 1400 m to 2600 m moves the 950-year temperature by 0.10 °C, and the tier stagger by under 0.01 °C. Decomposition shows why: the room lattice supplies 41 of the 50 °C present at 950 years but peaks around 150 years, because the nearest rooms outweigh the distant ones that would arrive late.
Limit. This is not a finding that the publication is wrong, and it is not a validation of anything. It states what one model class produces: under a monotonic digitized heat curve, a finite container, fixed fully saturated thermal properties and linear superposition, the second peak is not a container, row, array or repository-extent effect. Four explanations remain open and untested — how the total heat rate is applied to the container volume in the original model, how the published curve is interpolated between tabulated ages, how the symmetry fraction and per-container power are scaled together, and property evolution during resaturation, where the buffer conductivity in the source’s own fits is 21 % lower while unsaturated. Any of those could account for the difference. The work is parked pending that source audit, not concluded.
Software Example — Decay Heat × Spacing Interface
When a study gets run repeatedly, the workflow around it earns an interface. One internal browser tool takes a decay-heat curve and a spacing pair and launches a direct PFLOTRAN solve on a fixed 3 × 3 array — canister spacing 3–60 m, tunnel spacing 6–120 m, about 47 minutes per run — alongside bounded look-ups that reuse frozen result kernels in milliseconds for cases already qualified.
It is an interface for running the model, not a qualification of the layouts it will accept. Every combination it lets you enter is an exploratory run; any layout intended to support a decision still needs its own verified case, built and checked the way the five above were. Internal tool — a walkthrough is available on request.
Source Terms — Decay Heat Without a Proprietary Bottleneck
Every thermal case above starts with a decay-heat curve, and that curve decides the answer. There are two defensible ways to obtain one, and the choice is made per study rather than by habit. The first is a documented published curve: the array work in case 01 uses a tabulated specific decay heat in watts per tonne of heavy metal on a years-after-discharge clock, scaled from 45 to 47 GWd/t burnup, applied to a declared 1.844 tHM package loading at a declared 30-year cooling age — 1,968 W per package at emplacement, 408 W at 200 years. The second is a case-specific curve derived in the open: a full transport–depletion calculation in OpenMC 0.15.3 on ENDF/B-VIII.0 neutron data with the matching thermal depletion chain, normalized by fission Q value, burned at the reactor’s own specific power to its equilibrium discharge burnup and then decayed to the emplacement cooling age. For the gas-reactor fuel block in case 01 that route gives 25.1 W falling to 4.4 W per reference block between 30 and 230 years of cooling.
The point of stating this is practical: a licensed depletion code is not a precondition for every thermal screening study. Where a published curve exists for the fuel, burnup and cooling age in question, it is the cheaper and more traceable input. Where no published curve covers the case — a fuel form still at concept stage, or a cooling age outside a curve’s stated validity window — an open-source depletion route makes the study possible at all, provided the inputs and checks support it. Either way the assumptions are written down: fuel type and enrichment, burnup and its denominator, cooling age and its time origin, package loading, normalization basis and nuclear data evaluation. Two of the cases above changed materially once those were audited rather than assumed.
Limit. An open-source depletion calculation is not automatically equivalent to a licensed one — different nuclear data, chain, normalization and statistics give different answers, and no equivalence has been demonstrated here. A derived curve is not a validated curve: the statistics behind the depletion work above are formally unqualified, that benchmark line is parked pending source-convergence qualification and a missing reference dataset, and one earlier synthetic curve was retired precisely because it was never a depletion calculation. Published curves carry their own limits — validity windows are stated and enforced, which is why one case’s permissible horizon is 20 years rather than 200. Curve uncertainty is reported separately from package-loading uncertainty; the two are not combined into a single number.