1. Thermodynamic integration (primary)

The headline route is calphy non-equilibrium thermodynamic integration (NETI): each phase's absolute Gibbs free energy is fixed by reversibly switching the Hamiltonian, over a finite alchemical path λ ∈ [0, 1], between the MLIP and an analytically solvable reference —

• solid — a Frenkel–Ladd Einstein crystal (atoms harmonically tethered to their lattice sites, whose free energy is known in closed form);
• liquid — an Uhlenbeck–Ford (UF) fluid reference, realised in LAMMPS as pair_style hybrid/scaled ... grace ... ufm (the UF potential mixed in and scaled out over the switching path alongside the GRACE potential).

For each (potential, element) the campaign runs 10 single-temperature free-energy kernels — 5 solid + 5 liquid — on the reduced-temperature grid

T = [0.60, 0.75, 0.90, 1.00, 1.10] × Tm_exp

widened to [0.40 … 1.00] × Tm_exp for a handful of very soft, low-melting elements where the default grid sits too close to (or above) the true Tm to bracket a crossing. Cells are ~2000–2700 atoms; each kernel runs N_equil = 10,000 equilibration steps followed by N_switch = 25,000 switching steps, n_iter = 1 (no repeat averaging — see the switching-convergence study for why 25k is adequate).

Each phase's five (T, G) points are fit linearly in T; Tm is the temperature where the solid and liquid G(T) lines cross. Two quality gates run before the fit:

• melt-corrupted solid / freeze-corrupted liquid points are trimmed. A solid kernel that spontaneously melts mid-switch (or a liquid kernel that freezes) no longer samples the phase it was assigned to, so its (T, G) point is dropped before the linear fit — this is the same trimming shown as hollow/dashed points on the /melting page's per-element free-energy panel.
• no_sign_change rows are extrapolations. If the surviving points on both branches never actually bracket a crossing (both grids monotonically diverge or converge without swapping order), the reported Tm is a linear extrapolation outside the sampled range, not an interpolated crossing. These rows are flagged no_sign_change = true throughout the site — excluded from the clean ensemble by default (the "Exclude no-sign-change extrapolations" toggle on /melting) and always shown low-confidence (hollow parity markers, an on-panel warning banner) wherever they appear.

2. Switching-convergence study

Longer switching paths dissipate less irreversible work, so a natural question is whether the production N_switch = 25,000 setting is simply too fast to have converged. To check, Mo and Re were rerun at three switching lengths plus a Richardson extrapolation to N_switch → ∞:

Quadrupling the switching time from 25k to 100k steps recovers only +70 K (Mo) and +97 K (Re) — a small fraction of the ~300–650 K gap to experiment — and extrapolating all the way to infinite switching time (zero dissipation) still leaves both elements ~10 % (Mo) / ~18 % (Re) low. Production therefore stays at N_switch = 25,000: the residual gap to experiment is real model-plus-method error (the MLIP's PES and/or the NETI free-energy protocol itself), not an artefact of insufficiently slow switching. This is the evidence behind the TI ensemble's systematic low bias reported in §4.

3. Coexistence (secondary)

The independent cross-check builds an explicit two-phase solid–liquid cell (~4000 atoms) and brackets the melting point directly, rather than via free energies:

• At each candidate temperature, a short NVE run is performed at 11 strains bracketing the equilibrium volume.
• The solid fraction is measured per strain by adaptive common-neighbour analysis (CNA), distinguishing the ordered (solid) band from the disordered (liquid) band across the interface.
• The equilibrium coexistence point is located by interface-stationary extrapolation as P → 0 (the strain/pressure at which the solid–liquid interface neither advances nor retreats).
• The temperature itself is found by bracket bisection: if the interface grows solid, the bracket's lower bound rises; if it melts through, the upper bound falls. The search is capped at 28 iterations — most production runs report hitting this cap rather than a tight converged bracket, and are flagged on-page (converged: false badge) rather than silently reported as exact.
• Engine: LAMMPS pair_style grace (no reference potential needed — coexistence only requires the MLIP itself).

Coexistence is run for a subset of (potential, element) pairs — where an explicit two-phase cell could be constructed and stably equilibrated — not the full elemental grid TI covers.

4. The two methods disagree by design

TI and coexistence are kept separate, not averaged, because they fail in opposite regimes and agreeing with both would only hide that:

• TI is robust for high-melting refractory metals (large, well-separated solid/liquid free-energy branches) but its harmonic/UF reference free energies strain for very soft, low-Tm elements, where the widened [0.40…1.00]×Tm_exp grid is itself a symptom of the difficulty.
• Coexistence needs a mechanically well-defined solid–liquid interface and struggles precisely where the solid is barely stable at all — the same soft-element regime, but for a different reason (no clean interface to bracket, not a bad free-energy reference).

On the trusted-4 clean ensemble (sign-change crossings only, no-sign-change extrapolations excluded), the TI headline is:

• MAE ≈ 21.6 % vs experimental Tm
• MBE ≈ −21.4 % — i.e. the bias is almost the whole error: TI runs systematically low, not scattered around the true value
• 27 of 55 elements land within ±20 % of experiment

(The campaign report's ti_results_trusted.csv instead pools the no-sign-change extrapolations into its per-element means, giving MAE ≈ 24.1 %, MBE ≈ −23.1 %, 29 of 61 within ±20 % over that larger, mixed-quality pool — the regen tool's cross-check certifies against these pooled numbers.)

The honest reading is that a uniform ~20+ % low bias is a real, characterized systematic offset (confirmed by the switching-convergence study above), not noise to be averaged away. The more reassuring number sits alongside it: the inter-potential spread — how much the potentials agree with each other — is only ~16 %, tighter than the 21.6 % MAE against experiment. The models under-predict Tm in close agreement with one another; they are more precise than they are accurate. Where TI and coexistence both ran for the same (potential, element), treat the pair as two independent, differently-biased estimates of the same quantity — not as a single averaged number.

5. Data caveats

• 11 tm-null TI rows. A handful of (potential, element) TI runs produced no usable crossing at all (not even an extrapolated one) — these carry a null tm and are excluded from every metric, not folded in as zero or as a failure count.
• Coexistence pressures are in bar, not GPa/kbar like the rest of the site's stress-adjacent quantities — a unit-convention difference to watch when cross-referencing the downloadable coex iteration tables against other tracks.
• Kr and Xe tm_exp come from the CRC Handbook, not the same experimental-compilation source used for the metallic elements: Kr 115.78 K, Xe 161.40 K.
• Per-iteration converged flags are unreliable — only the top-level run flag should be trusted. The coexistence bisection trace records a converged field at each iteration, but only the top-level, final value (the 28-iteration-cap badge shown on /melting) is a meaningful convergence signal; per-iteration flags in the raw iteration table should not be read as "converged as of this iteration."
• λ web paths are downsampled to ~300 points for the on-page switching-hysteresis plot (the λ drill-in tab under a selected element); the raw, full-resolution ~25,000-point λ path per kernel ships in the download bundle, not on the page.

See the leaderboard → · ← Methodology overview