Vantasner Coherence Model
[CODE]
01 SUMMARY 02 MODEL 03 DETAIL 04 SOURCE

[01 / THE METHOD]

Start rough, then refine.

Build the cheapest version of the whole model, then spend compute only where the answer is most sensitive. A small team of AI agents made six models in a chain over four passes.

  • [PASSES]0→3from a few equations to electronic-structure calculations
  • [MODELS]6in a chain where each one's output is the next one's input
  • [TEAM]1+2a lead agent and two specialists
  • [LARGEST RUN]106lattice sites in one simulation
Three versions of the same small clevis bracket on a dark workbench, left to right: a rough clay model, a bent wire frame, and a machined metal part with bright edges. [Studio still life · one bracket in clay, wire and steel]
A studio still life used as a metaphor. The same part is made three times, and each version is sharper than the last. The models are built the same way.

[IN SHORT]

  1. 01

    Build the cheapest version of the whole model first so that every link of the chain exists, even if it is rough.

  2. 02

    Rank the unknowns by how far they move the answer and spend compute on the top of the list.

  3. 03

    Each pass swaps a guess for a calculation. The passes move from equations to fitted laws, to the electrons of the metal and finally to a million-site lattice.

  4. 04

    A lead agent and two specialist agents did the work. Every result is a file in a public repository, next to the script that made it.

[01.1 · COARSE TO FINE]

Each pass is sharper than the last

Pass 0 builds the whole chain from a few equations. Each later pass replaces the guesses the output depends on most with something computed. The passes are sized by compute rather than by calendar time. The last page sets out what the next pass would compute.

[FIG 01] The same model at four resolutions

  1. 0 [PASS 0 · NEGLIGIBLE]

    A few equations

    This pass builds the whole chain in its cheapest form. It solves protons coupled to the field exactly, then switches five candidate ingredients on and off using material scales that were looked up.

    [MODELS]01 Protons and the field02 Ingredients

    [COMPUTE]NumPy only

    g0_stabilityg0_modelg0_scales

  2. 1 [PASS 1 · LIGHT]

    Sweeps, laws and a ranking

    Sweep every unknown across its plausible range, write what comes out as closed-form laws, and rank the unknowns by how far they move the output.

    [MODELS]02 Ingredients

    [COMPUTE]NumPy only

    g1_rank

  3. 2 [PASS 2 · MODERATE]

    Compute the top of the ranking

    The top-ranked inputs are how stiff each proton well is and how it bends. Calculate them from the electrons of the metal, then refit the model with them.

    [MODELS]03 Computed wells04 Anharmonic oscillator

    [COMPUTE]164 core-minutes of GPAW

    g2_inputg2_refit

  4. 3 [PASS 3 · HEAVIER]

    Scale, layers and loading

    This pass adds a lattice of up to a million sites, nickel layers on palladium, strain and shear laws, and how the well changes with loading.

    [MODELS]05 A million sites06 Nickel on palladium

    [COMPUTE]About 2,870 core-minutes of GPAW for loading (the most expensive step) and 200 for the strain laws

    g3_latticeg3_skewg3_loading

Notice that the same bracket appears at every step and is drawn more sharply as the compute rises. Block height ranks the compute tier that the method names for each pass (negligible, light, moderate, heavier). Table C under The detail lists the costs the repository states.

[FIG 02] The ranking decides where the compute goes

Which unknowns move the answer most?

In pass 1 every unknown was swept across its plausible range and scored by how many decades (factors of ten) it moved the model's outputs. Screening and the shape of the well came out on top, and both can be computed from the electrons of the metal. Pass 2 computed them.

  • Electronic structure can compute it
  • Found by lookup
  • Neither
  1. Screening lengthHow strongly the metal electrons screen the pull between protons 5.0 Computed: force constants from GPAW · Model 3
  2. Quartic coefficientHow far the well departs from a parabola 2.0 Computed: quartic coefficients from GPAW · Model 3
  3. LossHow fast a proton vibration decays 2.0 Open: no measured decay rate found · Model 5
  4. Choice of vibration peakWhich measured vibration line is used 1.6 Replaced: levels of the computed well · Model 4
  5. Stiffness (gradient term)How strongly a displacement pulls on its neighbours 1.0 Replaced: computed neighbour couplings · Model 5
  6. LoadingHydrogen per metal atom 0.74 Computed: the loading law · g3_loading
  7. Field mode choiceWhich electromagnetic mode is coupled 0.28
  8. Electron countConduction electrons per site 0.10
  9. Lattice constantSpacing of the metal atoms 0.017
  10. Electron responseStrength of the electron-plasma coupling 0.0015
  11. Temperature 0

Decades the outputs move across the input's range (log10 of largest over smallest)

Notice how far screening stands out, with five decades against two for the next. Everything below the top six moves the outputs by less than 0.3 decades. Each bar is the largest move across seven outputs and three materials (NiH, TiH2, PdD). g1_rank/rank.json

[01.2 · THE CHAIN]

Six models in one chain

Each model takes what the one before it found and adds the next piece of physics. Two later calculations feed back into the wells. Select a model to see what goes in, what comes out and what it cost.

[FIG 03] The chain of models coloured by pass

PASS 0 PASS 0–1 PASS 2 PASS 3
[INPUT]Published theory
[INPUT]Measured scales
[FEEDS BACK] revised wells MODEL 01 Protons and the field Negligible MODEL 02 Ingredients Light MODEL 03 Computed wells Moderate MODEL 04 Anharmonic oscillator Moderate MODEL 05 A million sites Heavier MODEL 06 Nickel on palladium Heavier FEEDS MODEL 5 Loading law FEEDS MODEL 6 Strain and shear laws

[MODEL 01 · PASS 0 · NEGLIGIBLE]

Protons and the field

Compute tier: negligible

[IN]
The published coupling of protons to the light field (Gamberale and Modanese, arXiv:2205.14015), with every counter-rotating term kept.
[OUT]
The normal modes of protons and field together. When the same field enters the coupling and the diamagnetic term, the field alone gives no ordered state at any coupling strength. Another ingredient has to make the protons act together.
[COMPUTE]
NumPy. Closed forms checked by direct diagonalisation.

Model 01 in full g0_stability

[MODEL 02 · PASS 0–1 · NEGLIGIBLE TO LIGHT]

What could make them order

Compute tier: negligible to light

[IN]
Five candidate ingredients treated as switches and tested in all 32 combinations for NiH, TiH1.84 and PdD, using measured material scales that were looked up.
[OUT]
Two routes to order: a static one (attraction between protons plus an anharmonic well) and a pumped one. A ranking of every unknown by how far it moves the output, with screening at the top.
[COMPUTE]
NumPy. Sweeps and closed-form laws.

Model 02 in full g0_modelg0_scalesg1_rank

[MODEL 03 · PASS 2 · MODERATE]

The computed proton wells

Compute tier: moderate

[IN]
Electronic-structure calculations of NiH and PdD at full loading (GPAW, PBE), with each hydrogen nudged by 0.01 Å either way.
[OUT]
How stiff each well is, how far its shape departs from a parabola, and how straining the metal changes it.
[COMPUTE]
164 core-minutes of GPAW, including failed and stopped attempts.

Model 03 in full g2_input

[MODEL 04 · PASS 2 · MODERATE]

The anharmonic oscillator

Compute tier: moderate

[IN]
One hydrogen in its computed well, coupled to its neighbours by the computed force constants. It was run first on the coarser wells and then again on the corrected ones.
[OUT]
Its vibration levels, the stretch at which the static route orders, and the frequencies the pumped route would use. On the corrected wells the mode in which every hydrogen moves in step comes out at 85.6 meV for NiH, 62.3 for PdH and 39.8 for PdD, against measured peaks of 89, 55.8 and 37. The computed values are 3.8 percent below, 11.6 percent above and 7.6 percent above the measured ones. The static route orders at 8.1, 7.5 and 5.3 percent stretch at 0 K. These values are extrapolated from calculations at ±1 percent.
[COMPUTE]
NumPy. The levels come from direct diagonalisation in three dimensions.

Model 04 in full g2_refitg3_lattice/rerun

[MODEL 05 · PASS 3 · HEAVIER]

A million sites

Compute tier: heavier

[IN]
The computed wells and couplings on a lattice of up to 1,000,000 sites with loading, vacancy arrangement, strain, temperature and pump as settings. It was rerun on the corrected inputs.
[OUT]
Spatial patterns, and under the pump an oscillation at half the frequency at which the pump modulates the well. When pumped at 0 K on the corrected inputs, NiH moves in step across the whole 110,592-site box, with 99.8 percent of the motion in the uniform pattern. PdD stays in patches that move in step over about 1 to 2 nm.
[COMPUTE]
NumPy. The first main batch takes about a minute; a million sites for 32 steps takes 12 seconds. On one thread, the 110,592-site NiH rerun took under seven minutes and the 30 PdD cases about an hour.

Model 05 in full g3_latticeg3_lattice/rerun

[MODEL 06 · PASS 3 · HEAVIER]

Nickel on palladium

Compute tier: heavier

[IN]
Four nickel layers on six palladium layers, relaxed with the EMT potential in ASE. Ni at 3.52 Å on Pd at 3.89 Å gives about 9.5 percent mismatch. The stack is relaxed again with the Pd at the PdD spacing of 4.084 Å.
[OUT]
A repeating strain pattern, which the strain law turns into a map of the deuterium wells. At the PdD spacing the Pd layer at the interface splits into soft lines with stiffer islands between them. There is one line in each direction per repeat. When the map is run on the lattice, pumped motion gathers in the soft regions but does not repeat with the pattern. Every site lies outside the computed strain range, so the map is an extrapolation.
[COMPUTE]
About two CPU-minutes.

Model 06 in full g3_skewg3_lattice/rerun

[FEEDS BACK INTO MODEL 5 · PASS 3 · HEAVIER]

Loading law

Compute tier: heavier

[IN]
A 64-atom PdD cell at 500 eV with 3×3×3 k-points. One deuterium moves while the others and the 32 palladium atoms stay in place. The cell has every site filled, or 4 or 8 of the 32 sites empty, either spread out or bunched. The same settings were rechecked for NiH.
[OUT]
The well and the couplings to neighbours at full loading, replacing the coarser values behind Model 4. A vacancy softens the deuterium beside it and pushes it toward the vacancy, and the arrangement of the vacancies matters. A two-term law in lattice spacing and empty neighbours carries the trend with loading for the lattice.
[COMPUTE]
About 2,510 core-minutes of GPAW for the PdD loading set, including the full-loading well, and about 2,870 with the NiH check. This is the most expensive step so far.

g3_loading

[FEEDS BACK INTO MODEL 6 · PASS 3 · HEAVIER]

Strain and shear laws

Compute tier: heavier

[IN]
The two-atom PdD cell computed with GPAW, both sheared and squeezed or stretched as a thin film.
[OUT]
A stiffness law for any strain inside the computed range. The strains at the interface still lie outside that range, so the well map there remains an extrapolation.
[COMPUTE]
200 core-minutes of GPAW.

g3_skew/partA

[01.3 · WHO DID THE WORK]

A lead and two specialists

A small team of AI agents did the work. The lead plans the passes and checks every result. One specialist runs the electronic-structure calculations and the other builds the physics models. They coordinate through a shared project graph in which each piece of work is a node with a written goal.

[FIG 04] How one piece of work moves from goal to record

  1. 01 · GOALA node with a goal

    The lead writes what to compute, why, and what counts as done.

  2. 02 · WORKA specialist runs it

    The work is an electronic-structure calculation or a physics model, along with the script that makes it.

  3. 03 · DELIVERYDelivered as files

    Each delivery has a script, its output and a README. Calculations add a cost ledger that charges failed runs too, and SHA-256 checksums that tie each number to the file it came from.

  4. 04 · CHECKThe lead verifies

    Readout scripts rebuild the tables from the saved results, so a check does not need the costly run again. Anything wrong goes back as a correction.

  5. 05 · RECORDCommitted in public

    The result enters the coherencev1 repository, and the next pass starts from it.

Corrections travel the same path as results, so they stay in the record. Four of them are below.

Corrections are part of the record

[NUMERICAL SETTINGS] Found by the physics model and fixed by a new calculation

Coarse settings were not enough

Model 4 put the PdD vibration about 45 percent above the measured value. It pointed to the 64-atom cell behind its well, which had been computed with a single k-point (one sample of the electron states). The cause was the coarser numerical settings taken together. The cell was recomputed with a 500 eV cutoff instead of 400, 3×3×3 k-points and a tight force-convergence criterion.

PdD stiffness for all deuterium moving together: 1.03 eV per square ångström at the coarse settings, 0.202 at the tighter settings, 0.181 from the converged two-atom calculation. 64 atoms, 400 eV, 1 k-point 1.03 64 atoms, 500 eV, 3×3×3 k-points 0.202 2 atoms, k = 12 (converged) 0.181 00.40.81.2 eV/Ų, PdD, all D moving together

Notice the recomputed cell lands next to the converged reference.

It happened a second time. An agent building these report pages noticed that the NiH 64-atom well had used the same coarse settings. When rechecked at the tighter settings, the NiH neighbour coupling came out about three times larger (0.23 to 0.68 eV/Ų). The collective value then agreed with the two-atom cell (1.036 against 1.050).

g2_refitg3_loading

[TWO ROUTES] Diagnosed before the well shape was computed

Energy and force must agree

The stiffness of a well can be read from its energy curve or from its restoring force. With a loose stopping rule the two routes disagreed. The run was continued under a stricter rule, and the two routes now agree to 1.3 percent.

On-site PdD stiffness. Loose run: 0.598 from energy, 0.488 from force. Tightened run: 0.599 from energy, 0.591 from force, in eV per square ångström. LOOSE 0.488 0.598 TIGHT 0.591 0.599 0.450.500.550.600.65 force energy on-site stiffness, eV/Ų

Notice the gap closing from about a fifth to about a hundredth.

g3_loading/x1

[BLIND CHECK] A point the law never saw

The strain law predicts an unseen case

When the spacing is left free, the calculation prefers a PdD spacing of 4.134 Šrather than the 4.084 Šused as input. The strain law was fitted without that point. Its prediction of the stiffness there misses by 0.003 eV/Ų.

Collective PdD stiffness at 4.134 ångström: predicted by the strain law -0.058, computed directly -0.055 eV per square ångström. law −0.058 computed −0.055 −0.07−0.06−0.05−0.04 collective stiffness at 4.134 Å, eV/Ų

Notice how close the two marks sit on this scale.

g3_skew/partA

[UNCONVERGED INPUT] A striking result that did not survive

Stripes from a coupling three times too small

The first lattice simulation ran NiH with the neighbour coupling from the coarse 64-atom cell. That coupling was 0.23 eV/Ų, where the recomputed value is 0.68. When pumped at 0 K, the lattice broke into striped patches about 2 nm across. A million sites gave patches no larger than 110,592 did (2.16 against 2.01 nm). After the coupling was recomputed, the same 110,592-site box moved in step from edge to edge with 99.8 percent of the motion in the uniform pattern, and the stripes were gone. The rerun updated the NiH inputs together, so it does not single out the coupling.

The earlier result came from an unconverged input. It was caught because an agent building these report pages questioned the settings behind the NiH well (the first card above). The lattice page keeps the striped run as a record.

NiH in-step length at 110,592 sites, 0 K: 2.01 nanometres with the coarse coupling, 3.85 with the recomputed one, which is the size of the box. Coupling 0.23 eV/Ų, stripes 2.01 Coupling 0.68 eV/Ų, in step 3.85 whole box 01234 in-step length, nm, NiH, 110,592 sites, 0 K

Notice that the recomputed run reaches the dashed line, which is the most its box can show.

g3_lattice/coherenceg3_lattice/rerun/part2g3_loading

[01.4 · THE TOOLS]

Three tools

The work uses open-source scientific Python throughout. Each script, its output and its README sit together in one directory.

[FIG 05] Where each tool is used

Tool and what it does010203040506LLSL
NumPy and SciPyThe models themselves: normal modes, sweeps and laws, the oscillator's levels, the million-site lattice, and every readout.
GPAWDensity functional theory. It computes the electrons of the metal and from them the force on each hydrogen, which gives the stiffness and shape of its well.
ASE with EMTPlaces the metal atoms at the nickel-palladium interface using a fast built-in potential (effective-medium theory). ASE also builds the crystal cells GPAW runs on.
Columns 01 to 06 are the six models. LL is the loading law and SL the strain and shear laws. Notice that the costly GPAW step runs only where the ranking pointed.

[01.5 · THE DETAIL]

The detail

Three rules carry the method: how an unknown is ranked, how a stiffness is computed, and when a calculation counts as converged.

[1] Sensitivity of an unknown

Di = maxy, material log10 [ maxxi ∈ Ri |y| / minxi ∈ Ri |y| ]

xi is one unknown. It is scanned at 41 points across its plausible range Ri with every other input at its baseline. y runs over seven outputs: threshold loading, threshold density, mode frequency, static displacement, ordering energy, correlation length and ring-down time. The maximum is also taken over NiH, TiH2 and PdD. An output that is zero somewhere in the scan (at a phase boundary) has no finite ratio and is left out. Di is the bar length in Fig. 02.

[2] Two routes to the stiffness of a well

KF = − [ F(+h) − F(−h) ] / 2h     KE = [ E(+h) − 2E(0) + E(−h) ] / h2

h is the displacement of the hydrogen from the centre of its site. It is 0.01 Šin the two-atom cells and 0.03 Šin the 64-atom cell. F is the force on the displaced hydrogen along the displacement and E is the total energy. Both come from GPAW with the metal atoms held fixed. K is in eV/Ų. In the 64-atom cell the two routes must agree within 5 percent before the well is accepted.

[3] When a k-point mesh is accepted

| Kn − Kn′ | / | Kn′ | ≤ 0.05

Kn is the stiffness for all hydrogen moving together computed on an n×n×n mesh of k-points (sample points for the electron states). n′ is the next denser mesh. The accepted mesh is the coarsest one from which every following pair passes.

Table A. The ranking behind Fig. 02

UnknownWho can supply itDecades, DiCrosses the ordering thresholdWhat the later passes did
Screening length Electronic structure 5.0 yes Computed: force constants from GPAW
Quartic coefficient Electronic structure 2.0 no Computed: quartic coefficients from GPAW
Loss Lookup 2.0 no Open: no measured decay rate found
Choice of vibration peak Lookup 1.6 yes Replaced: levels of the computed well
Stiffness (gradient term) Electronic structure 1.0 no Replaced: computed neighbour couplings
Loading Lookup 0.74 yes Computed: the loading law
Field mode choice Neither 0.28 no
Electron count Electronic structure 0.10 no
Lattice constant Lookup 0.017 no
Electron response Electronic structure 0.0015 no
Temperature Lookup 0 no

Table B. Stiffness against k-point mesh

Mesh nNiH, K (eV/Ų)PdD, K (eV/Ų)
6 1.174 0.421
8 1.050 accepted 0.221
12 1.028 0.181 accepted
16 1.017 0.175

All hydrogen moving together in the two-atom cell, from GPAW PBE with a 500 eV plane-wave cutoff. g2_input/convergence.json

Stiffness against k-point mesh for NiH and PdD; the values are in Table B. 0 0.4 0.8 1.2 681216 k-point mesh n (n×n×n) K, eV/Ų NiH PdD
Notice that PdD needs a far denser mesh than NiH before it settles. Its value on the 6×6×6 mesh is more than twice the converged one. Large dots mark the accepted mesh.

Table C. Compute as the repository states it

DirectoryWorkStated computeNote
g2_inputGPAW wells, NiH and PdD164 core-minutesNiH 80, PdD 85, including failed and stopped attempts
g3_skew/partAGPAW strain and shear laws, PdD200 core-minutes
g3_loadingGPAW 64-atom wells: PdD at full and reduced loading, NiH checkabout 2,870 core-minutesabout 2,510 of it for the PdD loading set including the x = 1 well; the most expensive step
g3_skewInterface geometry, ASE with EMTabout two CPU-minutes
g3_latticeLattice simulation, main batchabout a minutea million sites for 32 steps takes 12 seconds
g3_lattice/rerunLattice rerun on the corrected inputs: 30 PdD cases and one 110,592-site NiH caseabout an hour (PdD) and under seven minutes (NiH)one thread; part2/updated/execution.json
g0_*, g1_rank, g2_refitField model, ingredients, ranking, oscillatorno figure statedNumPy and SciPy

[WHAT THIS MEANS FOR THE NEXT MODEL]

The chain starts at its cheapest link

Model 1 couples protons to the light field and solves the result exactly with NumPy alone. It asks whether the field alone can make the protons act together. Its answer sets what every later model has to add.