LENR lattice campaign 1 · DFT · September 2026 · all numbers machine-audited (17/17)

H in Ni(111): adsorption sites, strain response, and absorption barriers

A complete first-principles map of where hydrogen sits on and in nickel, what it costs to move it inside, and how lattice strain changes the answer — with the experimental design consequences spelled out.

1 · Scope and method

The question: for hydrogen at a Ni(111) surface, what are the stable sites, the barriers between them, and the strain dependence — i.e. the thermodynamic and kinetic basis for superstitial/coherent-domain loading scenarios. We computed adsorption energies for all four surface sites and two subsurface interstitial sites, the strain derivative of each, and two full minimum-energy paths (fcc→octa, fcc→tetra) with climbing-image NEB. Method: GPAW (PBE, plane waves, 500 eV cutoff, spin-polarized), 4×4 Ni(111) four-layer slab + dipole correction, k-meshes 3×3×1 (relaxation/NEB) and 4×4×1 (final single points, all quoted numbers). Adsorption energies are referenced to ½H₂ in the gas phase. Validation: the site ladder reproduces four independent published studies (Bhatia–Sholl, Noordhoek et al., Shirazi et al., Henkelman-group anchors) within ≤0.02 eV. Every quoted number re-derives from raw output through an independent audit script (17/17 checks).

2 · The structure: slab, sites, and routes

How to read this — Grey spheres: the 4-layer Ni slab (real relaxed coordinates, viewed from above at an angle; drag to rotate, scroll to zoom). Small spheres: hydrogen's energy minimum positions, one per site, colored by adsorption energy (green < 0, red > 0, scale spans ±0.6 eV). Dashed sticks: the two subsurface migration routes rendered from the path construction — orange stays in the fcc column (octa), purple shows the lateral displacement required for tetra. What to notice: all surface H minima sit above the topmost Ni plane; the octa and tetra minima sit below it, and the tetra minimum is laterally displaced from the fcc column — that displacement drives everything in section 5.

3 · Adsorption site ladder

The six local minima for H, with adsorption energies Eads relative to ½H₂(g) at the 4×4×1 mesh:

SiteGeometryEadsInterpretation
fccThreefold hollow, no atom directly beneath (ABC stacking)−0.519 eVGlobal minimum on clean Ni(111); anchor −0.50 ± 0.02 ✓ (Bhatia–Sholl, Noordhoek, Shirazi)
hcpThreefold hollow above a second-layer atom−0.503 eVQuasi-degenerate with fcc (16 meV); standard for (111) surfaces
bridgeTwofold bridge−0.377 eVMetastable ledge; agrees with Shirazi −0.35
topOn-top+0.037 eVBarely bound; agrees with Shirazi +0.05
octaSubsurface octahedral interstice, directly below the fcc hollow+0.178 eVEndothermic vs ½H₂ at zero coverage — the target state for subsurface loading
tetraSubsurface tetrahedral interstice, laterally displaced (see §5)+0.511 eVStrongly endothermic; also laterally offset from the fcc column

k-point bracket: k₄₄₁ vs k₆₆₁ differs by 17.6 meV on fcc → ±9 meV systematic on Eads. Coverage note: the 4×4 cell corresponds to 0.25 ML, which pushes subsurface states up relative to the dilute limit — direction and magnitude consistent with Shirazi's coverage dependence.

Result 1 — thermodynamic exclusion. At zero strain, both subsurface interstitials are endothermic (+0.18, +0.51 eV vs ½H₂) while the surface fcc site is exothermic (−0.52 eV). Clean Ni(111) at low coverage sequesters H at the surface and does not populate the subsurface. This reproduces, from first principles, the empirical requirement of extreme loading conditions before subsurface population appears.

4 · Strain response

How to read this — x-axis: tetragonal strain of the slab's lattice constant, ε ∈ [−2%, +2%] (each point uses the bare-slab reference at the same lattice constant, so the quantity is the strain-induced change of H site stability, not a thermal-expansion artifact). y-axis: Eads vs ½H₂. Slopes are the physics: dEads/dε = −28.0 (fcc), −52.1 (tetra), −49.8 (octa) meV/% — tensile strain stabilizes the subsurface interstitials ~2× faster than the surface site, because enlarging the interstitial volume lowers the H-induced strain energy stored in the surrounding Ni. What to look at: the octa curve crosses Eads = 0 between 0 and +2% — at +2% it is exothermic (−0.030 eV). That crossing is the loading switch.
Result 2 — strain reverses the thermodynamics. dEads/dε for subsurface sites is ≈ −50 meV/%, twice the surface value. A +2% tensile state makes the octa interstice exothermic: the lattice then gains energy by taking H in. Mechanistically, tensile strain enlarges the interstitial cage and relieves the H-Ni repulsion; this is also consistent with Shirazi's coverage trend (high coverage → local lattice expansion → subsurface favorable), reproduced here with an explicit external knob.

5 · Migration barriers: the two routes in

How to read this — Each curve is a climbing-image NEB band (6 images, x = image index: 0 = surface fcc endpoint, 5 = subsurface endpoint) after full relaxation to fmax ≤ 0.12 eV/Å. The maximum relative to the left endpoint is the forward activation barrier; relative to the right endpoint, the reverse barrier. Blue: fcc→octa (attempt-2 protocol, saddle refined by a 4×4×1 single point on the climbing image). Orange: fcc→tetra (attempt-3, converged band at 3×3×1 — see note below). What to look at: the blue peak sits +0.85 eV above fcc and the blue descent lands +0.18 above it — barrier_fwd 0.846, barrier_rev 0.149 eV (both k₄₄₁-refined). The orange curve peaks +1.07 eV above fcc and its right endpoint sits only 34 meV below the peak — the reverse barrier has essentially vanished.

Main road, fcc→octa. Ea = 0.846 eV forward, 0.149 eV reverse (k₄₄₁). Literature anchor range for this barrier is 0.60–0.84 eV (Henkelman-group NEB; Bhatia–Sholl) — we sit at the top edge, expected from 0.25 ML coverage in the cell. Kinetic consequence, quantified: with attempt frequency ν ≈ 10¹²⁻¹³ s⁻¹, the residence time of an untrapped octa H is τ = ν⁻¹exp(Erev/kBT) ≈ 30 ps at 300 K and ≈ 1–2 ps at 650 K. The subsurface state is kinetically evaporating on clean Ni.

Side road, fcc→tetra. Ea ≈ 1.07 eV (3×3×1 converged climbing image; the k₄₄₁ refinement was skipped after an infrastructure failure at exit — bounded by the octa band's k₃₃₁→k₄₄₁ shift of ≲60 meV), reverse barrier 0.034 eV: the tetra interstice is a shelf, not a store. Additionally, a structural finding: the tetrahedral interstice below the fcc hollow is laterally displaced by 1.44 Å (the octa interstice is exactly in-column, displacement 0.000 Å). A naive vertical path therefore drives H through the close-packed Ni layer — our first band literally hit +10 eV Pauli walls at the layer plane — and the MEP must dive in-column, traverse the interlayer gap, and descend into the displaced cage. The tetra route is harder for two independent reasons: a higher barrier and a misaligned door. The purple stick in the 3D view shows the reconstructed path.

Result 3 — kinetic exclusion. Thermodynamics excludes the subsurface (Result 1) and so does kinetics: reverse barriers of 0.15 eV (octa) and 0.03 eV (tetra) mean picosecond-to-nanosecond residence times. On clean, unstrained Ni(111), continuous loading does not accumulate a subsurface population — it only sustains a steady-state leak. Any experimental claim of subsurface accumulation must therefore identify the additional stabilizing ingredient: strain, traps, or coverage-induced modification.

6 · Consequences for NiH experiment design

Each computed quantity maps onto a protocol parameter and an observable. The chain is: number → mechanism → protocol → what you would see.

Protocol problemComputed basisDesign consequenceObservable signature
Getting H in (absorption)Ea = 0.85 eV (octa route): Arrhenius-attempt probability exp(−0.85/kBT) ≈ 10⁻⁸ at 300 K vs ≈ 10⁻⁴ at 700 K; pressure/electrochemical potential shifts the driving force, not the barrierLoad hot — 600–700 K — where thermal activation plus chemical-potential driving can move H over the hill; room-temperature loading of clean Ni cannot populate the subsurface at any practical feed rateUptake rate with ~Arrhenius T-dependence; activation energy ≈ 0.85 eV (or ≈ 1.07 eV if the tetra route dominated — it does not)
Keeping H in (retention)Erev = 0.149 eV → τ ≈ 30 ps (300 K); unassisted retention is zero on experimental timescalesContinuous H₂ feed (needle-valve protocols) is mandatory and exactly compensates the leak; steady-state subsurface coverage is set by feed rate × residence time. For true accumulation you need the stabilization (strain, below) or traps (vacancies/voids — binding energies 0.4–0.6 eV in the literature, to be computed next)Post-loading desorption transients on a ps-leak system look like instant loss at T; a trapped population would instead show staged, activated release (trapping energies set the release temperatures)
Making inside favorable (the switch)dEads/dε ≈ −50 meV/% for subsurface vs −28 meV/% surface; exothermicity flip at ε ≈ +2%. Mechanism: interstitial volume enlargement. Practical strain sources: external tensile strain, or internal — the FM→PM magnetoelastic expansion at the Curie point, tunable via alloying (NiCu, NiFe), plus thermal expansion itselfUse tensile-strained foils/foams or a Curie-matched alloy so that the magnetic transition delivers the expansion at the operating temperature; expect the sign of uptake enthalpy to flip as the material crosses the thresholdA sharp uptake knee / exothermic burst when strain state crosses ≈ +2% (or the alloy equivalent); absence in the unstrained control
Coupled feedbackH lowers the Ni moment and (literature) the Curie temperature; the Curie transition changes the lattice; loading therefore feeds back on the switchTreat loading fraction, strain state, and magnetic state as one coupled system during analysis — not as independent knobsHysteresis between loading and de-loading; uptake features tracking TC(loading) rather than fixed TC(0)
The decisive experiment this campaign points at: a Curie-matched Ni alloy (TC tuned to the intended operating temperature by Cu/Fe/Cr fraction) loaded at 600–700 K under continuous H₂, against an unstrained pure-Ni control. Prediction from these results: the alloy shows favorable (exothermic) subsurface uptake once the magnetoelastic expansion crosses the ≈ +2%-equivalent threshold; the control never does; the crossover tracks TC(composition), not temperature alone. Cheap, falsifiable, and it tests the strain/Curie switch picture directly. The tetra shelf is predicted to accumulate nothing on clean Ni — its absence in any dataset is a consistency check, not a surprise.

7 · Next computational rungs (proposed)

(1) Paramagnetic/disordered-local-moment site ladder at experimental lattice constants vs T — quantifies whether the Curie transition alone delivers the ≈+2%-equivalent stabilization. (2) Vacancy and void trapping energies for H — the quantitative "unless trapped" term with release temperatures. (3) H-induced Curie shift and the loading↔magnetic feedback loop. (4) QE cross-validation of the two barrier numbers with the independent code.

8 · Tooling

ToolTypeRole hereStatus
GPAW 26.7Plane-wave/real-space DFTAll campaign numbers (PBE-PW, spin-polarized, dipole-corrected slabs)Validated — 4 literature anchors reproduced
ASE 3.29Atomistics frameworkNEB/climbing-image machinery, FIRE optimization, structure I/O, image-parallel orchestrationValidated — campaign infrastructure
Quantum ESPRESSO + SSSP pseudosIndependent plane-wave DFTCross-validation of barriers/site energies with a second code — guards against implementation artifactsInstalled, staged for campaign 2
LAMMPSClassical MDμs-scale diffusion and many-H loading configurations once a Ni–H potential is fitted to this campaign's DFT dataInstalled, awaiting fitted potential
SymPy / SciPySymbolic & numeric oraclesIndependent re-derivation of every quoted claim (mutation audit)Active — 17/17
Simulations: GPAW PBE-PW 500 eV, spin-polarized; 4×4 Ni(111), 4 layers; CI-NEB. Site ladder and strain ladder at k₄₄₁; barriers as annotated. Technical report: RESULT.md · audit script and raw JSONs in the same repository · commit 403d01f+ · generated 2026-09-22.