Method-equivalence validation: the double-well on one engine

This page records a reproducibility check that every sampling strategy and every code path in PyRETIS reproduces the same crossing probability and rate constant on one simple system, and that the result matches an absolute analytical truth rather than only agreeing with itself. The system is a 1D double well sampled on the internal engine; the cases and the driver script live in examples/validation/methods/ and are launched by examples/validation/run_validation.py.

The point of the test is that every sampling strategy, both RNG generations (the legacy MT19937 pin and the shipping PCG64 default), the multi-worker pool, and both rate estimators (the matched crossing-probability report and WHAM) all work over the identical potential, so if their rate estimates agree – and land on the analytical Kramers rate – the sampling machinery is reproducible across strategies, generators and estimators. This complements the cross-engine validation, which checks engine equivalence on a shared force field.

The shared system

Every method case models the same 1D double well, so they must give the same rate within statistical error:

  • Potential: the quartic double well \(V = a x^4 - b x^2 + c\) with a = 1 and b = 2 (so the barrier is \(\Delta V = b^2 / (4 a) = 1\)), c = 0.

  • Dynamics: inertial Langevin at temperature T = 0.07 with a time step dt = 0.025 and friction gamma = 0.3 (reduced units), unit mass. With these numbers \(\beta\,\Delta V = 14.3\), so escape over the barrier is a genuine rare event.

  • Order parameter: the particle position, with eight interfaces spanning \([-0.99 \ldots 1.0]\).

The system parameters are not hard-coded into the analysis – they are read back out of each case’s input TOML (read_system_params), so the analytical reference below tracks the actual runs, and every method case is cross-checked to confirm it really models the same well.

The analytical reference – a known truth

The double well is simple enough to have a closed-form escape rate from Kramers’ theory, so the suite checks the simulations against an absolute truth, not only against each other. analytical_double_well_rate() evaluates it in three increasingly complete forms:

Table 63 Analytical double-well escape rate (reduced units).

Estimator

Rate

Meaning

k_TST

2.81e-07

transition-state theory, no recrossing (upper bound)

k_Kramers (spatial)

2.61e-07

moderate/high-friction Kramers prefactor

k_Kramers-MM (truth)

2.52e-07

  • Mel’nikov–Meshkov turnover factor (headline)

With \(\beta\,\Delta V = 14.3\) and reduced energy loss \(\delta = 11.4\) the turnover factor is \(\Upsilon = 0.97\) – the system sits firmly in the spatial-diffusion regime, so the Kramers result is accurate to a few percent. Converged cases should land on 2.5e-07 within their statistical error. A case can agree with every other case yet sit several sigma from the analytical rate – a shared systematic that the self-consistency check alone cannot see.

The methods

The same double well is sampled with every available strategy and checked through both estimators. Every run goes through the single pyretis run entry point; the cases differ in their RNG pin (*_mt19937 requests the legacy generator explicitly, *_pcg64 pins the shipping default, *_default leaves the config unpinned – byte-identical to its *_pcg64 sibling) and in their analysis kind: matched (the point-matching product of per-ensemble crossing probabilities, reported by pyretis analyse) or wham (the WHAM stitch of the same per-ensemble output, e.g. infswap_wham).

Table 64 Method cases on the shared double well.

Case

Strategy / what it adds

sh_mt19937

standard shooting, pinned to the legacy MT19937 generator

ss_mt19937

stone skipping, legacy MT19937 pin

wt_mt19937

web throwing, legacy MT19937 pin

wf_mt19937

wire fencing, legacy MT19937 pin

wf_ha_mt19937

wire fencing with high acceptance, MT19937 pin

wf_cap_*

wire fencing with a capped window ([tis] interface_cap)

wt_sour_*

web throwing with a shifted source sub-interface ([tis] interface_sour)

relshoot

RETIS with relative per-ensemble shoot frequencies ([retis] relative_shoots)

*_pcg64

the same samplers pinned to the PCG64 generator (rgen = "pcg64") – the A3.4 RNG-migration validation set

infswap_wham

infinite swapping analysed with WHAM – the cross-check of the second estimator

*_default

the same samplers with NO RNG pin (the shipping default): byte-identical twins of their *_pcg64 siblings, a routing/determinism guard

sh_w3

the default config with a multi-worker pool of three workers (PYRETIS_WORKERS=3) – the parallel-worker route

bias_gaussian

biased shooting-point selection (Gaussian selector toward the barrier) – the move biases the shooting point but preserves the path ensemble, so it must reproduce the same rate (bias_gaussian_default is its byte-twin, bias_uniform the uniform-selector identity check)

The wf_cap_* / wt_sour_* / relshoot cases (group params) vary non-default TIS knobs that change the sampling but not the physics, so each must still reach the same analytical rate. The wf_convention group (wf_mt19937/wf_ha_mt19937, infswap_wham, wf_default) cross-checks wire fencing across the estimators. On the scheduler’s per-ensemble-output route the wf occupancy is HA-weighted (the compute_weight crossing count drives the swap, as in the infinite-swapping sampler), so the output writer applies the WHAM Cxy/HA unweighting (Weight = compute_weight / frac) to recover the per-ensemble crossing probability – without it the rate over-counts by ~two orders of magnitude; ss_default carries the same treatment for stone skipping. The suite config’s select key (group tags such as wf_convention or params) runs a chosen subset without dropping the rest.

The matched-kind rows (analysed with the standard PyRETIS crossing-probability report) and the wham-kind rows (analysed with WHAM over the identical per-ensemble output) are a direct cross-check of the two estimators. The *_pcg64 cases are the A3.4 step 2 go/no-go: they run the identical samplers with the PCG64 generator instead of the legacy MT19937 and must reproduce the MT19937 rate within statistical error before the default generator is flipped. The *_default twins must stay byte-identical to their *_pcg64 equivalents – a routing/determinism guard.

Note

ss_default (stone skipping under the worker pool) earlier crashed the internal engine’s streaming dump (FileNotFoundError on ss_shoot.xyz) because the move re-pointed a persistent path’s phase point at the transient dump file. That is fixed (stone skipping now dumps a copy), and the case is enabled like the others.

Running it and reading the output

The suite is not a unit test and not a tutorial, and it is not run per change or per release – see when to run it. Short developer runs can be launched manually. The full multi-seed campaign does not run in GitLab CI: it needs on the order of 36 hours, far beyond the one-hour job timeout, so it belongs on a workstation or cluster node you control, typically driven by campaign_run.sh. CI carries only the method_validation_smoke launch check – a single case at ten cycles – which proves the driver still runs end to end but deliberately produces no validation evidence. All modes use a per-run config such as validation.toml – an ordered list of cases with their target cycle counts plus a per-machine seed and a reverse_list toggle. The first positional argument selects the action:

cd examples/validation

# usage helper (also printed when no action is given)
python run_validation.py

# show the recorded results -- READ-ONLY: runs nothing, writes nothing
python run_validation.py status
python run_validation.py status sh_mt19937  # ... for one case

# (re)analyse the runs already on disk and update the results table
python run_validation.py analyze
python run_validation.py analyze sh_mt19937  # ... for one case

# run every case listed in ./validation.toml, then analyse
python run_validation.py run
python run_validation.py run sh_mt19937      # ... for one case

One-off overrides, usable with run:

python run_validation.py run --config machineB.toml  # different per-machine config
python run_validation.py run --cycles 20000          # override every case's target
python run_validation.py run --seed 2                # this machine's seed
python run_validation.py run --jobs 8                # internal cases in parallel

--jobs N runs that many internal-engine (methods/) cases at once; they are independent single-process runs, so this only uses more cores and does not change any result (each case has its own directory and seed). --cycles is a target total, not an increment – a case with existing output is continued up to that count, otherwise it starts fresh from seed. python run_validation.py --write-config validation.toml writes a fresh config template with every case enabled.

Each case is analysed as soon as it finishes – its rate prints on the go – and a combined inventory, convergence plot and comparison follow at the end. The analysis prints a “Rate vs analytical Kramers reference” table (rate, k/k_ref, |d|/sigma, agree?) next to the suite-mean table. The persistent summaries land in validation_results.json (one entry per case: rate, relative error, cycles, number of independent runs, agreement flag, plus provenance – when, at which git commit, and whether the tree was dirty) and the human-readable validation_results.html (the same table with the convergence rate_vs_cycles.png embedded). A CHECK in the comparison almost always means not converged yet – raise that case’s cycles and re-analyse.

See examples/validation/README.rst for the full driver reference, including the two-machine combine-for-statistics workflow.

When to run it

Seldom, and only for substantial changes. This campaign is not the routine correctness gate and must not be treated as one. Day-to-day confidence that a modification still produces the correct behaviour comes from the ordinary suites:

  • the unit and integration tests (test-easy.sh), which pin the algorithms against hand-computed and synthetic oracles;

  • the example suites (test-heavy.sh), which compare complete runs byte-for-byte against committed reference output on every engine;

  • the tutorial smoke run, which proves every documented workflow still launches.

Those catch a regression far faster, far more cheaply, and far more precisely than a 36-hour statistical campaign can: a golden comparison fails on the first differing digit, whereas a rate estimate has to out-run its own error bars before a change becomes visible at all. A campaign that merely reproduces what the goldens already pin has told you nothing new.

Re-run the campaign when a change can move the sampled distribution itself, so that the goldens are expected to shift and can no longer serve as the reference:

  • a new or modified shooting move, swap move, or acceptance rule;

  • a change to the random-number generator, or to the order in which random numbers are drawn;

  • a change to the scheduler or the coordinator’s ensemble bookkeeping;

  • a new or modified rate / crossing-probability estimator, or a change to the high-acceptance reweighting;

  • a deliberate re-blessing of the reference output, where a campaign is what distinguishes “the numbers moved for the named reason” from “the numbers moved and we no longer know whether they are right”.

Do not re-run it for refactors that keep the goldens byte-identical, nor for I/O, CLI, packaging, documentation, engine plumbing, new examples, or test-infrastructure work. If test-heavy.sh reproduces every committed reference byte-for-byte, the sampled distribution did not move, and the recorded campaign still describes the current code.

How the evidence is recorded

The campaign runs off CI, on whatever machines you have. campaign_run.sh takes the seed, the reverse_list flag and a wall-clock budget as arguments, so two or more machines given different seeds sample independently and combine for statistics, and it continues each case to an ever-increasing cycle target; validation_campaign.toml lists the thirteen cases at 20,000 target cycles. The recorded result is the tracked previous_results/validation_results.json, whose per-case rows carry the rate, relative error, cycle count, agreement flag, and the provenance of the host and commit that produced them. That file is the standing evidence: it stays valid until one of the changes listed above invalidates it, which is why the campaign needs to run rarely rather than per release.

What CI contributes is only method_validation_smoke: one case at ten cycles, run on a throwaway copy of examples/validation so it cannot touch the recorded table. It proves the driver still executes end to end after a refactor and nothing more – ten cycles cannot resolve a rate, and the job publishes no artifacts precisely so its numbers can never be mistaken for evidence.

A row is evidence only when validation_results.json records the analytical-rate agreement and confidence information at a full cycle count; job completion alone is not a scientific pass, and neither is a green method_validation_smoke.