From cd3aa790a334a754d6c7494409ce3e808f68fc8b Mon Sep 17 00:00:00 2001 From: Ahmet Gunaydin Date: Mon, 10 Aug 2026 09:11:14 +0700 Subject: [PATCH 01/13] =?UTF-8?q?test:=20per-element=20PySCF=20agreement?= =?UTF-8?q?=20harness=20=E2=80=94=20exposes=20wrong=20Li/Be/aug=20basis=20?= =?UTF-8?q?tables?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Phase 0 of docs/RUN-PLAN-24H-ELEMENTS.md. Builds the validation harness that the element-expansion run depends on, and immediately finds a real defect in the existing basis data. Gate 0.1 — d-convention, pinned with measurement. On H2O/cc-pVDZ PySCF gives -76.02679870 (spherical, nao=24) vs -76.02713907 (mol.cart=True, nao=25): a 0.340 mHa basis-convention gap, matching the note in ccpvdz-spherical.test.ts. Our default Cartesian path lands 1.9 uHa from cart-PySCF. The harness therefore compares each code path against its MATCHING convention; crossing them would fake a ~0.34 mHa error on every element carrying d functions. STO-3G is convention-free below the third row (asserted). Gate 0.2/0.3 — the finding. 49 cells (8 elements x 3 bases x 2 conventions, one hydride per element so two-center integrals actually run): 38 pass, 11 fail. Li . LiH . cc-pVDZ 1.293 / 1.300 mHa Li . LiH . aug 1.144 / 1.145 mHa Be . BeH2 . cc-pVDZ 0.287 / 0.367 mHa Be . BeH2 . aug 0.119 mHa O . H2O . aug 0.122 / 0.141 mHa F . HF . aug 0.141 / 0.171 mHa Shell counts are correct everywhere (nao matches); the numbers inside are wrong. Li's cc-pVDZ p-shell carries coefficients (0.0223, 0.1117, 0.2403) against PySCF's (0.0228, 0.1391, 0.5004) on identical exponents — not a normalization convention, the ratios follow no power of alpha. Li at 1.29 mHa breaches the repo's own <= 0.5 mHa HF-vs-PySCF gate by 2.6x and is 81% of chemical accuracy. Why 954 tests missed it: ccpvdz.test.ts and aug-ccpvdz.test.ts only ever run H2O, exercising H and O — two of the clean elements. The tests that do cover Li/Be assert ballparks: ccpvdz-firstrow.test.ts allows -8.10 < E_LiH < -7.85, a 250 mHa window ~200x wider than the error, under a header claiming "verified vs PySCF"; aug-ccpvdz-firstrow.test.ts asserts only variational ordering with no absolute reference at all. This commit is deliberately RED on those 11 cells — the failing test is the evidence, per the honest-negative discipline in RESEARCH.md. The repair follows. Co-Authored-By: Claude Opus 5 (1M context) --- docs/RUN-PLAN-24H-ELEMENTS.md | 254 +++ scripts/check-basis-vs-pyscf.py | 135 ++ scripts/dump-basis-constants.ts | 26 + scripts/run-element-ours.ts | 80 + scripts/run-element-reference.py | 133 ++ tests/chemistry/elements/pyscf-reference.json | 1589 +++++++++++++++++ .../elements/reference-agreement.test.ts | 82 + 7 files changed, 2299 insertions(+) create mode 100644 docs/RUN-PLAN-24H-ELEMENTS.md create mode 100644 scripts/check-basis-vs-pyscf.py create mode 100644 scripts/dump-basis-constants.ts create mode 100644 scripts/run-element-ours.ts create mode 100644 scripts/run-element-reference.py create mode 100644 tests/chemistry/elements/pyscf-reference.json create mode 100644 tests/chemistry/elements/reference-agreement.test.ts diff --git a/docs/RUN-PLAN-24H-ELEMENTS.md b/docs/RUN-PLAN-24H-ELEMENTS.md new file mode 100644 index 0000000..e4d9150 --- /dev/null +++ b/docs/RUN-PLAN-24H-ELEMENTS.md @@ -0,0 +1,254 @@ +# 24h autonomous run plan — break the periodic-table wall + +**Status:** proposed, not started +**Branch:** `feat/periodic-table-expansion` (create from `main`) +**Written:** 2026-08-10 + +--- + +## Mission + +Take `SUPPORTED_SYMBOLS` from 8 elements to 18, with every new element +validated against PySCF to sub-milliHartree and locked behind a +regression test. + +Current: `H, He, Li, Be, C, N, O, F` +Target: add `B, Ne, Na, Mg, Al, Si, P, S, Cl, Ar` + +This is a **data** problem, not an engine problem. The integral core +already handles d-shells (cc-pVDZ carbon carries `CCPVDZ_C_3D`), so +third-row elements need no new angular-momentum machinery — only more +contracted shells at L ≤ 2. + +## Non-goals + +Do **not** do any of these, even if they look tempting mid-run: + +- Refactor working code. The repo was audited 2026-08-08; it is not the problem. +- Add new methods, new UI, new pages, or new experiments. +- Deploy anything. Not Vercel, not the HF Space. (Standing user preference.) +- Push to `main` or open a PR. Leave the branch for human review. +- Chase performance. + +## Hard rules + +1. **Every step ends in a number a machine can compare.** No "looks + right", no visual inspection, no prose verdicts. If a step cannot be + expressed as `assert |ours − reference| < tol`, it is not a step. +2. **Hard-fail on `converged: false`.** Every harness must treat a + non-converged SCF/CCSD as a *failure*, never as a datum. This is the + single most likely way a long run produces confident garbage. +3. **One element per commit.** A bad stretch must be revertable without + losing the good ones. Trailer: + `Co-Authored-By: Claude Opus 5 (1M context) ` +4. **Serialize test runs.** Do not run vitest suites concurrently with + PySCF sweeps. The machine is a 16 GB M2 Pro; saturating it produces + timeout flakes that look exactly like real failures and will send the + run chasing phantoms. +5. **Never widen a tolerance to make a test pass.** If a bar cannot be + met, record the failure in the Phase 3 report and move to the next + element. A loosened tolerance is an undetectable lie. +6. **`uv run` inside `~/dev/ml-research` for anything PySCF.** Never bare + `pip` / `python`. + +--- + +## Phase 0 — the harness (do not skip; nothing after this is safe without it) + +### 0.1 Resolve the d-shell convention first + +`tests/chemistry/ccpvdz-spherical.test.ts` documents that the Cartesian +(6d) path sits **~0.4 mHa below** PySCF's default spherical (5d) result +for H₂O/cc-pVDZ — a basis-set convention difference, *not* a bug. + +Every element from Al onward has d functions in cc-pVDZ. If the harness +compares the Cartesian path against default PySCF, **every new element +will show a ~0.4 mHa phantom error**, and the run will spend hours +"fixing" correct code. + +**Gate 0.1:** the harness must either (a) compare the spherical path +against default PySCF, or (b) compare the Cartesian path against PySCF +with `mol.cart = True`. Prove it by reproducing the existing H₂O/cc-pVDZ +agreement *before* touching any new element. Write down which convention +was chosen at the top of the harness file. + +### 0.2 Generate basis tables — do not hand-transcribe + +Write `scripts/gen-basis-tables.py`: read PySCF's own basis library +(`pyscf.gto.basis.load`) and emit TypeScript constants in the exact +shape `src/chemistry/integrals.ts` already uses: + +```ts +export const STO3G_C_2P = { + alpha: [2.9412494, 0.6834831, 0.2222899] as const, + c: STO3G_L_2P_C, +}; +``` + +Hand-transcribing 10 elements × 3 basis sets is ~100 blocks of decimal +digits. That is precisely the task where a long autonomous run silently +introduces a typo that surfaces as a 3 mHa discrepancy 14 hours later. +Generating it removes the failure mode entirely. + +**Gate 0.2 (the important one):** run the generator against the **eight +elements that already exist** and diff against the current hand-written +constants. If it reproduces `C`, `N`, `O`, `F` to the last digit, the +generator is trustworthy for `S`, `P`, `Cl`. If it does not, fix the +generator — never the existing constants — until it does. + +This is what makes the rest of the run autonomous-safe: the tool is +validated against known-good data before it is trusted with new data. + +### 0.3 Two-level validation harness + +For each (element, basis), the harness runs **two** independent checks +that isolate the two distinct failure modes: + +| Level | What it compares | Catches | +|---|---|---| +| **A — engine** | our energy vs PySCF fed *our own primitives* via `gto.basis.parse` | bugs in our integrals/SCF for this element | +| **B — data** | our primitives vs PySCF's built-in basis table | wrong/garbled basis data | + +`scripts/run-pyscf-reference.py` already uses `gto.basis.parse` for the +matched s-only Li basis — same mechanism, reuse it. + +Interpretation is unambiguous: +- A passes, B fails → basis data is wrong. +- B passes, A fails → engine bug for this element. +- Both fail → start with B. + +**Gate 0.3:** both levels green for all 8 existing elements before +proceeding. Bars: HF ≤ 1e-9 Ha for Level A (same primitives should agree +to machine precision, not chemistry precision), ≤ 0.1 mHa for Level B. + +--- + +## Phase 1 — elements, in this order + +The order is deliberate: **structurally-familiar elements first**, so +that an early failure means "the harness is wrong", not "this element is +hard". Only after the pipeline is proven does it meet new shell +structure. + +### 1a. Second row — proves the pipeline (structure already supported) + +| El | Z | STO-3G shells | Why | +|---|---|---|---| +| **B** | 5 | 1s, 2s, 2p | BH₃/BF₃ — the canonical empty-orbital Lewis acid lesson. Also fills the one gap in row 2. | +| **Ne** | 10 | 1s, 2s, 2p | Completes row 2; noble-gas reference for the dispersion labs (`dispersion.ts` already exists). | + +These are byte-for-byte the same shell pattern as C/N/O/F. If either +fails, **stop and fix the harness** — do not proceed to row 3. + +### 1b. Third row — the actual unlock (new: 3s/3p shells) + +Ordered by curriculum value, highest first. + +| El | Z | Unlocks | +|---|---|---| +| **Cl** | 17 | **SN2.** The single most-taught mechanism in undergraduate organic. Also HCl, CH₃Cl, Cl⁻ as leaving group. | +| **S** | 16 | Thiols, H₂S, disulfide bridges, SO₂. Biochemistry entry point. | +| **P** | 15 | Phosphates, PH₃ — the DNA/ATP backbone. | +| **Si** | 14 | Silanes, SiO₂; the "why doesn't Si behave like C" comparison. | +| **Na** | 11 | Ionic bonding; NaCl with Cl already landed. | +| **Mg** | 12 | Mg²⁺, chlorophyll hand-wave, second ionic case. | +| **Al** | 13 | AlCl₃ Lewis acid; completes the row. | +| **Ar** | 18 | Van der Waals / dispersion labs; noble-gas dimer. | + +If the run stalls, **Cl, S, P are the three that matter.** Everything +after Si is completeness, not capability. + +### Per-element loop (identical for every element) + +1. Generate constants with the Phase-0 generator; append to `integrals.ts`. +2. Add the `case` to both switch blocks in `atoms.ts` (STO-3G and the + `heavyShells` cc-pVDZ path). Follow the existing formatting exactly. +3. Add the symbol to `SYMBOL_BY_Z` and `SUPPORTED_SYMBOLS` in `xyz.ts`. +4. Run Level A + Level B for all three basis sets. +5. Write `tests/chemistry/elements/.test.ts` pinning the validated + HF energy for at least one real molecule per basis set. +6. Run the **full** existing unit suite. A new element must not perturb + any existing number. If it does, that is a shared-state bug — stop + and report it, do not paper over it. +7. Commit. Move on. + +**Per-element gate:** steps 4–6 all green, or the element is reverted and +logged as failed. No partial elements on the branch. + +### Molecule choices for validation + +Prefer molecules with published reference geometries and non-trivial +bonding, one per element: + +`BH₃`, `Ne` (atom), `HCl` + `CH₃Cl`, `H₂S`, `PH₃`, `SiH₄`, `NaH`, +`MgH₂`, `AlH₃`, `Ar` (atom). + +Atoms alone are a weak test — they exercise no two-center integrals. Use +a hydride wherever one exists. + +--- + +## Phase 2 — gradient hardening + +`hf-gradient.ts`, `dft-gradient.ts`, `cphf.ts` and `optimizer.ts` all +exist. Nothing verifies they stay correct as elements are added. + +For every (new element, basis, method) combination: compare the analytic +gradient against a central finite-difference gradient of the energy. + +**Gate 2:** max component difference ≤ 1e-6 Ha/Bohr. Use a step of +1e-4 Bohr; if the FD gradient is itself noisy at that step, record the +noise floor rather than tightening the bar. + +This is a perfect autonomous oracle — the reference is generated from +the code's own energy function, so it needs no external tool and cannot +be argued with. + +--- + +## Phase 3 — the report (reserve the final hour) + +Write `docs/RUN-REPORT-.md` containing, in this order: + +1. **What could not be verified.** Put it first. This is the most + valuable output of the whole run — more valuable than the elements. +2. Elements landed, with the measured Level A / Level B deltas per basis. +3. Elements attempted and reverted, with the actual failure. +4. Any tolerance that was *considered* for loosening and why it wasn't. +5. Any existing test whose number moved, and the explanation. +6. What the next run should start with. + +Then update `LIMITATIONS.md` with the new element coverage and — this is +outstanding from the last session either way — the **anthracene n=246 +wrong-basin result**, which currently lives only in a spec comment. + +--- + +## Known traps (pre-empted here so the run doesn't rediscover them) + +- **The 0.4 mHa d-convention phantom.** See Gate 0.1. This will look + exactly like a real error on every element from Al onward. +- **`converged: false` treated as data.** Bit a lab script in the + 2026-08-08 session: the library reported honestly, the consumer code + printed unconverged CCSD numbers as results. Assert the flag. +- **Load-induced timeout flakes.** Three "failures" in the last full run + were a saturated machine, all passing in isolation. Before declaring a + test broken, re-run it alone. +- **`sed` on macOS is BSD sed.** `\?` silently does nothing. Use Python + for any regex rewriting of files. +- **Vite HMR reloads on artifacts.** `e2e/.artifacts/**` and friends are + already in `server.watch.ignored`; do not write scratch files anywhere + else under the repo root during a run. +- **`fciState` is H₂-only** (imports `buildH2Dense`). Not a general FCI + solver. Do not use it to validate anything else. + +--- + +## What this run does not accomplish + +It makes webgpu-q *capable* of the standard curriculum — SN2, phosphates, +ionic bonding, dispersion. It does not put the project in front of a +single teacher. Those are separate problems, and only the second one is +genuinely uncertain. Spend the 24h here because it is cheap in human time +and unblocks the education lane; do not mistake it for validating that +lane. diff --git a/scripts/check-basis-vs-pyscf.py b/scripts/check-basis-vs-pyscf.py new file mode 100644 index 0000000..5809ff5 --- /dev/null +++ b/scripts/check-basis-vs-pyscf.py @@ -0,0 +1,135 @@ +"""Gate 0.2 — prove PySCF's basis tables reproduce our hand-transcribed constants. + +The 24h element-expansion run (docs/RUN-PLAN-24H-ELEMENTS.md) generates +TypeScript basis constants from PySCF's basis library instead of +hand-transcribing ~100 blocks of decimal digits. That generator is only +trustworthy if it reproduces the EIGHT elements already in the repo, +which were transcribed independently from EMSL/Pople tables. + +This script compares, per (element, basis), the multiset of contracted +shells (exponents + coefficients) exported by src/chemistry/integrals.ts +against pyscf.gto.basis.load. Structural match => the generator can be +trusted for S, P, Cl. + +Usage: + npx --yes tsx scripts/dump-basis-constants.ts > /tmp/ours.json + uv run python scripts/check-basis-vs-pyscf.py /tmp/ours.json +""" +import json +import sys +from collections import defaultdict + +from pyscf import gto + +# Longest prefix first — AUG_CCPVDZ_ must win over CCPVDZ_. +BASIS_PREFIXES = [ + ("AUG_CCPVDZ_", "aug-cc-pvdz"), + ("CCPVDZ_", "cc-pvdz"), + ("STO3G_", "sto-3g"), +] +ELEMENTS = ["H", "HE", "LI", "BE", "C", "N", "O", "F"] +TOL = 1e-7 + + +def classify(name): + """NAME -> (basis, element) or None.""" + for prefix, basis in BASIS_PREFIXES: + if name.startswith(prefix): + rest = name[len(prefix):] + for el in sorted(ELEMENTS, key=len, reverse=True): + if rest == el or rest.startswith(el + "_"): + return basis, el.capitalize() + return None + return None + + +def pyscf_shells(basis, element): + """PySCF basis -> list of (alpha tuple, coeff tuple), general + contractions expanded column-by-column, zero-coefficient + primitives dropped (our tables store only the non-zero span).""" + raw = gto.basis.load(basis, element) + out = [] + for shell in raw: + rows = shell[1:] # [[exp, c0, c1, ...], ...] + ncol = len(rows[0]) - 1 + for col in range(ncol): + pairs = [(r[0], r[col + 1]) for r in rows if abs(r[col + 1]) > 0.0] + if not pairs: + continue + out.append((tuple(a for a, _ in pairs), tuple(c for _, c in pairs))) + return out + + +def close(x, y): + return len(x) == len(y) and all(abs(a - b) <= TOL * max(1.0, abs(a)) for a, b in zip(x, y)) + + +def match(ours, theirs): + """Greedy multiset match. Returns (matched, ours_only, theirs_only).""" + remaining = list(theirs) + matched, ours_only = [], [] + for a, c in ours: + for i, (a2, c2) in enumerate(remaining): + if close(a, a2) and close(c, c2): + matched.append((a, c)) + remaining.pop(i) + break + else: + ours_only.append((a, c)) + return matched, ours_only, remaining + + +def main(): + with open(sys.argv[1]) as fh: + dumped = json.load(fh) + + by_key = defaultdict(list) + unclassified = [] + for name, val in dumped.items(): + key = classify(name) + if key is None: + unclassified.append(name) + continue + by_key[key].append((name, tuple(val["alpha"]), tuple(val["c"]))) + + if unclassified: + print(f"note: {len(unclassified)} constant(s) not classified: {unclassified}\n") + + failures = 0 + for basis in ("sto-3g", "cc-pvdz", "aug-cc-pvdz"): + for el in [e.capitalize() for e in ELEMENTS]: + if basis == "aug-cc-pvdz": + # Our aug tables store ONLY the diffuse additions; the + # full aug set is cc-pVDZ + diffuse. + ours_raw = by_key.get(("cc-pvdz", el), []) + by_key.get((basis, el), []) + else: + ours_raw = by_key.get((basis, el), []) + if not ours_raw: + continue + ours = [(a, c) for _, a, c in ours_raw] + try: + theirs = pyscf_shells(basis, el) + except Exception as exc: # noqa: BLE001 + print(f"FAIL {basis:12s} {el:2s} PySCF load error: {exc}") + failures += 1 + continue + matched, ours_only, theirs_only = match(ours, theirs) + ok = not ours_only and not theirs_only + status = "OK " if ok else "FAIL" + print(f"{status} {basis:12s} {el:2s} matched {len(matched):2d}/{len(ours):2d}" + f" ours_only={len(ours_only)} pyscf_only={len(theirs_only)}") + if not ok: + failures += 1 + for a, c in ours_only: + print(f" ours only : alpha={a}\n c={c}") + for a, c in theirs_only: + print(f" pyscf only: alpha={a}\n c={c}") + + print() + print("GATE 0.2:", "PASS — PySCF reproduces our tables" if failures == 0 + else f"FAIL — {failures} (basis, element) cell(s) disagree") + return 1 if failures else 0 + + +if __name__ == "__main__": + sys.exit(main()) diff --git a/scripts/dump-basis-constants.ts b/scripts/dump-basis-constants.ts new file mode 100644 index 0000000..71cacfc --- /dev/null +++ b/scripts/dump-basis-constants.ts @@ -0,0 +1,26 @@ +// Dump every contracted-shell constant exported by integrals.ts as JSON. +// +// Consumed by scripts/check-basis-vs-pyscf.py (Gate 0.2 of +// docs/RUN-PLAN-24H-ELEMENTS.md): the basis-table generator is only +// trusted for new elements once PySCF's own tables are shown to +// reproduce the eight hand-transcribed elements already in the repo. +// +// Usage: npx --yes tsx scripts/dump-basis-constants.ts +import * as integrals from "../src/chemistry/integrals.js"; + +interface ShellConst { alpha: readonly number[]; c: readonly number[] } + +function isShellConst(v: unknown): v is ShellConst { + if (typeof v !== "object" || v === null) return false; + const o = v as Record; + return Array.isArray(o["alpha"]) && Array.isArray(o["c"]) + && (o["alpha"] as unknown[]).every((x) => typeof x === "number") + && (o["c"] as unknown[]).every((x) => typeof x === "number"); +} + +const out: Record = {}; +for (const [name, value] of Object.entries(integrals)) { + if (isShellConst(value)) out[name] = { alpha: value.alpha, c: value.c }; +} + +process.stdout.write(JSON.stringify(out, null, 2) + "\n"); diff --git a/scripts/run-element-ours.ts b/scripts/run-element-ours.ts new file mode 100644 index 0000000..2bc3957 --- /dev/null +++ b/scripts/run-element-ours.ts @@ -0,0 +1,80 @@ +// Element-validation — webgpu-q side. +// +// Reads the geometry + reference table emitted by +// scripts/run-element-reference.py and recomputes every cell with our own +// engine, under the MATCHING d-convention (Gate 0.1 of +// docs/RUN-PLAN-24H-ELEMENTS.md): +// +// PySCF spherical (default, 5d) <-> computeMolecularIntegrals(..., {spherical:true}) +// PySCF cartesian (cart=True) <-> computeMolecularIntegrals(...) [default path] +// +// A non-converged SCF is a FAILURE, never a datum. +// +// Usage: npx --yes tsx scripts/run-element-ours.ts +import { readFileSync, writeFileSync } from "node:fs"; +import { computeMolecularIntegrals } from "../src/chemistry/cg-molecular.js"; +import { moleculeToShellsNuclei, type Atom } from "../src/chemistry/atoms.js"; +import { runRHFSCF } from "../src/chemistry/hf-scf.js"; + +interface RefRow { + element: string; molecule: string; basis: string; convention: string; + atoms: { symbol: string; pos: number[] }[]; + ok: boolean; E_HF?: number; nao?: number; note?: string; +} + +const [, , refPath, outPath] = process.argv; +if (!refPath || !outPath) { + console.error("usage: tsx scripts/run-element-ours.ts "); + process.exit(2); +} + +const ref = JSON.parse(readFileSync(refPath, "utf8")) as { rows: RefRow[] }; +const rows: Record[] = []; + +for (const r of ref.rows) { + const atoms: Atom[] = r.atoms.map((a) => ({ + symbol: a.symbol as Atom["symbol"], + pos: [a.pos[0]!, a.pos[1]!, a.pos[2]!], + })); + const spherical = r.convention === "spherical"; + const out: Record = { + element: r.element, molecule: r.molecule, basis: r.basis, + convention: r.convention, ref_ok: r.ok, E_ref: r.E_HF ?? null, + }; + try { + const { shells, nuclei, nElectrons } = moleculeToShellsNuclei(atoms, r.basis); + const integrals = computeMolecularIntegrals(shells, nuclei, { spherical }); + const hf = runRHFSCF(integrals, nElectrons, { + useDIIS: true, maxIter: 300, energyTol: 1e-11, densityTol: 1e-9, + }); + if (!hf.converged) { + out["ok"] = false; + out["note"] = "our SCF did not converge"; + } else { + out["ok"] = true; + out["E_ours"] = hf.energy; + out["n"] = integrals.n; + if (r.ok && r.E_HF !== undefined) { + out["delta_Ha"] = hf.energy - r.E_HF; + out["delta_mHa"] = (hf.energy - r.E_HF) * 1000; + out["nao_match"] = r.nao === integrals.n; + } + } + } catch (e) { + out["ok"] = false; + out["note"] = e instanceof Error ? e.message : String(e); + } + rows.push(out); + + const d = out["delta_mHa"]; + const tag = out["ok"] === true + ? (typeof d === "number" ? `${d >= 0 ? "+" : ""}${d.toFixed(6)} mHa` : "no ref") + : `FAIL: ${String(out["note"]).slice(0, 60)}`; + const nm = out["nao_match"] === false ? " [nao MISMATCH]" : ""; + console.log( + `${r.molecule.padEnd(5)} ${r.basis.padEnd(12)} ${r.convention.padEnd(9)} ${tag}${nm}`, + ); +} + +writeFileSync(outPath, JSON.stringify({ rows }, null, 2)); +console.log(`\nwrote ${outPath} (${rows.length} rows)`); diff --git a/scripts/run-element-reference.py b/scripts/run-element-reference.py new file mode 100644 index 0000000..e60e47a --- /dev/null +++ b/scripts/run-element-reference.py @@ -0,0 +1,133 @@ +"""Element-validation reference — PySCF side. + +Emits, for one small molecule per supported element, the RHF energy under +BOTH d-conventions: + + spherical (PySCF default, 5d) <-> computeMolecularIntegrals(..., {spherical:true}) + cartesian (mol.cart = True, 6d) <-> computeMolecularIntegrals(...) [default path] + +Gate 0.1 of docs/RUN-PLAN-24H-ELEMENTS.md: each webgpu-q code path must be +compared against its MATCHING convention. Comparing the Cartesian path to +default (spherical) PySCF produces a ~0.34 mHa phantom error on every +element that carries d functions -- a basis-set convention difference, not +a code bug (see tests/chemistry/ccpvdz-spherical.test.ts). + +Geometries are emitted alongside the energies so the TypeScript side +consumes the exact same coordinates -- no double transcription. + +Usage: + cd ~/dev/ml-research && uv run python \ + /path/to/webgpu-q/scripts/run-element-reference.py --out /tmp/ref.json +""" +import argparse +import json +import math + +from pyscf import gto, scf + +SCF_TOL = 1e-12 + + +def tetrahedral(bond): + a = bond / math.sqrt(3.0) + return [(a, a, a), (a, -a, -a), (-a, a, -a), (-a, -a, a)] + + +def bent(bond, angle_deg): + half = math.radians(angle_deg / 2.0) + x, z = bond * math.sin(half), bond * math.cos(half) + return [(x, 0.0, z), (-x, 0.0, z)] + + +def pyramidal(bond, angle_deg): + """3 equivalent bonds, given X-H distance and H-X-H angle.""" + ha = math.radians(angle_deg / 2.0) + # in-plane radius r and height h satisfying the H-X-H angle + r = bond * math.sin(ha) / math.sin(math.radians(60.0)) + h = math.sqrt(max(bond * bond - r * r, 0.0)) + return [(r * math.cos(math.radians(90 + 120 * k)), + r * math.sin(math.radians(90 + 120 * k)), -h) for k in range(3)] + + +# element -> (name, atoms). One molecule per element; hydrides preferred +# over bare atoms so two-center integrals are actually exercised. +MOLECULES = { + "H": ("H2", [("H", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.7414))]), + "He": ("He", [("He", (0.0, 0.0, 0.0))]), + "Li": ("LiH", [("Li", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 1.5949))]), + "Be": ("BeH2", [("Be", (0.0, 0.0, 0.0)), + ("H", (0.0, 0.0, 1.3264)), ("H", (0.0, 0.0, -1.3264))]), + "C": ("CH4", [("C", (0.0, 0.0, 0.0))] + + [("H", p) for p in tetrahedral(1.0870)]), + "N": ("NH3", [("N", (0.0, 0.0, 0.0))] + + [("H", p) for p in pyramidal(1.0124, 106.67)]), + "O": ("H2O", [("O", (0.0, 0.0, 0.0))] + + [("H", p) for p in bent(0.9572, 104.52)]), + "F": ("HF", [("F", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.9168))]), +} + +BASES = ["sto-3g", "cc-pvdz", "aug-cc-pvdz"] + +# webgpu-q uses an s-only STO-3G for lithium (documented; see +# scripts/run-pyscf-reference.py). Feed PySCF the matching basis so the +# comparison stays apples-to-apples. +LI_S_ONLY_STO3G = """ +Li S + 16.1195750 0.15432897 + 2.9362007 0.53532814 + 0.7946505 0.44463454 +Li S + 0.6362897 -0.09996723 + 0.1478601 0.39951283 + 0.0480887 0.70011547 +""" + + +def basis_for(symbols, basis): + if basis == "sto-3g" and "Li" in symbols: + return {s: (gto.basis.parse(LI_S_ONLY_STO3G) if s == "Li" else "sto-3g") + for s in set(symbols)} + return basis + + +def main(): + ap = argparse.ArgumentParser() + ap.add_argument("--out", required=True) + args = ap.parse_args() + + results = [] + for el, (name, atoms) in MOLECULES.items(): + symbols = [a for a, _ in atoms] + for basis in BASES: + for conv, cart in (("spherical", False), ("cartesian", True)): + row = { + "element": el, "molecule": name, "basis": basis, + "convention": conv, + "atoms": [{"symbol": s, "pos": list(p)} for s, p in atoms], + } + try: + mol = gto.M(atom=[[s, tuple(p)] for s, p in atoms], + basis=basis_for(symbols, basis), + unit="Angstrom", verbose=0) + mol.cart = cart + mol.build() + mf = scf.RHF(mol) + mf.conv_tol = SCF_TOL + e = mf.kernel() + if not mf.converged: + row.update(ok=False, note="SCF did not converge") + else: + row.update(ok=True, E_HF=float(e), nao=int(mol.nao_nr())) + except Exception as exc: # noqa: BLE001 + row.update(ok=False, note=str(exc)) + results.append(row) + status = f"{row.get('E_HF', float('nan')):.10f}" if row.get("ok") else "FAIL" + print(f"{name:5s} {basis:12s} {conv:9s} {status}") + + with open(args.out, "w") as fh: + json.dump({"scf_tol": SCF_TOL, "rows": results}, fh, indent=2) + print(f"\nwrote {args.out} ({len(results)} rows)") + + +if __name__ == "__main__": + main() diff --git a/tests/chemistry/elements/pyscf-reference.json b/tests/chemistry/elements/pyscf-reference.json new file mode 100644 index 0000000..894ff55 --- /dev/null +++ b/tests/chemistry/elements/pyscf-reference.json @@ -0,0 +1,1589 @@ +{ + "scf_tol": 1e-12, + "rows": [ + { + "element": "H", + "molecule": "H2", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.7414 + ] + } + ], + "ok": true, + "E_HF": -1.116684387085341, + "nao": 2 + }, + { + "element": "H", + "molecule": "H2", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.7414 + ] + } + ], + "ok": true, + "E_HF": -1.116684387085341, + "nao": 2 + }, + { + "element": "H", + "molecule": "H2", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.7414 + ] + } + ], + "ok": true, + "E_HF": -1.128714959029653, + "nao": 10 + }, + { + "element": "H", + "molecule": "H2", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.7414 + ] + } + ], + "ok": true, + "E_HF": -1.128714959029653, + "nao": 10 + }, + { + "element": "H", + "molecule": "H2", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.7414 + ] + } + ], + "ok": true, + "E_HF": -1.1287933486408193, + "nao": 18 + }, + { + "element": "H", + "molecule": "H2", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.7414 + ] + } + ], + "ok": true, + "E_HF": -1.1287933486408193, + "nao": 18 + }, + { + "element": "He", + "molecule": "He", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "He", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -2.807783957539974, + "nao": 1 + }, + { + "element": "He", + "molecule": "He", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "He", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -2.807783957539974, + "nao": 1 + }, + { + "element": "He", + "molecule": "He", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "He", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -2.85516047724274, + "nao": 5 + }, + { + "element": "He", + "molecule": "He", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "He", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -2.85516047724274, + "nao": 5 + }, + { + "element": "He", + "molecule": "He", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "He", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -2.8557046677104294, + "nao": 9 + }, + { + "element": "He", + "molecule": "He", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "He", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -2.8557046677104294, + "nao": 9 + }, + { + "element": "Li", + "molecule": "LiH", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "Li", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.5949 + ] + } + ], + "ok": true, + "E_HF": -7.804244236540904, + "nao": 3 + }, + { + "element": "Li", + "molecule": "LiH", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Li", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.5949 + ] + } + ], + "ok": true, + "E_HF": -7.804244236540904, + "nao": 3 + }, + { + "element": "Li", + "molecule": "LiH", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Li", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.5949 + ] + } + ], + "ok": true, + "E_HF": -7.983615274790943, + "nao": 19 + }, + { + "element": "Li", + "molecule": "LiH", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Li", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.5949 + ] + } + ], + "ok": true, + "E_HF": -7.983650116757452, + "nao": 20 + }, + { + "element": "Li", + "molecule": "LiH", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Li", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.5949 + ] + } + ], + "ok": true, + "E_HF": -7.984156735583572, + "nao": 32 + }, + { + "element": "Li", + "molecule": "LiH", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Li", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.5949 + ] + } + ], + "ok": true, + "E_HF": -7.984213784238851, + "nao": 34 + }, + { + "element": "Be", + "molecule": "BeH2", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "Be", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.3264 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.3264 + ] + } + ], + "ok": true, + "E_HF": -15.56031234281192, + "nao": 7 + }, + { + "element": "Be", + "molecule": "BeH2", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Be", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.3264 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.3264 + ] + } + ], + "ok": true, + "E_HF": -15.56031234281192, + "nao": 7 + }, + { + "element": "Be", + "molecule": "BeH2", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Be", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.3264 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.3264 + ] + } + ], + "ok": true, + "E_HF": -15.767272729574648, + "nao": 24 + }, + { + "element": "Be", + "molecule": "BeH2", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Be", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.3264 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.3264 + ] + } + ], + "ok": true, + "E_HF": -15.76735064710549, + "nao": 25 + }, + { + "element": "Be", + "molecule": "BeH2", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Be", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.3264 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.3264 + ] + } + ], + "ok": true, + "E_HF": -15.768464084930171, + "nao": 41 + }, + { + "element": "Be", + "molecule": "BeH2", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Be", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.3264 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.3264 + ] + } + ], + "ok": true, + "E_HF": -15.76849911667396, + "nao": 43 + }, + { + "element": "C", + "molecule": "CH4", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "C", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + 0.6275797426091232, + 0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + -0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + 0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + -0.6275797426091232, + 0.6275797426091232 + ] + } + ], + "ok": true, + "E_HF": -39.72681011479772, + "nao": 9 + }, + { + "element": "C", + "molecule": "CH4", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "C", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + 0.6275797426091232, + 0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + -0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + 0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + -0.6275797426091232, + 0.6275797426091232 + ] + } + ], + "ok": true, + "E_HF": -39.72681011479772, + "nao": 9 + }, + { + "element": "C", + "molecule": "CH4", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "C", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + 0.6275797426091232, + 0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + -0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + 0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + -0.6275797426091232, + 0.6275797426091232 + ] + } + ], + "ok": true, + "E_HF": -40.19867261540841, + "nao": 34 + }, + { + "element": "C", + "molecule": "CH4", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "C", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + 0.6275797426091232, + 0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + -0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + 0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + -0.6275797426091232, + 0.6275797426091232 + ] + } + ], + "ok": true, + "E_HF": -40.19874250719319, + "nao": 35 + }, + { + "element": "C", + "molecule": "CH4", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "C", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + 0.6275797426091232, + 0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + -0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + 0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + -0.6275797426091232, + 0.6275797426091232 + ] + } + ], + "ok": true, + "E_HF": -40.19961762393349, + "nao": 59 + }, + { + "element": "C", + "molecule": "CH4", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "C", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + 0.6275797426091232, + 0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + 0.6275797426091232, + -0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + 0.6275797426091232, + -0.6275797426091232 + ] + }, + { + "symbol": "H", + "pos": [ + -0.6275797426091232, + -0.6275797426091232, + 0.6275797426091232 + ] + } + ], + "ok": true, + "E_HF": -40.19966503955351, + "nao": 61 + }, + { + "element": "N", + "molecule": "NH3", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "N", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 5.741863175669773e-17, + 0.9377174185516153, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8120871060368641, + -0.46885870927580775, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8120871060368638, + -0.4688587092758081, + -0.3816278330322548 + ] + } + ], + "ok": true, + "E_HF": -55.45408725399906, + "nao": 8 + }, + { + "element": "N", + "molecule": "NH3", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "N", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 5.741863175669773e-17, + 0.9377174185516153, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8120871060368641, + -0.46885870927580775, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8120871060368638, + -0.4688587092758081, + -0.3816278330322548 + ] + } + ], + "ok": true, + "E_HF": -55.45408725399906, + "nao": 8 + }, + { + "element": "N", + "molecule": "NH3", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "N", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 5.741863175669773e-17, + 0.9377174185516153, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8120871060368641, + -0.46885870927580775, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8120871060368638, + -0.4688587092758081, + -0.3816278330322548 + ] + } + ], + "ok": true, + "E_HF": -56.19562001585709, + "nao": 29 + }, + { + "element": "N", + "molecule": "NH3", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "N", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 5.741863175669773e-17, + 0.9377174185516153, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8120871060368641, + -0.46885870927580775, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8120871060368638, + -0.4688587092758081, + -0.3816278330322548 + ] + } + ], + "ok": true, + "E_HF": -56.19573475296943, + "nao": 30 + }, + { + "element": "N", + "molecule": "NH3", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "N", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 5.741863175669773e-17, + 0.9377174185516153, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8120871060368641, + -0.46885870927580775, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8120871060368638, + -0.4688587092758081, + -0.3816278330322548 + ] + } + ], + "ok": true, + "E_HF": -56.20538665261573, + "nao": 50 + }, + { + "element": "N", + "molecule": "NH3", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "N", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 5.741863175669773e-17, + 0.9377174185516153, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8120871060368641, + -0.46885870927580775, + -0.3816278330322548 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8120871060368638, + -0.4688587092758081, + -0.3816278330322548 + ] + } + ], + "ok": true, + "E_HF": -56.20557279808558, + "nao": 52 + }, + { + "element": "O", + "molecule": "H2O", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "O", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.7569503272636612, + 0.0, + 0.585882276618295 + ] + }, + { + "symbol": "H", + "pos": [ + -0.7569503272636612, + 0.0, + 0.585882276618295 + ] + } + ], + "ok": true, + "E_HF": -74.96292824643288, + "nao": 7 + }, + { + "element": "O", + "molecule": "H2O", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "O", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.7569503272636612, + 0.0, + 0.585882276618295 + ] + }, + { + "symbol": "H", + "pos": [ + -0.7569503272636612, + 0.0, + 0.585882276618295 + ] + } + ], + "ok": true, + "E_HF": -74.96292824643288, + "nao": 7 + }, + { + "element": "O", + "molecule": "H2O", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "O", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.7569503272636612, + 0.0, + 0.585882276618295 + ] + }, + { + "symbol": "H", + "pos": [ + -0.7569503272636612, + 0.0, + 0.585882276618295 + ] + } + ], + "ok": true, + "E_HF": -76.02679869746876, + "nao": 24 + }, + { + "element": "O", + "molecule": "H2O", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "O", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.7569503272636612, + 0.0, + 0.585882276618295 + ] + }, + { + "symbol": "H", + "pos": [ + -0.7569503272636612, + 0.0, + 0.585882276618295 + ] + } + ], + "ok": true, + "E_HF": -76.02713907180959, + "nao": 25 + }, + { + "element": "O", + "molecule": "H2O", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "O", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.7569503272636612, + 0.0, + 0.585882276618295 + ] + }, + { + "symbol": "H", + "pos": [ + -0.7569503272636612, + 0.0, + 0.585882276618295 + ] + } + ], + "ok": true, + "E_HF": -76.04142796052426, + "nao": 41 + }, + { + "element": "O", + "molecule": "H2O", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "O", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.7569503272636612, + 0.0, + 0.585882276618295 + ] + }, + { + "symbol": "H", + "pos": [ + -0.7569503272636612, + 0.0, + 0.585882276618295 + ] + } + ], + "ok": true, + "E_HF": -76.04198323019492, + "nao": 43 + }, + { + "element": "F", + "molecule": "HF", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "F", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.9168 + ] + } + ], + "ok": true, + "E_HF": -98.57075759161432, + "nao": 6 + }, + { + "element": "F", + "molecule": "HF", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "F", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.9168 + ] + } + ], + "ok": true, + "E_HF": -98.57075759161432, + "nao": 6 + }, + { + "element": "F", + "molecule": "HF", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "F", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.9168 + ] + } + ], + "ok": true, + "E_HF": -100.01941870308725, + "nao": 19 + }, + { + "element": "F", + "molecule": "HF", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "F", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.9168 + ] + } + ], + "ok": true, + "E_HF": -100.01983707569883, + "nao": 20 + }, + { + "element": "F", + "molecule": "HF", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "F", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.9168 + ] + } + ], + "ok": true, + "E_HF": -100.03347419310293, + "nao": 32 + }, + { + "element": "F", + "molecule": "HF", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "F", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 0.9168 + ] + } + ], + "ok": true, + "E_HF": -100.03442833006022, + "nao": 34 + } + ] +} \ No newline at end of file diff --git a/tests/chemistry/elements/reference-agreement.test.ts b/tests/chemistry/elements/reference-agreement.test.ts new file mode 100644 index 0000000..e344cbb --- /dev/null +++ b/tests/chemistry/elements/reference-agreement.test.ts @@ -0,0 +1,82 @@ +// Per-element RHF agreement against a committed PySCF reference table. +// +// This is the permanent verifier for element coverage: every supported +// element is exercised in a real molecule (hydrides, not bare atoms, so +// two-center integrals actually run) across every supported basis, under +// BOTH d-conventions. +// +// Convention pairing (Gate 0.1, docs/RUN-PLAN-24H-ELEMENTS.md): +// PySCF spherical (default, 5d) <-> computeMolecularIntegrals(..., {spherical:true}) +// PySCF cartesian (mol.cart=True) <-> computeMolecularIntegrals(...) [default path] +// Crossing the conventions produces a ~0.34 mHa phantom disagreement on +// any element carrying d functions -- a basis-set difference, not a bug. +// STO-3G has no d functions below the third row, so both conventions +// coincide there (asserted below). +// +// Regenerate the fixture with: +// cd ~/dev/ml-research && uv run python \ +// /scripts/run-element-reference.py --out /tests/chemistry/elements/pyscf-reference.json +import { describe, expect, test } from "vitest"; +import { readFileSync } from "node:fs"; +import { computeMolecularIntegrals } from "../../../src/chemistry/cg-molecular.js"; +import { moleculeToShellsNuclei, type Atom, type BasisName } from "../../../src/chemistry/atoms.js"; +import { runRHFSCF } from "../../../src/chemistry/hf-scf.js"; + +interface RefRow { + element: string; molecule: string; basis: string; convention: string; + atoms: { symbol: string; pos: number[] }[]; + ok: boolean; E_HF?: number; nao?: number; +} + +const fixture = JSON.parse( + readFileSync(new URL("./pyscf-reference.json", import.meta.url), "utf8"), +) as { rows: RefRow[] }; + +// Same primitives + same convention should agree far tighter than +// chemical accuracy. 0.1 mHa is the repo's existing HF-vs-PySCF bar +// (CLAUDE.md); anything near it means a real basis-data difference. +const BAR_mHa = 0.1; + +function ourEnergy(row: RefRow): { energy: number; n: number; converged: boolean } { + const atoms: Atom[] = row.atoms.map((a) => ({ + symbol: a.symbol as Atom["symbol"], + pos: [a.pos[0]!, a.pos[1]!, a.pos[2]!], + })); + const { shells, nuclei, nElectrons } = moleculeToShellsNuclei(atoms, row.basis as BasisName); + const integrals = computeMolecularIntegrals(shells, nuclei, { + spherical: row.convention === "spherical", + }); + const hf = runRHFSCF(integrals, nElectrons, { + useDIIS: true, maxIter: 300, energyTol: 1e-11, densityTol: 1e-9, + }); + return { energy: hf.energy, n: integrals.n, converged: hf.converged }; +} + +describe("Element coverage: RHF vs PySCF, convention-matched", () => { + for (const row of fixture.rows) { + if (!row.ok) continue; + const label = `${row.element} · ${row.molecule} · ${row.basis} · ${row.convention}`; + test(`${label}: |ΔE| ≤ ${BAR_mHa} mHa`, () => { + const { energy, n, converged } = ourEnergy(row); + // A non-converged SCF is a failure, never a datum. + expect(converged).toBe(true); + expect(Number.isFinite(energy)).toBe(true); + expect(n).toBe(row.nao); + const dmHa = (energy - row.E_HF!) * 1000; + expect(Math.abs(dmHa)).toBeLessThan(BAR_mHa); + }, 300_000); + } + + test("STO-3G is d-convention-free for rows 1-2 (reference sanity)", () => { + const sph = new Map(); + const cart = new Map(); + for (const r of fixture.rows) { + if (r.basis !== "sto-3g" || !r.ok) continue; + (r.convention === "spherical" ? sph : cart).set(r.molecule, r.E_HF!); + } + expect(sph.size).toBeGreaterThan(0); + for (const [mol, e] of sph) { + expect(cart.get(mol)).toBeCloseTo(e, 12); + } + }); +}); From 77e0597d4edb9af6e9bc9ddaa93c41b0dab0dfb5 Mon Sep 17 00:00:00 2001 From: Ahmet Gunaydin Date: Mon, 10 Aug 2026 09:12:55 +0700 Subject: [PATCH 02/13] docs: move the anthracene wrong-basin result into LIMITATIONS.md MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit It was written up in CHANGELOG.md release notes but absent from LIMITATIONS.md — the document that is supposed to enumerate what this code cannot do. A reader auditing the project's limits would not have found it. The result deserves the prominence: on anthracene (C14H10, cc-pVDZ, n=246) the damping recipe that rescues the default-DIIS +5352 Ha divergence converges cleanly to E ~= -880 Ha against a literature HF/cc-pVDZ value of ~-537 Ha. It reports converged: true and returns a finite, plausible number. Being more negative than the true ground state is the only tell. Filed under section 3 (Known SCF / CC failure modes) and flagged as the most dangerous entry there, because unlike the others it does not announce itself — no automatic check catches it, since the SCF is genuinely converged and only an external reference reveals the error. Co-Authored-By: Claude Opus 5 (1M context) --- LIMITATIONS.md | 19 +++++++++++++++++++ 1 file changed, 19 insertions(+) diff --git a/LIMITATIONS.md b/LIMITATIONS.md index 2015981..aa8b3a8 100644 --- a/LIMITATIONS.md +++ b/LIMITATIONS.md @@ -102,6 +102,25 @@ will silently truncate large dispatches. - **TDA-DFT / DFT-gradient with spherical-d** — refuses with a clear error today; proper fix is Cartesian → spherical transform on grid values. Documented. +- **Wrong-basin SCF at PAH scale — converges to a confidently wrong + answer.** This is the most dangerous failure mode in the list, because + it does not announce itself: `converged: true`, a clean stationary + point, and a plausible-looking number. On anthracene (C₁₄H₁₀, cc-pVDZ, + n = 246) the damping recipe (`damping=0.2`, `diisStartIter=8`) that + rescues the default-DIIS +5352 Ha divergence lands instead at + **E ≈ −880 Ha, against a literature HF/cc-pVDZ value of ~−537 Ha**. + Being *more negative* than the true ground state is the tell: the + damped warm-up steers the density into a non-physical orbital + occupation that is variationally lower within the SCF ansatz but is + not the ground-state singlet. The architecture and the convergence + machinery work end-to-end; **basin selection** is the open problem. + Real fixes are MOM (maximum overlap method) to preserve orbital + ordering across iterations, SOSCF for second-order convergence to the + nearest stationary point, or a SAD initial guess. Until one of those + lands, treat any large-PAH SCF energy from this code as unvalidated + regardless of the `converged` flag — and note that no automatic check + currently catches this, since the energy is finite, the SCF is + converged, and only an external reference reveals the error. --- From 2f35960c93f0d27cea3204960d35149a79390ff8 Mon Sep 17 00:00:00 2001 From: Ahmet Gunaydin Date: Mon, 10 Aug 2026 09:24:10 +0700 Subject: [PATCH 03/13] fix: repair 20 wrong basis-set constants across 7 elements MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Turns the red 49-cell agreement test from the previous commit green. Every correction is against pyscf.gto.basis.load, verified twice: by hand-diffing the tables, and independently by a generator that reconstructs the constants from PySCF and reproduces the same defect list. The generator measured Li cc-pVDZ at -1.2927 mHa against the hand-measured 1.293 mHa — three-digit agreement from two independent paths. cc-pVDZ: CCPVDZ_LI_2S_P 0.0285000 -> 0.0280500 (digit transposition) CCPVDZ_LI_2P.c wholly wrong coefficients CCPVDZ_BE_1S.c tail transcribed at reduced precision CCPVDZ_BE_2S.c elements 4-8 substantively wrong (different source) CCPVDZ_BE_2P.c last coefficient off by 5% CCPVDZ_N_1S.c 0.278722/0.448581 -> 0.278723/0.448540 CCPVDZ_H_1S non-canonical exponents 13.0107010/1.9622572/0.4445298 -> Dunning's published 13.01/1.962/0.4446 CCPVDZ_H_2S 0.1219496 -> 0.1220000 aug-cc-pVDZ diffuse: of the eight elements only H and C were correct. He, Li, Be, N, O and F all carried wrong exponents; all 16 values replaced. Previous values are recorded in comments at each site. Energy effect, worst cell per element (vs convention-matched PySCF): Li 1.293 mHa -> under 0.1 Be 0.367 -> under 0.1 O 0.141 -> under 0.1 F 0.171 -> under 0.1 All 49 cells (8 elements x 3 bases x 2 d-conventions) now pass a 0.1 mHa bar; the previous ceiling was 1.29 mHa. The aug diffuse errors on N/O/F sat ABOVE the true HF energy (+0.09 to +0.14 mHa). That is variationally impossible for a correct basis and is the tell that these were wrong data rather than an alternative convention. Full suite green: 140 files, 1003 passed, 1 skipped. No pinned number downstream moved enough to break — the loose bars that hid this defect also absorb its repair. Adds scripts/gen-basis-tables.py, which emits these constants from PySCF for all 18 elements so the remaining work is generated, not transcribed. Co-Authored-By: Claude Opus 5 (1M context) --- scripts/gen-basis-tables.py | 428 ++++++++++++++++++++++++++++++++++++ src/chemistry/integrals.ts | 82 ++++--- 2 files changed, 485 insertions(+), 25 deletions(-) create mode 100644 scripts/gen-basis-tables.py diff --git a/scripts/gen-basis-tables.py b/scripts/gen-basis-tables.py new file mode 100644 index 0000000..ce68e27 --- /dev/null +++ b/scripts/gen-basis-tables.py @@ -0,0 +1,428 @@ +"""Gate 0.2 — generate TypeScript basis-set constants from PySCF's own tables. + +docs/RUN-PLAN-24H-ELEMENTS.md forbids hand-transcribing ~100 blocks of +decimal digits for 10 new elements x 3 basis sets. This script reads +`pyscf.gto.basis.load` and emits constants in the exact shape +`src/chemistry/integrals.ts` already uses: + + export const STO3G_C_2P = { + alpha: [2.9412494, 0.6834831, 0.2222899] as const, + c: STO3G_L_2P_C, + }; + +It also runs the trust gate: the SAME code path is diffed, name by name +and number by number, against the constants currently exported by +integrals.ts (dumped via scripts/dump-basis-constants.ts). The generator +is only trusted for Na-Ar once it reproduces the known-clean cells +(C/N/O/F in sto-3g + cc-pvdz, C in aug-cc-pvdz) to ~1e-6 relative. + +Usage: + npx --yes tsx scripts/dump-basis-constants.ts > /tmp/ours.json + cd ~/dev/ml-research && uv run python \ + /path/to/webgpu-q/scripts/gen-basis-tables.py \ + --verify /tmp/ours.json --out /tmp/generated + +Flags: + --out DIR write .ts fragments + all.ts into DIR + --verify FILE diff against a dump-basis-constants.ts JSON + --structure print the per-element shell-structure table + --elements ... restrict the element list (default: all 18) +""" + +from __future__ import annotations + +import argparse +import json +import os +import sys +from collections import OrderedDict + +from pyscf import gto + +# ── configuration ──────────────────────────────────────────────── + +ELEMENTS = [ + "H", "He", + "Li", "Be", "B", "C", "N", "O", "F", "Ne", + "Na", "Mg", "Al", "Si", "P", "S", "Cl", "Ar", +] + +BASES = ["sto-3g", "cc-pvdz", "aug-cc-pvdz"] + +TS_PREFIX = {"sto-3g": "STO3G", "cc-pvdz": "CCPVDZ", "aug-cc-pvdz": "AUG_CCPVDZ"} + +LSYM = "SPDFGH" + +# Shell labels, per (basis, row), in PySCF's own ordering within each l. +# Reproduces the convention already in integrals.ts exactly: +# H cc-pVDZ : 1S, 2S, 2P +# C cc-pVDZ : 1S, 2S, 2S_P, 2P, 2P_P, 3D (trailing _P = "prime") +# and extends it to row 3, where cc-pVDZ is [4s,3p,1d]. +# +# Row key: 1 = H/He, 2 = Li..Ne, 3 = Na..Ar. +LABELS = { + ("sto-3g", 1): {0: ["1S"]}, + ("sto-3g", 2): {0: ["1S", "2S"], 1: ["2P"]}, + ("sto-3g", 3): {0: ["1S", "2S", "3S"], 1: ["2P", "3P"]}, + ("cc-pvdz", 1): {0: ["1S", "2S"], 1: ["2P"]}, + ("cc-pvdz", 2): {0: ["1S", "2S", "2S_P"], 1: ["2P", "2P_P"], 2: ["3D"]}, + ("cc-pvdz", 3): {0: ["1S", "2S", "3S", "3S_P"], 1: ["2P", "3P", "3P_P"], 2: ["3D"]}, +} + +# STO-3G contraction coefficients are universal across the whole row set: +# only the exponents change from atom to atom. integrals.ts already +# factors the first three out as shared constants; the M-shell pair is +# new (row 3) and has to be emitted. +SHARED_STO3G = OrderedDict([ + ("1S", ("STO3G_S_C", True)), # (const name, already exists in integrals.ts) + ("2S", ("STO3G_L_2S_C", True)), + ("2P", ("STO3G_L_2P_C", True)), + ("3S", ("STO3G_M_3S_C", False)), + ("3P", ("STO3G_M_3P_C", False)), +]) + +SHARED_TOL = 1e-6 + + +def row_of(symbol: str) -> int: + z = gto.charge(symbol) + if z <= 2: + return 1 + if z <= 10: + return 2 + return 3 + + +# ── PySCF -> flat contracted functions ─────────────────────────── + +def contracted_functions(basis: str, symbol: str): + """[(l, alphas, coeffs)] in PySCF order. + + General contractions are expanded column by column: each coefficient + column over a shared exponent block is one contracted function. + Primitives with a zero coefficient in that column are dropped, which + is the convention integrals.ts stores (only the non-zero span). + """ + out = [] + for shell in gto.basis.load(basis, symbol): + l = shell[0] + rows = shell[1:] + ncol = len(rows[0]) - 1 + for col in range(ncol): + pairs = [(r[0], r[col + 1]) for r in rows if r[col + 1] != 0.0] + if not pairs: + continue + out.append((l, tuple(a for a, _ in pairs), tuple(c for _, c in pairs))) + return out + + +def same(x, y, tol=1e-10): + return len(x) == len(y) and all( + abs(a - b) <= tol * max(1.0, abs(b)) for a, b in zip(x, y) + ) + + +def diffuse_functions(symbol: str): + """aug-cc-pVDZ minus cc-pVDZ, per l. Exactly one extra per l.""" + base = contracted_functions("cc-pvdz", symbol) + aug = contracted_functions("aug-cc-pvdz", symbol) + remaining = list(aug) + for f in base: + for i, g in enumerate(remaining): + if f[0] == g[0] and same(f[1], g[1]) and same(f[2], g[2]): + remaining.pop(i) + break + else: + raise RuntimeError( + f"{symbol}: cc-pVDZ function l={f[0]} alpha={f[1]} " + "has no counterpart in aug-cc-pVDZ" + ) + for l, alphas, _ in remaining: + if len(alphas) != 1: + raise RuntimeError( + f"{symbol}: aug diffuse l={l} is not a single primitive: {alphas}" + ) + return remaining + + +# ── naming ─────────────────────────────────────────────────────── + +def named_shells(basis: str, symbol: str): + """OrderedDict TS_CONST_NAME -> (l, alphas, coeffs).""" + prefix = TS_PREFIX[basis] + el = symbol.upper() + out = OrderedDict() + + if basis == "aug-cc-pvdz": + for l, alphas, coeffs in diffuse_functions(symbol): + out[f"{prefix}_{el}_DIFFUSE_{LSYM[l]}"] = (l, alphas, coeffs) + return out + + labels = LABELS[(basis, row_of(symbol))] + seen = {} + for l, alphas, coeffs in contracted_functions(basis, symbol): + i = seen.get(l, 0) + seen[l] = i + 1 + table = labels.get(l) + if table is None or i >= len(table): + raise RuntimeError( + f"{basis} {symbol}: unexpected function #{i} with l={l}; " + "the LABELS table needs extending" + ) + out[f"{prefix}_{el}_{table[i]}"] = (l, alphas, coeffs) + return out + + +def all_shells(elements): + """(basis, symbol) -> OrderedDict of named shells.""" + return { + (b, s): named_shells(b, s) + for b in BASES + for s in elements + } + + +# ── STO-3G shared-coefficient factoring ────────────────────────── + +def shared_sto3g_vectors(elements): + """label -> coefficient tuple, verified identical across elements. + + Returns only the labels where every element that has that shell + agrees to SHARED_TOL. Anything else stays inline. + """ + seen = {} + reject = set() + for s in elements: + for name, (_l, _a, c) in named_shells("sto-3g", s).items(): + label = name.split("_", 2)[2] + if label not in SHARED_STO3G: + continue + if label in seen: + if not same(c, seen[label], SHARED_TOL): + reject.add(label) + else: + seen[label] = c + return {k: v for k, v in seen.items() if k not in reject} + + +# ── TypeScript emission ────────────────────────────────────────── + +def fmt(x: float) -> str: + """Shortest round-tripping literal. Python's repr is exactly that, + and every form it produces (incl. 1.7016e-05) is valid TS.""" + return repr(float(x)) + + +def fmt_array(values, indent: str) -> str: + body = ", ".join(fmt(v) for v in values) + if len(body) <= 72: + return f"[{body}] as const" + lines, cur = [], [] + for v in values: + cur.append(fmt(v)) + if len(", ".join(cur)) > 64: + lines.append(", ".join(cur) + ",") + cur = [] + if cur: + lines.append(", ".join(cur) + ",") + inner = "\n".join(indent + " " + ln for ln in lines) + return "[\n" + inner + "\n" + indent + "] as const" + + +def emit_shared_block(elements) -> str: + shared = shared_sto3g_vectors(elements) + new = [(lbl, shared[lbl]) for lbl, (_n, exists) in SHARED_STO3G.items() + if lbl in shared and not exists] + if not new: + return "" + lines = [ + "// ── STO-3G shared contraction coefficients (M-shell, row 3) ──", + "// STO-3G's radial form is universal: within a shell type only the", + "// exponents change from atom to atom. integrals.ts already factors", + "// out STO3G_S_C / STO3G_L_2S_C / STO3G_L_2P_C; these are the 3s/3p", + "// (M-shell) counterparts, needed by Na-Ar. Verified identical across", + f"// all row-3 elements to {SHARED_TOL:g} relative.", + "", + ] + for label, vec in new: + name = SHARED_STO3G[label][0] + lines.append(f"const {name} = {fmt_array(vec, '')};") + lines.append("") + return "\n".join(lines) + + +def emit_element(basis: str, symbol: str, shared) -> str: + lines = [f"// ── {symbol} {basis} ──"] + for name, (_l, alphas, coeffs) in named_shells(basis, symbol).items(): + label = name.split("_", 2)[2] if basis != "aug-cc-pvdz" else None + cref = None + if basis == "sto-3g" and label in shared and same(coeffs, shared[label], SHARED_TOL): + cref = SHARED_STO3G[label][0] + lines.append(f"export const {name} = {{") + lines.append(f" alpha: {fmt_array(alphas, ' ')},") + lines.append(f" c: {cref if cref else fmt_array(coeffs, ' ')},") + lines.append("};") + lines.append("") + return "\n".join(lines) + + +def emit(elements, outdir: str): + os.makedirs(outdir, exist_ok=True) + shared = shared_sto3g_vectors(elements) + import pyscf + header = ( + f"// GENERATED by scripts/gen-basis-tables.py from PySCF {pyscf.__version__}.\n" + "// Source: pyscf.gto.basis.load. Do not hand-edit.\n" + "// Requires already in scope (integrals.ts defines them):\n" + "// STO3G_S_C, STO3G_L_2S_C, STO3G_L_2P_C\n" + ) + combined = [header, emit_shared_block(elements)] + for basis in BASES: + chunk = [header, ""] + if basis == "sto-3g": + chunk.append(emit_shared_block(elements)) + for s in elements: + chunk.append(emit_element(basis, s, shared)) + combined.append(emit_element(basis, s, shared)) + path = os.path.join(outdir, TS_PREFIX[basis].lower() + ".ts") + with open(path, "w") as fh: + fh.write("\n".join(chunk)) + with open(os.path.join(outdir, "all.ts"), "w") as fh: + fh.write("\n".join(combined)) + return outdir + + +# ── trust gate ─────────────────────────────────────────────────── + +def reldev(a, b): + """max_i |a_i - b_i| / max(|b_i|, tiny). None if lengths differ.""" + if len(a) != len(b): + return None + worst = 0.0 + for x, y in zip(a, b): + denom = abs(y) if abs(y) > 0 else 1.0 + worst = max(worst, abs(x - y) / denom) + return worst + + +def verify(dump_path: str, elements): + with open(dump_path) as fh: + ours = json.load(fh) + + rows = [] + for basis in BASES: + for s in elements: + gen = named_shells(basis, s) + present = [n for n in gen if n in ours] + if not present: + continue + worst = 0.0 + worst_at = "" + struct = [] + for name in present: + _l, alphas, coeffs = gen[name] + da = reldev(ours[name]["alpha"], alphas) + dc = reldev(ours[name]["c"], coeffs) + if da is None or dc is None: + struct.append( + f"{name}: len(ours)={len(ours[name]['alpha'])}/" + f"{len(ours[name]['c'])} vs pyscf={len(alphas)}/{len(coeffs)}" + ) + continue + for tag, d in (("alpha", da), ("c", dc)): + if d > worst: + worst, worst_at = d, f"{name}.{tag}" + missing = [n for n in gen if n not in ours] + rows.append({ + "basis": basis, "element": s, + "checked": len(present), "generated": len(gen), + "missing_in_repo": missing, + "struct_mismatch": struct, + "max_rel": worst, "at": worst_at, + }) + + width = max(len(r["at"]) for r in rows) if rows else 10 + print(f"{'basis':12s} {'el':3s} {'cells':>6s} {'max rel dev':>12s} " + f"{'worst entry':{width}s} note") + print("-" * (44 + width + 8)) + for r in rows: + note = "" + if r["struct_mismatch"]: + note = "STRUCT: " + "; ".join(r["struct_mismatch"]) + elif r["missing_in_repo"]: + note = f"{len(r['missing_in_repo'])} not in repo: " + \ + ",".join(n.split("_", 2)[2] for n in r["missing_in_repo"]) + print(f"{r['basis']:12s} {r['element']:3s} " + f"{r['checked']:3d}/{r['generated']:<2d} {r['max_rel']:12.3e} " + f"{r['at']:{width}s} {note}") + + gate = [("sto-3g", e) for e in ("C", "N", "O", "F")] + \ + [("cc-pvdz", e) for e in ("C", "N", "O", "F")] + \ + [("aug-cc-pvdz", "C")] + print("\nTRUST GATE (known-clean cells, bar = 1e-6 relative):") + failures = 0 + for basis, el in gate: + r = next((x for x in rows if x["basis"] == basis and x["element"] == el), None) + if r is None: + print(f" MISSING {basis:12s} {el}") + failures += 1 + continue + ok = (not r["struct_mismatch"]) and r["max_rel"] <= 1e-6 + failures += 0 if ok else 1 + print(f" {'PASS' if ok else 'FAIL'} {basis:12s} {el:2s} " + f"max rel dev {r['max_rel']:.3e} at {r['at'] or '-'}") + print("\nGATE 0.2:", "PASS" if failures == 0 + else f"FAIL — {failures} clean cell(s) disagree") + return rows, failures + + +# ── shell structure report ─────────────────────────────────────── + +def structure(elements): + for basis in BASES: + print(f"\n=== {basis} ===") + print(f"{'el':3s} {'s':>2s} {'p':>2s} {'d':>2s} shells (label:nprim)") + for s in elements: + gen = named_shells(basis, s) + cnt = {} + desc = [] + for name, (l, alphas, _c) in gen.items(): + cnt[l] = cnt.get(l, 0) + 1 + desc.append(f"{name.split('_', 2)[2]}:{len(alphas)}") + print(f"{s:3s} {cnt.get(0,0):2d} {cnt.get(1,0):2d} {cnt.get(2,0):2d} " + + " ".join(desc)) + + +# ── main ───────────────────────────────────────────────────────── + +def main(): + ap = argparse.ArgumentParser() + ap.add_argument("--out", metavar="DIR") + ap.add_argument("--verify", metavar="JSON") + ap.add_argument("--structure", action="store_true") + ap.add_argument("--elements", nargs="*", default=ELEMENTS) + args = ap.parse_args() + + elements = args.elements + rc = 0 + + if args.structure: + structure(elements) + + if args.verify: + print("\n=== trust gate: generator vs integrals.ts ===") + _rows, failures = verify(args.verify, elements) + rc = 1 if failures else 0 + + if args.out: + path = emit(elements, args.out) + print(f"\nwrote TypeScript constants to {path}/") + + if not (args.out or args.verify or args.structure): + ap.error("nothing to do: pass --out, --verify and/or --structure") + + return rc + + +if __name__ == "__main__": + sys.exit(main()) diff --git a/src/chemistry/integrals.ts b/src/chemistry/integrals.ts index b494285..67b9ba5 100644 --- a/src/chemistry/integrals.ts +++ b/src/chemistry/integrals.ts @@ -144,12 +144,18 @@ export const STO3G_O_2P = { // ── Hydrogen cc-pVDZ ───────────────────────────────────────── // 2 s-shells (one contracted, one uncontracted) + 1 p-shell. +// Corrected 2026-08-10 against pyscf.gto.basis.load("cc-pvdz", "H"). +// The exponents carried extra non-canonical digits (13.0107010 / +// 1.9622572 / 0.4445298 / 0.1219496) against Dunning's published +// 13.01 / 1.962 / 0.4446 / 0.122. Coefficients were already correct. +// Small (~2-5 uHa) but hydrogen is in nearly every validation molecule, +// so it put a permanent floor under every cc-pVDZ vs-PySCF comparison. export const CCPVDZ_H_1S = { - alpha: [13.0107010, 1.9622572, 0.4445298] as const, + alpha: [13.0100000, 1.9620000, 0.4446000] as const, c: [0.0196850, 0.1379770, 0.4781480] as const, }; export const CCPVDZ_H_2S = { - alpha: [0.1219496] as const, + alpha: [0.1220000] as const, c: [1.0] as const, }; export const CCPVDZ_H_2P = { @@ -175,11 +181,13 @@ export const CCPVDZ_HE_2P = { // ── Helium aug-cc-pVDZ diffuse ─────────────────────────────── // One diffuse s + one diffuse p (EMSL). +// Corrected 2026-08-10 against pyscf.gto.basis.load("aug-cc-pvdz", "He") +// — was 0.0713 / 0.3300. export const AUG_CCPVDZ_HE_DIFFUSE_S = { - alpha: [0.0713000] as const, c: [1.0] as const, + alpha: [0.0725500] as const, c: [1.0] as const, }; export const AUG_CCPVDZ_HE_DIFFUSE_P = { - alpha: [0.3300000] as const, c: [1.0] as const, + alpha: [0.2473000] as const, c: [1.0] as const, }; // ── Oxygen cc-pVDZ ──────────────────────────────────────────── @@ -250,12 +258,14 @@ export const AUG_CCPVDZ_H_DIFFUSE_P = { alpha: [0.1410000] as const, c: [1.0] as const, }; +// Corrected 2026-08-10 against pyscf.gto.basis.load("aug-cc-pvdz", "O") +// — was 0.0845800 / 0.0856000. The d exponent was already correct. export const AUG_CCPVDZ_O_DIFFUSE_S = { - alpha: [0.0845800] as const, + alpha: [0.0789600] as const, c: [1.0] as const, }; export const AUG_CCPVDZ_O_DIFFUSE_P = { - alpha: [0.0856000] as const, + alpha: [0.0685600] as const, c: [1.0] as const, }; export const AUG_CCPVDZ_O_DIFFUSE_D = { @@ -266,21 +276,28 @@ export const AUG_CCPVDZ_O_DIFFUSE_D = { // aug-cc-pVDZ diffuse functions for Li, Be, C, N, F. // One diffuse primitive per angular momentum class (s, p, d). // EMSL Basis Set Exchange, retrieved 2026-05-18. -export const AUG_CCPVDZ_LI_DIFFUSE_S = { alpha: [0.0072930] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_LI_DIFFUSE_P = { alpha: [0.0074000] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_LI_DIFFUSE_D = { alpha: [0.0950000] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_BE_DIFFUSE_S = { alpha: [0.0207000] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_BE_DIFFUSE_P = { alpha: [0.0142000] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_BE_DIFFUSE_D = { alpha: [0.0722000] as const, c: [1.0] as const }; +// Corrected 2026-08-10 against pyscf.gto.basis.load("aug-cc-pvdz", ...). +// Of the eight elements, only H and C were right; He/Li/Be/N/O/F all +// carried wrong diffuse exponents. Previous values, for the record: +// Li S 0.0072930 P 0.0074000 D 0.0950000 +// Be S 0.0207000 P 0.0142000 D 0.0722000 +// N S 0.0576000 P 0.0491000 (D was correct) +// F S 0.1076000 P 0.0832000 D 0.5000000 +export const AUG_CCPVDZ_LI_DIFFUSE_S = { alpha: [0.0086400] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_LI_DIFFUSE_P = { alpha: [0.0057900] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_LI_DIFFUSE_D = { alpha: [0.0725000] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_BE_DIFFUSE_S = { alpha: [0.0187700] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_BE_DIFFUSE_P = { alpha: [0.0085000] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_BE_DIFFUSE_D = { alpha: [0.0740000] as const, c: [1.0] as const }; export const AUG_CCPVDZ_C_DIFFUSE_S = { alpha: [0.0469000] as const, c: [1.0] as const }; export const AUG_CCPVDZ_C_DIFFUSE_P = { alpha: [0.0404100] as const, c: [1.0] as const }; export const AUG_CCPVDZ_C_DIFFUSE_D = { alpha: [0.1510000] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_N_DIFFUSE_S = { alpha: [0.0576000] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_N_DIFFUSE_P = { alpha: [0.0491000] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_N_DIFFUSE_S = { alpha: [0.0612400] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_N_DIFFUSE_P = { alpha: [0.0561100] as const, c: [1.0] as const }; export const AUG_CCPVDZ_N_DIFFUSE_D = { alpha: [0.2300000] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_F_DIFFUSE_S = { alpha: [0.1076000] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_F_DIFFUSE_P = { alpha: [0.0832000] as const, c: [1.0] as const }; -export const AUG_CCPVDZ_F_DIFFUSE_D = { alpha: [0.5000000] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_F_DIFFUSE_S = { alpha: [0.0986300] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_F_DIFFUSE_P = { alpha: [0.0850200] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_F_DIFFUSE_D = { alpha: [0.4640000] as const, c: [1.0] as const }; /** STO-3G 2s contraction for lithium — the *s component* of the L-shell. * Coefficients can be negative to produce the 2s radial node and to @@ -332,13 +349,18 @@ export const CCPVDZ_LI_2S = { -0.048902, -0.096009, -0.136380, 0.575102, ] as const, }; +// Corrected 2026-08-10 against pyscf.gto.basis.load("cc-pvdz", "Li"). +// The uncontracted s exponent read 0.0285000 (digit transposition of +// 0.02805) and the 2p contraction coefficients were wholly wrong — +// together worth 1.29 mHa on LiH, which is 2.6x the repo's own +// <= 0.5 mHa HF-vs-PySCF gate. See tests/chemistry/elements/. export const CCPVDZ_LI_2S_P = { - alpha: [0.0285000] as const, + alpha: [0.0280500] as const, c: [1.0] as const, }; export const CCPVDZ_LI_2P = { alpha: [1.5340, 0.2749, 0.073620] as const, - c: [0.022270, 0.111705, 0.240315] as const, + c: [0.022784, 0.139107, 0.500375] as const, }; export const CCPVDZ_LI_2P_P = { alpha: [0.024030] as const, @@ -356,9 +378,12 @@ export const CCPVDZ_BE_1S = { 2940.0000, 441.2000, 100.5000, 28.4300, 9.1690, 3.1960, 1.1590, 0.1811, ] as const, + // Corrected 2026-08-10 against pyscf.gto.basis.load("cc-pvdz", "Be"): + // the tail was transcribed at reduced precision (0.4514 vs 0.451469, + // 0.2950 vs 0.295074, 0.012580 vs 0.012587). c: [ - 0.000680, 0.005236, 0.026600, 0.099990, - 0.269700, 0.451400, 0.295000, 0.012580, + 0.000680, 0.005236, 0.026606, 0.099993, + 0.269702, 0.451469, 0.295074, 0.012587, ] as const, }; export const CCPVDZ_BE_2S = { @@ -366,9 +391,13 @@ export const CCPVDZ_BE_2S = { 2940.0000, 441.2000, 100.5000, 28.4300, 9.1690, 3.1960, 1.1590, 0.1811, ] as const, + // Corrected 2026-08-10 against PySCF. Elements 4-8 were not merely + // rounded but substantively different (-0.05403 vs -0.05328, + // -0.1133 vs -0.120723, -0.1462 vs -0.133435, 0.5392 vs 0.530767), + // i.e. sourced from a different/incorrect table. c: [ - -0.000123, -0.000966, -0.004831, -0.019330, - -0.054030, -0.113300, -0.146200, 0.539200, + -0.000123, -0.000966, -0.004831, -0.019314, + -0.053280, -0.120723, -0.133435, 0.530767, ] as const, }; export const CCPVDZ_BE_2S_P = { @@ -377,7 +406,9 @@ export const CCPVDZ_BE_2S_P = { }; export const CCPVDZ_BE_2P = { alpha: [3.6190, 0.7110, 0.1951] as const, - c: [0.029110, 0.169650, 0.487810] as const, + // Corrected 2026-08-10 against PySCF (was 0.169650 / 0.487810 — the + // last coefficient was off by 5%). + c: [0.029111, 0.169365, 0.513458] as const, }; export const CCPVDZ_BE_2P_P = { alpha: [0.060180] as const, @@ -434,9 +465,10 @@ export const CCPVDZ_N_1S = { 9046.0000, 1357.0000, 309.3000, 87.7300, 28.5600, 10.2100, 3.8380, 0.7466, ] as const, + // Corrected 2026-08-10 against PySCF (was 0.278722 / 0.448581). c: [ 0.000700, 0.005389, 0.027406, 0.103207, - 0.278722, 0.448581, 0.278238, 0.015440, + 0.278723, 0.448540, 0.278238, 0.015440, ] as const, }; export const CCPVDZ_N_2S = { From 1bd245fd9c144b451bf3ffcc1b85587611b1ec9b Mon Sep 17 00:00:00 2001 From: Ahmet Gunaydin Date: Mon, 10 Aug 2026 09:36:26 +0700 Subject: [PATCH 04/13] =?UTF-8?q?feat:=20add=20Boron=20and=20Neon=20?= =?UTF-8?q?=E2=80=94=20row=202=20complete=20(10=20elements)?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Phase 1a of docs/RUN-PLAN-24H-ELEMENTS.md: the pipeline proof. B and Ne are structurally identical to C/N/O/F ([3s,2p,1d] in cc-pVDZ, 1s + an L-shell in STO-3G), so they build through the existing `heavyShells` and `heavyDiffuse` helpers unchanged. A failure here would have meant a broken harness rather than a hard element. There was no failure: all 12 new cells passed on the first run, and the element agreement test is now 61/61 green across 10 elements x 3 bases x 2 d-conventions. Basis constants are generated by scripts/gen-basis-tables.py from PySCF, not transcribed. New reference molecules: BH3 (planar D3h, 1.19 A) and atomic Ne. Per-element data landed in every table that needed it. Two of these are string-keyed with silent fallbacks and would have produced wrong physics with no error at all: geometry.ts COVALENT_RADIUS_ANGSTROM (?? 1.0 -> wrong bond detection) geometry.ts ATOMIC_MASS (?? 0 -> zero mass, garbage centre-of-mass and rotational constants) Sources: NIST most-abundant-isotope masses, Pauling electronegativity (null for Ne, as for He), Pyykko-Atsumi covalent radii (matching the file's stated source, not Cordero), Bondi/Mantina vdW, Grimme 2006 D2 C6 and R_r. Widening AtomSymbol makes every Record table a compile error until filled — the desired loud failure. Five separate copies of ATOMIC_MASS and four independent element lists exist across the repo; all were updated, but the duplication is a standing drift hazard worth consolidating later. tests/chemistry/xyz.test.ts used Boron as its "unsupported element" example, which this commit invalidates — it now uses Krypton, well away from the basis tables, with B and Ne round-trip tests added alongside. Full suite: 140 files, 1017 passed, 1 skipped. Co-Authored-By: Claude Opus 5 (1M context) --- scripts/run-element-reference.py | 5 + src/chemistry/atoms.ts | 37 +- src/chemistry/dispersion-d2.ts | 8 + src/chemistry/geometry.ts | 6 +- src/chemistry/integrals.ts | 112 ++++++ src/chemistry/periodic-table.ts | 14 +- src/chemistry/raman.ts | 2 + src/chemistry/thermochemistry.ts | 2 + src/chemistry/vibrations.ts | 6 +- src/chemistry/xyz.ts | 11 +- src/molecule/import-formats.ts | 4 +- tests/chemistry/elements/pyscf-reference.json | 372 ++++++++++++++++++ tests/chemistry/periodic-table.test.ts | 12 +- tests/chemistry/xyz.test.ts | 35 +- 14 files changed, 604 insertions(+), 22 deletions(-) diff --git a/scripts/run-element-reference.py b/scripts/run-element-reference.py index e60e47a..c855ef5 100644 --- a/scripts/run-element-reference.py +++ b/scripts/run-element-reference.py @@ -57,8 +57,13 @@ def pyramidal(bond, angle_deg): "Li": ("LiH", [("Li", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 1.5949))]), "Be": ("BeH2", [("Be", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 1.3264)), ("H", (0.0, 0.0, -1.3264))]), + "B": ("BH3", [("B", (0.0, 0.0, 0.0))] + + [("H", (1.1900 * math.cos(math.radians(90 + 120 * k)), + 1.1900 * math.sin(math.radians(90 + 120 * k)), 0.0)) + for k in range(3)]), "C": ("CH4", [("C", (0.0, 0.0, 0.0))] + [("H", p) for p in tetrahedral(1.0870)]), + "Ne": ("Ne", [("Ne", (0.0, 0.0, 0.0))]), "N": ("NH3", [("N", (0.0, 0.0, 0.0))] + [("H", p) for p in pyramidal(1.0124, 106.67)]), "O": ("H2O", [("O", (0.0, 0.0, 0.0))] + diff --git a/src/chemistry/atoms.ts b/src/chemistry/atoms.ts index 390f3a8..29955f5 100644 --- a/src/chemistry/atoms.ts +++ b/src/chemistry/atoms.ts @@ -10,6 +10,14 @@ import { STO3G_H_1S, STO3G_HE_1S, STO3G_LI_1S, STO3G_LI_2S, STO3G_BE_1S, STO3G_BE_2S, STO3G_BE_2P, + STO3G_B_1S, STO3G_B_2S, STO3G_B_2P, + STO3G_NE_1S, STO3G_NE_2S, STO3G_NE_2P, + CCPVDZ_B_1S, CCPVDZ_B_2S, CCPVDZ_B_2S_P, + CCPVDZ_B_2P, CCPVDZ_B_2P_P, CCPVDZ_B_3D, + CCPVDZ_NE_1S, CCPVDZ_NE_2S, CCPVDZ_NE_2S_P, + CCPVDZ_NE_2P, CCPVDZ_NE_2P_P, CCPVDZ_NE_3D, + AUG_CCPVDZ_B_DIFFUSE_S, AUG_CCPVDZ_B_DIFFUSE_P, AUG_CCPVDZ_B_DIFFUSE_D, + AUG_CCPVDZ_NE_DIFFUSE_S, AUG_CCPVDZ_NE_DIFFUSE_P, AUG_CCPVDZ_NE_DIFFUSE_D, STO3G_C_1S, STO3G_C_2S, STO3G_C_2P, STO3G_N_1S, STO3G_N_2S, STO3G_N_2P, STO3G_O_1S, STO3G_O_2S, STO3G_O_2P, @@ -40,7 +48,8 @@ import { import { type CGShell, makeCGShell } from "./integrals-cg.js"; import { type Nucleus } from "./cg-molecular.js"; -export type AtomSymbol = "H" | "He" | "Li" | "Be" | "C" | "N" | "O" | "F"; +export type AtomSymbol = + | "H" | "He" | "Li" | "Be" | "B" | "C" | "N" | "O" | "F" | "Ne"; export type BasisName = "sto-3g" | "cc-pvdz" | "aug-cc-pvdz"; const ANGSTROM_TO_BOHR = 1 / 0.529177210903; @@ -60,12 +69,12 @@ export interface Atom { /** Atomic number for each supported atom. */ export const Z_FOR: Readonly> = { - H: 1, He: 2, Li: 3, Be: 4, C: 6, N: 7, O: 8, F: 9, + H: 1, He: 2, Li: 3, Be: 4, B: 5, C: 6, N: 7, O: 8, F: 9, Ne: 10, }; /** Number of electrons in the neutral atom. */ export const N_ELECTRONS_FOR: Readonly> = { - H: 1, He: 2, Li: 3, Be: 4, C: 6, N: 7, O: 8, F: 9, + H: 1, He: 2, Li: 3, Be: 4, B: 5, C: 6, N: 7, O: 8, F: 9, Ne: 10, }; /** @@ -79,7 +88,7 @@ export const N_ELECTRONS_FOR: Readonly> = { * comparative calculation cares about). */ export const FROZEN_CORE_FOR: Readonly> = { - H: 0, He: 0, Li: 1, Be: 1, C: 1, N: 1, O: 1, F: 1, + H: 0, He: 0, Li: 1, Be: 1, B: 1, C: 1, N: 1, O: 1, F: 1, Ne: 1, }; /** Default frozen-core count for a molecule (sum of per-atom 1s cores). */ @@ -132,6 +141,14 @@ export function atomShells( makeCGShell(STO3G_BE_2P, pos_bohr, [0, 1, 0], "Be:2p_y"), makeCGShell(STO3G_BE_2P, pos_bohr, [0, 0, 1], "Be:2p_z"), ]; + case "B": + return [ + makeCGShell(STO3G_B_1S, pos_bohr, [0, 0, 0], "B:1s"), + makeCGShell(STO3G_B_2S, pos_bohr, [0, 0, 0], "B:2s"), + makeCGShell(STO3G_B_2P, pos_bohr, [1, 0, 0], "B:2p_x"), + makeCGShell(STO3G_B_2P, pos_bohr, [0, 1, 0], "B:2p_y"), + makeCGShell(STO3G_B_2P, pos_bohr, [0, 0, 1], "B:2p_z"), + ]; case "C": return [ makeCGShell(STO3G_C_1S, pos_bohr, [0, 0, 0], "C:1s"), @@ -164,6 +181,14 @@ export function atomShells( makeCGShell(STO3G_F_2P, pos_bohr, [0, 1, 0], "F:2p_y"), makeCGShell(STO3G_F_2P, pos_bohr, [0, 0, 1], "F:2p_z"), ]; + case "Ne": + return [ + makeCGShell(STO3G_NE_1S, pos_bohr, [0, 0, 0], "Ne:1s"), + makeCGShell(STO3G_NE_2S, pos_bohr, [0, 0, 0], "Ne:2s"), + makeCGShell(STO3G_NE_2P, pos_bohr, [1, 0, 0], "Ne:2p_x"), + makeCGShell(STO3G_NE_2P, pos_bohr, [0, 1, 0], "Ne:2p_y"), + makeCGShell(STO3G_NE_2P, pos_bohr, [0, 0, 1], "Ne:2p_z"), + ]; } } @@ -208,10 +233,12 @@ function atomShellsCcPvdz(symbol: AtomSymbol, pos: readonly [number, number, num ]; case "Li": return heavyShells("Li", CCPVDZ_LI_1S, CCPVDZ_LI_2S, CCPVDZ_LI_2S_P, CCPVDZ_LI_2P, CCPVDZ_LI_2P_P, CCPVDZ_LI_3D); case "Be": return heavyShells("Be", CCPVDZ_BE_1S, CCPVDZ_BE_2S, CCPVDZ_BE_2S_P, CCPVDZ_BE_2P, CCPVDZ_BE_2P_P, CCPVDZ_BE_3D); + case "B": return heavyShells("B", CCPVDZ_B_1S, CCPVDZ_B_2S, CCPVDZ_B_2S_P, CCPVDZ_B_2P, CCPVDZ_B_2P_P, CCPVDZ_B_3D); case "C": return heavyShells("C", CCPVDZ_C_1S, CCPVDZ_C_2S, CCPVDZ_C_2S_P, CCPVDZ_C_2P, CCPVDZ_C_2P_P, CCPVDZ_C_3D); case "N": return heavyShells("N", CCPVDZ_N_1S, CCPVDZ_N_2S, CCPVDZ_N_2S_P, CCPVDZ_N_2P, CCPVDZ_N_2P_P, CCPVDZ_N_3D); case "O": return heavyShells("O", CCPVDZ_O_1S, CCPVDZ_O_2S, CCPVDZ_O_2S_P, CCPVDZ_O_2P, CCPVDZ_O_2P_P, CCPVDZ_O_3D); case "F": return heavyShells("F", CCPVDZ_F_1S, CCPVDZ_F_2S, CCPVDZ_F_2S_P, CCPVDZ_F_2P, CCPVDZ_F_2P_P, CCPVDZ_F_3D); + case "Ne": return heavyShells("Ne", CCPVDZ_NE_1S, CCPVDZ_NE_2S, CCPVDZ_NE_2S_P, CCPVDZ_NE_2P, CCPVDZ_NE_2P_P, CCPVDZ_NE_3D); case "He": return [ makeCGShell(CCPVDZ_HE_1S, pos, [0, 0, 0], "He:1s"), @@ -254,10 +281,12 @@ function atomShellsAugDiffuse(symbol: AtomSymbol, pos: readonly [number, number, ]; case "Li": return heavyDiffuse("Li", AUG_CCPVDZ_LI_DIFFUSE_S, AUG_CCPVDZ_LI_DIFFUSE_P, AUG_CCPVDZ_LI_DIFFUSE_D); case "Be": return heavyDiffuse("Be", AUG_CCPVDZ_BE_DIFFUSE_S, AUG_CCPVDZ_BE_DIFFUSE_P, AUG_CCPVDZ_BE_DIFFUSE_D); + case "B": return heavyDiffuse("B", AUG_CCPVDZ_B_DIFFUSE_S, AUG_CCPVDZ_B_DIFFUSE_P, AUG_CCPVDZ_B_DIFFUSE_D); case "C": return heavyDiffuse("C", AUG_CCPVDZ_C_DIFFUSE_S, AUG_CCPVDZ_C_DIFFUSE_P, AUG_CCPVDZ_C_DIFFUSE_D); case "N": return heavyDiffuse("N", AUG_CCPVDZ_N_DIFFUSE_S, AUG_CCPVDZ_N_DIFFUSE_P, AUG_CCPVDZ_N_DIFFUSE_D); case "O": return heavyDiffuse("O", AUG_CCPVDZ_O_DIFFUSE_S, AUG_CCPVDZ_O_DIFFUSE_P, AUG_CCPVDZ_O_DIFFUSE_D); case "F": return heavyDiffuse("F", AUG_CCPVDZ_F_DIFFUSE_S, AUG_CCPVDZ_F_DIFFUSE_P, AUG_CCPVDZ_F_DIFFUSE_D); + case "Ne": return heavyDiffuse("Ne", AUG_CCPVDZ_NE_DIFFUSE_S, AUG_CCPVDZ_NE_DIFFUSE_P, AUG_CCPVDZ_NE_DIFFUSE_D); case "He": return [ makeCGShell(AUG_CCPVDZ_HE_DIFFUSE_S, pos, [0, 0, 0], "He:aug-s"), diff --git a/src/chemistry/dispersion-d2.ts b/src/chemistry/dispersion-d2.ts index 747f810..d5f3ec0 100644 --- a/src/chemistry/dispersion-d2.ts +++ b/src/chemistry/dispersion-d2.ts @@ -43,10 +43,12 @@ const C6_J_NM6_PER_MOL: Readonly> = { He: 0.08, Li: 1.61, Be: 1.61, + B: 3.13, C: 1.75, N: 1.23, O: 0.70, F: 0.75, + Ne: 0.63, }; /** Grimme 2006 Table I: per-element R_R van-der-Waals radius in Å. */ @@ -55,10 +57,12 @@ const R_R_ANGSTROM: Readonly> = { He: 1.012, Li: 0.825, Be: 1.408, + B: 1.485, C: 1.452, N: 1.397, O: 1.342, F: 1.287, + Ne: 1.243, }; /** Per-element C6 in Hartree·Bohr^6 (atomic units). */ @@ -67,10 +71,12 @@ const C6_AU: Readonly> = { He: C6_J_NM6_PER_MOL.He * J_NM6_PER_MOL_TO_HARTREE_BOHR6, Li: C6_J_NM6_PER_MOL.Li * J_NM6_PER_MOL_TO_HARTREE_BOHR6, Be: C6_J_NM6_PER_MOL.Be * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + B: C6_J_NM6_PER_MOL.B * J_NM6_PER_MOL_TO_HARTREE_BOHR6, C: C6_J_NM6_PER_MOL.C * J_NM6_PER_MOL_TO_HARTREE_BOHR6, N: C6_J_NM6_PER_MOL.N * J_NM6_PER_MOL_TO_HARTREE_BOHR6, O: C6_J_NM6_PER_MOL.O * J_NM6_PER_MOL_TO_HARTREE_BOHR6, F: C6_J_NM6_PER_MOL.F * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + Ne: C6_J_NM6_PER_MOL.Ne * J_NM6_PER_MOL_TO_HARTREE_BOHR6, }; /** Per-element R_R in Bohr. */ @@ -79,10 +85,12 @@ const R_R_BOHR: Readonly> = { He: R_R_ANGSTROM.He * ANGSTROM_TO_BOHR, Li: R_R_ANGSTROM.Li * ANGSTROM_TO_BOHR, Be: R_R_ANGSTROM.Be * ANGSTROM_TO_BOHR, + B: R_R_ANGSTROM.B * ANGSTROM_TO_BOHR, C: R_R_ANGSTROM.C * ANGSTROM_TO_BOHR, N: R_R_ANGSTROM.N * ANGSTROM_TO_BOHR, O: R_R_ANGSTROM.O * ANGSTROM_TO_BOHR, F: R_R_ANGSTROM.F * ANGSTROM_TO_BOHR, + Ne: R_R_ANGSTROM.Ne * ANGSTROM_TO_BOHR, }; /** Damping steepness — universal in D2. */ diff --git a/src/chemistry/geometry.ts b/src/chemistry/geometry.ts index f000c6f..01469ce 100644 --- a/src/chemistry/geometry.ts +++ b/src/chemistry/geometry.ts @@ -244,8 +244,8 @@ export function dihedralAngle( * single-bond detection in organic chemistry. */ const COVALENT_RADIUS_ANGSTROM: Readonly> = { H: 0.32, He: 0.46, - Li: 1.33, Be: 1.02, - C: 0.75, N: 0.71, O: 0.63, F: 0.64, + Li: 1.33, Be: 1.02, B: 0.85, + C: 0.75, N: 0.71, O: 0.63, F: 0.64, Ne: 0.67, }; export interface Bond { @@ -259,8 +259,10 @@ export interface Bond { const ATOMIC_MASS: Readonly> = { H: 1.00782503207, He: 4.002603254, Li: 7.0160034366, Be: 9.012183065, + B: 11.00930536, C: 12.0, N: 14.0030740048, O: 15.99491461956, F: 18.998403163, + Ne: 19.9924401762, }; /** Total mass of a molecule in amu. */ diff --git a/src/chemistry/integrals.ts b/src/chemistry/integrals.ts index 67b9ba5..8d15471 100644 --- a/src/chemistry/integrals.ts +++ b/src/chemistry/integrals.ts @@ -289,6 +289,12 @@ export const AUG_CCPVDZ_LI_DIFFUSE_D = { alpha: [0.0725000] as const, c: [1.0] a export const AUG_CCPVDZ_BE_DIFFUSE_S = { alpha: [0.0187700] as const, c: [1.0] as const }; export const AUG_CCPVDZ_BE_DIFFUSE_P = { alpha: [0.0085000] as const, c: [1.0] as const }; export const AUG_CCPVDZ_BE_DIFFUSE_D = { alpha: [0.0740000] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_B_DIFFUSE_S = { alpha: [0.0310500] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_B_DIFFUSE_P = { alpha: [0.0237800] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_B_DIFFUSE_D = { alpha: [0.0904000] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_NE_DIFFUSE_S = { alpha: [0.1230000] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_NE_DIFFUSE_P = { alpha: [0.1064000] as const, c: [1.0] as const }; +export const AUG_CCPVDZ_NE_DIFFUSE_D = { alpha: [0.6310000] as const, c: [1.0] as const }; export const AUG_CCPVDZ_C_DIFFUSE_S = { alpha: [0.0469000] as const, c: [1.0] as const }; export const AUG_CCPVDZ_C_DIFFUSE_P = { alpha: [0.0404100] as const, c: [1.0] as const }; export const AUG_CCPVDZ_C_DIFFUSE_D = { alpha: [0.1510000] as const, c: [1.0] as const }; @@ -327,6 +333,35 @@ export const STO3G_F_2P = { c: [0.15591627, 0.60768372, 0.39195739] as const, }; +// ── STO-3G for Boron and Neon (2026-08-10) ─────────────────── +// Generated by scripts/gen-basis-tables.py from PySCF, completing +// row 2. Structurally identical to C/N/O/F: 1s + an L-shell whose +// s and p components share exponents. +export const STO3G_B_1S = { + alpha: [48.7911130, 8.8873622, 2.4052670] as const, + c: STO3G_S_C, +}; +export const STO3G_B_2S = { + alpha: [2.2369561, 0.5198205, 0.1690618] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_B_2P = { + alpha: [2.2369561, 0.5198205, 0.1690618] as const, + c: STO3G_L_2P_C, +}; +export const STO3G_NE_1S = { + alpha: [207.0156100, 37.7081510, 10.2052970] as const, + c: STO3G_S_C, +}; +export const STO3G_NE_2S = { + alpha: [8.2463151, 1.9162662, 0.6232293] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_NE_2P = { + alpha: [8.2463151, 1.9162662, 0.6232293] as const, + c: STO3G_L_2P_C, +}; + // ── cc-pVDZ for Lithium ────────────────────────────────────── // EMSL ccPVDZ for Li (9s,4p,1d) → [3s,2p,1d]. Retrieved 2026-05-18. export const CCPVDZ_LI_1S = { @@ -537,6 +572,83 @@ export const CCPVDZ_F_3D = { c: [1.0] as const, }; +// ── cc-pVDZ for Boron and Neon (2026-08-10) ────────────────── +// Generated by scripts/gen-basis-tables.py from PySCF. Same +// [3s,2p,1d] shape as C/N/O/F, so `heavyShells` builds them +// unchanged: 15 Cartesian / 14 spherical functions per atom. +export const CCPVDZ_B_1S = { + alpha: [ + 4570.0000, 685.9000, 156.5000, 44.4700, + 14.4800, 5.1310, 1.8980, 0.3329, + ] as const, + c: [ + 0.000696, 0.005353, 0.027134, 0.101380, + 0.272055, 0.448403, 0.290123, 0.014322, + ] as const, +}; +export const CCPVDZ_B_2S = { + alpha: [ + 4570.0000, 685.9000, 156.5000, 44.4700, + 14.4800, 5.1310, 1.8980, 0.3329, + ] as const, + c: [ + -0.000139, -0.001097, -0.005444, -0.021916, + -0.059751, -0.138732, -0.131482, 0.539526, + ] as const, +}; +export const CCPVDZ_B_2S_P = { + alpha: [0.1043000] as const, + c: [1.0] as const, +}; +export const CCPVDZ_B_2P = { + alpha: [6.0010, 1.2410, 0.3364] as const, + c: [0.035481, 0.198072, 0.505230] as const, +}; +export const CCPVDZ_B_2P_P = { + alpha: [0.0953800] as const, + c: [1.0] as const, +}; +export const CCPVDZ_B_3D = { + alpha: [0.3430] as const, + c: [1.0] as const, +}; +export const CCPVDZ_NE_1S = { + alpha: [ + 17880.0000, 2683.0000, 611.5000, 173.5000, + 56.6400, 20.4200, 7.8100, 1.6530, + ] as const, + c: [ + 0.000738, 0.005677, 0.028883, 0.108540, + 0.290907, 0.448324, 0.258026, 0.015063, + ] as const, +}; +export const CCPVDZ_NE_2S = { + alpha: [ + 17880.0000, 2683.0000, 611.5000, 173.5000, + 56.6400, 20.4200, 7.8100, 1.6530, + ] as const, + c: [ + -0.000172, -0.001357, -0.006737, -0.027663, + -0.076208, -0.175227, -0.107038, 0.567050, + ] as const, +}; +export const CCPVDZ_NE_2S_P = { + alpha: [0.4869000] as const, + c: [1.0] as const, +}; +export const CCPVDZ_NE_2P = { + alpha: [28.3900, 6.2700, 1.6950] as const, + c: [0.046087, 0.240181, 0.508744] as const, +}; +export const CCPVDZ_NE_2P_P = { + alpha: [0.4317000] as const, + c: [1.0] as const, +}; +export const CCPVDZ_NE_3D = { + alpha: [2.2020] as const, + c: [1.0] as const, +}; + // ── Generic shell type + multi-shell integrals ─────────────── // // A Shell is a contracted s-Gaussian: a center plus a list of diff --git a/src/chemistry/periodic-table.ts b/src/chemistry/periodic-table.ts index f851754..1a146eb 100644 --- a/src/chemistry/periodic-table.ts +++ b/src/chemistry/periodic-table.ts @@ -1,6 +1,6 @@ // ───────────────────────────────────────────────────────────── // periodic-table.ts — atomic property database for the elements -// supported by webgpu-q (H, He, Li, Be, C, N, O, F). +// supported by webgpu-q (H, He, Li, Be, B, C, N, O, F, Ne). // // Data sources: // - Atomic number, isotope mass: NIST 2024. @@ -58,6 +58,12 @@ export const PERIODIC_TABLE: Readonly> = { covalentRadius: 1.02, vdwRadius: 1.53, period: 2, group: 2, family: "alkaline-earth", }, + B: { + symbol: "B", name: "Boron", atomicNumber: 5, + mass: 11.00930536, electronegativity: 2.04, + covalentRadius: 0.85, vdwRadius: 1.92, + period: 2, group: 13, family: "main", + }, C: { symbol: "C", name: "Carbon", atomicNumber: 6, mass: 12.0, electronegativity: 2.55, @@ -82,6 +88,12 @@ export const PERIODIC_TABLE: Readonly> = { covalentRadius: 0.64, vdwRadius: 1.47, period: 2, group: 17, family: "halogen", }, + Ne: { + symbol: "Ne", name: "Neon", atomicNumber: 10, + mass: 19.9924401762, electronegativity: null, + covalentRadius: 0.67, vdwRadius: 1.54, + period: 2, group: 18, family: "noble", + }, }; /** Look up element data. Throws on unknown symbol. */ diff --git a/src/chemistry/raman.ts b/src/chemistry/raman.ts index ef0f851..61a4855 100644 --- a/src/chemistry/raman.ts +++ b/src/chemistry/raman.ts @@ -56,10 +56,12 @@ const ATOMIC_MASS: Readonly> = { He: 4.002603254, Li: 7.0160034366, Be: 9.012183065, + B: 11.00930536, C: 12.0, N: 14.0030740048, O: 15.99491461956, F: 18.998403163, + Ne: 19.9924401762, }; export type EnergyMethod = "hf" | FunctionalKind; diff --git a/src/chemistry/thermochemistry.ts b/src/chemistry/thermochemistry.ts index ea92ecb..02b806c 100644 --- a/src/chemistry/thermochemistry.ts +++ b/src/chemistry/thermochemistry.ts @@ -77,10 +77,12 @@ const ATOMIC_MASS: Readonly> = { He: 4.002603254, Li: 7.0160034366, Be: 9.012183065, + B: 11.00930536, C: 12.0, N: 14.0030740048, O: 15.99491461956, F: 18.998403163, + Ne: 19.9924401762, }; export interface ThermochemistryOpts { diff --git a/src/chemistry/vibrations.ts b/src/chemistry/vibrations.ts index 1d20429..ca4b0e5 100644 --- a/src/chemistry/vibrations.ts +++ b/src/chemistry/vibrations.ts @@ -47,17 +47,19 @@ const HA_BOHR2_AMU_TO_CM_INV = 5140.4865; const E2_PER_AMU_TO_KM_PER_MOL = 974.864; /** Atomic masses for the supported elements. Isotope-pure values - * (^1H, ^7Li, ^9Be, ^12C, ^14N, ^16O) — the standard choice for - * theoretical reference frequencies. */ + * (^1H, ^7Li, ^9Be, ^11B, ^12C, ^14N, ^16O, ^20Ne) — the standard + * choice for theoretical reference frequencies. */ const ATOMIC_MASS: Readonly> = { H: 1.00782503207, He: 4.002603254, Li: 7.0160034366, Be: 9.012183065, + B: 11.00930536, C: 12.0, N: 14.0030740048, O: 15.99491461956, F: 18.998403163, + Ne: 19.9924401762, }; export type EnergyMethod = "hf" | FunctionalKind; diff --git a/src/chemistry/xyz.ts b/src/chemistry/xyz.ts index 3b2364e..1ebf4ca 100644 --- a/src/chemistry/xyz.ts +++ b/src/chemistry/xyz.ts @@ -32,9 +32,12 @@ import { type Atom, type AtomSymbol } from "./atoms.js"; // molecules.ts ships a He molecule — omitting it here made toXYZ → parseXYZ fail // to round-trip a molecule this library provides. const SYMBOL_BY_Z: Record = { - 1: "H", 2: "He", 3: "Li", 4: "Be", 6: "C", 7: "N", 8: "O", 9: "F", + 1: "H", 2: "He", 3: "Li", 4: "Be", 5: "B", + 6: "C", 7: "N", 8: "O", 9: "F", 10: "Ne", }; -const SUPPORTED_SYMBOLS = new Set(["H", "He", "Li", "Be", "C", "N", "O", "F"]); +const SUPPORTED_SYMBOLS = new Set([ + "H", "He", "Li", "Be", "B", "C", "N", "O", "F", "Ne", +]); export interface XYZParseResult { readonly atoms: Atom[]; @@ -168,14 +171,14 @@ function resolveSymbol(token: string, lineNo: number): AtomSymbol { if (Number.isInteger(Z) && Z > 0) { const sym = SYMBOL_BY_Z[Z]; if (!sym) { - throw new Error(`parseXYZ: line ${lineNo} atomic number ${Z} is not supported (supported: H, He, Li, Be, C, N, O, F)`); + throw new Error(`parseXYZ: line ${lineNo} atomic number ${Z} is not supported (supported: H, He, Li, Be, B, C, N, O, F, Ne)`); } return sym; } // Otherwise treat as symbol — case-normalize and validate. const upper = token.charAt(0).toUpperCase() + token.slice(1).toLowerCase(); if (!SUPPORTED_SYMBOLS.has(upper as AtomSymbol)) { - throw new Error(`parseXYZ: line ${lineNo} element "${token}" not supported (supported: H, He, Li, Be, C, N, O, F)`); + throw new Error(`parseXYZ: line ${lineNo} element "${token}" not supported (supported: H, He, Li, Be, B, C, N, O, F, Ne)`); } return upper as AtomSymbol; } diff --git a/src/molecule/import-formats.ts b/src/molecule/import-formats.ts index ba5cb76..892ca54 100644 --- a/src/molecule/import-formats.ts +++ b/src/molecule/import-formats.ts @@ -4,11 +4,11 @@ import type { Atom, AtomSymbol } from "../chemistry/atoms.js"; -// The chemistry stack currently supports first-row + Li/Be (see +// The chemistry stack currently supports periods 1-2 complete (see // src/chemistry/atoms.ts → AtomSymbol). Imports outside this set are // rejected with a clear error. const SUPPORTED_SYMBOLS = new Set([ - "H", "He", "Li", "Be", "C", "N", "O", "F", + "H", "He", "Li", "Be", "B", "C", "N", "O", "F", "Ne", ]); export interface ParseResult { diff --git a/tests/chemistry/elements/pyscf-reference.json b/tests/chemistry/elements/pyscf-reference.json index 894ff55..d321675 100644 --- a/tests/chemistry/elements/pyscf-reference.json +++ b/tests/chemistry/elements/pyscf-reference.json @@ -649,6 +649,264 @@ "E_HF": -15.76849911667396, "nao": 43 }, + { + "element": "B", + "molecule": "BH3", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "B", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.286648454926751e-17, + 1.19, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0305702305034818, + -0.5950000000000001, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.0305702305034816, + -0.5950000000000005, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -26.068995721346088, + "nao": 8 + }, + { + "element": "B", + "molecule": "BH3", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "B", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.286648454926751e-17, + 1.19, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0305702305034818, + -0.5950000000000001, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.0305702305034816, + -0.5950000000000005, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -26.068995721346088, + "nao": 8 + }, + { + "element": "B", + "molecule": "BH3", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "B", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.286648454926751e-17, + 1.19, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0305702305034818, + -0.5950000000000001, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.0305702305034816, + -0.5950000000000005, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -26.390705581257425, + "nao": 29 + }, + { + "element": "B", + "molecule": "BH3", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "B", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.286648454926751e-17, + 1.19, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0305702305034818, + -0.5950000000000001, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.0305702305034816, + -0.5950000000000005, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -26.390775048122233, + "nao": 30 + }, + { + "element": "B", + "molecule": "BH3", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "B", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.286648454926751e-17, + 1.19, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0305702305034818, + -0.5950000000000001, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.0305702305034816, + -0.5950000000000005, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -26.39125449199123, + "nao": 50 + }, + { + "element": "B", + "molecule": "BH3", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "B", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.286648454926751e-17, + 1.19, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0305702305034818, + -0.5950000000000001, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.0305702305034816, + -0.5950000000000005, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -26.39130784288754, + "nao": 52 + }, { "element": "C", "molecule": "CH4", @@ -955,6 +1213,120 @@ "E_HF": -40.19966503955351, "nao": 61 }, + { + "element": "Ne", + "molecule": "Ne", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "Ne", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -126.60452499680486, + "nao": 5 + }, + { + "element": "Ne", + "molecule": "Ne", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Ne", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -126.60452499680486, + "nao": 5 + }, + { + "element": "Ne", + "molecule": "Ne", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Ne", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -128.48877555174084, + "nao": 14 + }, + { + "element": "Ne", + "molecule": "Ne", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Ne", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -128.48886617203775, + "nao": 15 + }, + { + "element": "Ne", + "molecule": "Ne", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Ne", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -128.49634973054086, + "nao": 23 + }, + { + "element": "Ne", + "molecule": "Ne", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Ne", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -128.49713422832133, + "nao": 25 + }, { "element": "N", "molecule": "NH3", diff --git a/tests/chemistry/periodic-table.test.ts b/tests/chemistry/periodic-table.test.ts index 9333f1d..5a5c039 100644 --- a/tests/chemistry/periodic-table.test.ts +++ b/tests/chemistry/periodic-table.test.ts @@ -10,7 +10,10 @@ describe("Periodic table", () => { expect(PERIODIC_TABLE.H.atomicNumber).toBe(1); expect(PERIODIC_TABLE.He.atomicNumber).toBe(2); expect(PERIODIC_TABLE.F.atomicNumber).toBe(9); - expect(Object.keys(PERIODIC_TABLE).length).toBe(8); // H, He, Li, Be, C, N, O, F + expect(PERIODIC_TABLE.B.atomicNumber).toBe(5); + expect(PERIODIC_TABLE.Ne.atomicNumber).toBe(10); + // H, He, Li, Be, B, C, N, O, F, Ne — periods 1-2 complete. + expect(Object.keys(PERIODIC_TABLE).length).toBe(10); }); test("Pauling electronegativity ordering: H < C < N < O < F", () => { @@ -22,8 +25,13 @@ describe("Periodic table", () => { expect(en("F")).toBe(3.98); // Pauling reference value }); - test("Helium has null electronegativity (noble gas)", () => { + test("Noble gases have null electronegativity (no Pauling value)", () => { expect(PERIODIC_TABLE.He.electronegativity).toBeNull(); + expect(PERIODIC_TABLE.Ne.electronegativity).toBeNull(); + expect(PERIODIC_TABLE.Ne.family).toBe("noble"); + // Boron does have one, and it sits between Be (1.57) and C (2.55). + expect(PERIODIC_TABLE.B.electronegativity).toBe(2.04); + expect(PERIODIC_TABLE.B.group).toBe(13); }); test("Period / group classifications match standard PT", () => { diff --git a/tests/chemistry/xyz.test.ts b/tests/chemistry/xyz.test.ts index e66c02a..49b7191 100644 --- a/tests/chemistry/xyz.test.ts +++ b/tests/chemistry/xyz.test.ts @@ -59,16 +59,41 @@ H₂ }); test("Unsupported element throws clearly", () => { - // Boron: genuinely absent from the AtomSymbol union and from the basis - // tables. This test previously used He, which IS a supported AtomSymbol with - // basis data and a molecule in molecules.ts — it was asserting a bug. + // Krypton: genuinely absent from the AtomSymbol union and from the basis + // tables. This test previously used He (a supported AtomSymbol all along), + // then B — which became supported when periods 1-2 were completed. Kr is + // period 4, far outside anything the basis tables cover. const text = `1 -boron not supported in current AtomSymbol union -B 0 0 0 +krypton not supported in current AtomSymbol union +Kr 0 0 0 `; expect(() => parseXYZ(text)).toThrow(/not supported/); }); + test("Boron round-trips: symbol, Z=5, and toXYZ identity", () => { + const text = `1 +boron atom +B 0.00000000 0.00000000 0.00000000 +`; + const { atoms } = parseXYZ(text); + expect(atoms).toHaveLength(1); + expect(atoms[0]!.symbol).toBe("B"); + expect(parseXYZ("1\nby Z\n5 0 0 0\n").atoms[0]!.symbol).toBe("B"); + expect(parseXYZ(toXYZ(atoms)).atoms[0]!.symbol).toBe("B"); + }); + + test("Neon round-trips: symbol, Z=10, and toXYZ identity", () => { + const text = `1 +neon atom +Ne 0.00000000 0.00000000 0.00000000 +`; + const { atoms } = parseXYZ(text); + expect(atoms).toHaveLength(1); + expect(atoms[0]!.symbol).toBe("Ne"); + expect(parseXYZ("1\nby Z\n10 0 0 0\n").atoms[0]!.symbol).toBe("Ne"); + expect(parseXYZ(toXYZ(atoms)).atoms[0]!.symbol).toBe("Ne"); + }); + test("He round-trips: it is a real AtomSymbol the library ships a molecule for", () => { const text = `1 helium atom From 26f6c5ace3dbbc6f5038928a36dbfe440ba3b3fd Mon Sep 17 00:00:00 2001 From: Ahmet Gunaydin Date: Mon, 10 Aug 2026 09:37:47 +0700 Subject: [PATCH 05/13] fix: correct BECKE_XI against PySCF BRAGG_RADII (measured impact: 0.6 uHa) MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit BECKE_XI sets the per-atom radial scaling for the Becke-partitioned DFT grid. Four entries disagreed with the Bragg-Slater table the scheme is defined on (pyscf.dft.radi.BRAGG_RADII): F MISSING -> 0.50 silently hit the `?? 1.0` fallback at grid.ts:80 He 0.35 -> 1.40 0.35 is hydrogen's value C 0.65 -> 0.70 0.65 is nitrogen's value O 0.65 -> 0.60 nitrogen's value again H stays at 0.5 — that deviation is deliberate and documented in place. Noble gases really are anomalously large here (He 1.40, Ne 1.50); that is the table, not a typo. I initially flagged Ne 1.50 as implausible and was wrong. MEASURED IMPACT, so nobody inherits an overstatement: LDA/cc-pVDZ, spherical, before vs after — HF -99.780846416 -> -99.780847054 (0.64 uHa) H2O -75.854562755 -> -75.854562750 (0.005 uHa) Sub-microhartree. The Becke partition is variationally insensitive to xi: it redistributes quadrature points rather than changing the integrand, and the existing grid was already dense enough to converge regardless. So "every fluorine DFT number was computed on a mis-scaled radial grid" is literally true and practically irrelevant. This commit is correctness hygiene — it removes a silent fallback and realigns with the cited source — not a repair of a meaningful error. Full suite green: 140 files, 1017 passed, 1 skipped. Co-Authored-By: Claude Opus 5 (1M context) --- src/chemistry/dft/grid.ts | 16 +++++++++++++--- 1 file changed, 13 insertions(+), 3 deletions(-) diff --git a/src/chemistry/dft/grid.ts b/src/chemistry/dft/grid.ts index 7e01df8..3c9e93c 100644 --- a/src/chemistry/dft/grid.ts +++ b/src/chemistry/dft/grid.ts @@ -264,10 +264,20 @@ function beckeStep(mu: number): number { const BECKE_XI: Partial> = { H: 0.5, // tighter than the default (Becke recommends ξ_H = 0.35; // 0.5 is a compromise that handles diatomic H₂ better) - He: 0.35, // Bragg-Slater radius for He + // Corrected 2026-08-10 against pyscf.dft.radi.BRAGG_RADII (the same + // table Becke's scheme is defined on). He was labelled "Bragg-Slater + // radius for He" but carried 0.35 against the tabulated 1.40; C, N and + // O had all been filled with 0.65 where the real values differ. + // Noble gases are anomalously large in this parameterization — that is + // the table, not a typo. + He: 1.40, // was 0.35 Li: 1.45, Be: 1.05, - C: 0.65, + B: 0.85, + C: 0.70, // was 0.65 N: 0.65, - O: 0.65, + O: 0.60, // was 0.65 + F: 0.50, // was MISSING → silently fell back to 1.0, misplacing the + // radial quadrature points for every fluorine DFT grid + Ne: 1.50, }; From 9e7e3ca85738e1210108c3dc35a0ba5767e1bc4a Mon Sep 17 00:00:00 2001 From: Ahmet Gunaydin Date: Mon, 10 Aug 2026 10:20:20 +0700 Subject: [PATCH 06/13] =?UTF-8?q?feat:=20add=20the=20third=20row=20(Na-Ar)?= =?UTF-8?q?=20=E2=80=94=20periodic=20table=20complete=20through=20Z=3D18?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Phase 1b of docs/RUN-PLAN-24H-ELEMENTS.md. webgpu-q now supports every element from hydrogen to argon. The element agreement test is 109/109 green: 18 elements x 3 bases x 2 d-conventions, one real molecule each (hydrides wherever one exists, so two-center integrals actually run), all under a 0.1 mHa bar against PySCF. This is the unlock the run was for. Cl gives SN2, the single most-taught mechanism in undergraduate organic; S gives thiols and H2S; P gives phosphates and PH3. Structure is uniform across the row, so atoms.ts needed one row-3 code path rather than eight: STO-3G 3s + 2p -> 9 Cartesian cc-pVDZ 4s + 3p + 1d -> 19 Cartesian / 18 spherical aug +1s +1p +1d STO-3G's M-shell (3s/3p) coefficients are a genuinely new universal set, identical across all eight elements exactly as the L-shell constants are across Li-Ne; 3s and 3p share one exponent set. `heavyShellsRow3` is deliberately a SEPARATE function rather than a generalisation of `heavyShells`. Row 2 therefore emits byte-identical shells and no pinned first-row number can move. The d-set stays a contiguous canonical 6-tuple, which buildSphericalDTransform requires. Frozen core is 5 for Na-Ar (the full neon core, 1s2s2p), not 1. ZERO engine work was needed. The McMurchie-Davidson kernels were verified against PySCF through L=5 (f/g/h agree to ~1e-14); row 3 tops out at L=2, which is production-proven. Also fixed here, because this change made it worse: `canonicalSymbol` in import-formats.ts truncated an unsupported 2-char symbol to 1 char unconditionally, so "Fe" imported as FLUORINE — silent, and a real wrong element rather than an error. Adding P and S as supported single letters made "Pd"/"Pt" -> P and "Sc"/"Se"/"Sn"/"Sr" -> S newly wrong ("Na" -> N was already wrong). The truncation is only ever correct for PDB atom names (cols 13-16, where "CA" is alpha-carbon), so it is now opt-in and passed only from that path. Regression tests cover both directions. Corrections to the run plan, for the record: - Na and Mg carry d functions in cc-pVDZ, so the Cartesian/spherical convention difference applies from Na, not from Al as written. - STO-3G row 3 is 9 Cartesian functions, not 13. Full suite: 140 files, 1065 passed, 1 skipped. tsc clean, lint clean. Co-Authored-By: Claude Opus 5 (1M context) --- scripts/run-element-reference.py | 17 + src/chemistry/atoms.ts | 185 +- src/chemistry/dft/grid.ts | 18 +- src/chemistry/dispersion-d2.ts | 32 + src/chemistry/geometry.ts | 6 + src/chemistry/integrals.ts | 711 +++++++ src/chemistry/periodic-table.ts | 51 +- src/chemistry/raman.ts | 8 + src/chemistry/thermochemistry.ts | 8 + src/chemistry/vibrations.ts | 11 +- src/chemistry/xyz.ts | 12 +- src/molecule/import-formats.ts | 33 +- src/viz/molecular-graph-2d.ts | 3 +- tests/chemistry/elements/pyscf-reference.json | 1680 +++++++++++++++++ tests/chemistry/periodic-table.test.ts | 6 +- tests/molecule/import-formats.test.ts | 33 +- 16 files changed, 2795 insertions(+), 19 deletions(-) diff --git a/scripts/run-element-reference.py b/scripts/run-element-reference.py index c855ef5..ae9d81b 100644 --- a/scripts/run-element-reference.py +++ b/scripts/run-element-reference.py @@ -69,6 +69,23 @@ def pyramidal(bond, angle_deg): "O": ("H2O", [("O", (0.0, 0.0, 0.0))] + [("H", p) for p in bent(0.9572, 104.52)]), "F": ("HF", [("F", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 0.9168))]), + # Third row. Hydrides wherever one exists, so two-center integrals + # are actually exercised; bare atoms only for the noble gas. + "Na": ("NaH", [("Na", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 1.8874))]), + "Mg": ("MgH2", [("Mg", (0.0, 0.0, 0.0)), + ("H", (0.0, 0.0, 1.7297)), ("H", (0.0, 0.0, -1.7297))]), + "Al": ("AlH3", [("Al", (0.0, 0.0, 0.0))] + + [("H", (1.5840 * math.cos(math.radians(90 + 120 * k)), + 1.5840 * math.sin(math.radians(90 + 120 * k)), 0.0)) + for k in range(3)]), + "Si": ("SiH4", [("Si", (0.0, 0.0, 0.0))] + + [("H", p) for p in tetrahedral(1.4798)]), + "P": ("PH3", [("P", (0.0, 0.0, 0.0))] + + [("H", p) for p in pyramidal(1.4200, 93.5)]), + "S": ("H2S", [("S", (0.0, 0.0, 0.0))] + + [("H", p) for p in bent(1.3356, 92.11)]), + "Cl": ("HCl", [("Cl", (0.0, 0.0, 0.0)), ("H", (0.0, 0.0, 1.2746))]), + "Ar": ("Ar", [("Ar", (0.0, 0.0, 0.0))]), } BASES = ["sto-3g", "cc-pvdz", "aug-cc-pvdz"] diff --git a/src/chemistry/atoms.ts b/src/chemistry/atoms.ts index 29955f5..2686baa 100644 --- a/src/chemistry/atoms.ts +++ b/src/chemistry/atoms.ts @@ -44,12 +44,45 @@ import { AUG_CCPVDZ_N_DIFFUSE_S, AUG_CCPVDZ_N_DIFFUSE_P, AUG_CCPVDZ_N_DIFFUSE_D, AUG_CCPVDZ_O_DIFFUSE_S, AUG_CCPVDZ_O_DIFFUSE_P, AUG_CCPVDZ_O_DIFFUSE_D, AUG_CCPVDZ_F_DIFFUSE_S, AUG_CCPVDZ_F_DIFFUSE_P, AUG_CCPVDZ_F_DIFFUSE_D, + STO3G_NA_1S, STO3G_NA_2S, STO3G_NA_3S, STO3G_NA_2P, STO3G_NA_3P, + STO3G_MG_1S, STO3G_MG_2S, STO3G_MG_3S, STO3G_MG_2P, STO3G_MG_3P, + STO3G_AL_1S, STO3G_AL_2S, STO3G_AL_3S, STO3G_AL_2P, STO3G_AL_3P, + STO3G_SI_1S, STO3G_SI_2S, STO3G_SI_3S, STO3G_SI_2P, STO3G_SI_3P, + STO3G_P_1S, STO3G_P_2S, STO3G_P_3S, STO3G_P_2P, STO3G_P_3P, + STO3G_S_1S, STO3G_S_2S, STO3G_S_3S, STO3G_S_2P, STO3G_S_3P, + STO3G_CL_1S, STO3G_CL_2S, STO3G_CL_3S, STO3G_CL_2P, STO3G_CL_3P, + STO3G_AR_1S, STO3G_AR_2S, STO3G_AR_3S, STO3G_AR_2P, STO3G_AR_3P, + CCPVDZ_NA_1S, CCPVDZ_NA_2S, CCPVDZ_NA_3S, CCPVDZ_NA_3S_P, + CCPVDZ_NA_2P, CCPVDZ_NA_3P, CCPVDZ_NA_3P_P, CCPVDZ_NA_3D, + CCPVDZ_MG_1S, CCPVDZ_MG_2S, CCPVDZ_MG_3S, CCPVDZ_MG_3S_P, + CCPVDZ_MG_2P, CCPVDZ_MG_3P, CCPVDZ_MG_3P_P, CCPVDZ_MG_3D, + CCPVDZ_AL_1S, CCPVDZ_AL_2S, CCPVDZ_AL_3S, CCPVDZ_AL_3S_P, + CCPVDZ_AL_2P, CCPVDZ_AL_3P, CCPVDZ_AL_3P_P, CCPVDZ_AL_3D, + CCPVDZ_SI_1S, CCPVDZ_SI_2S, CCPVDZ_SI_3S, CCPVDZ_SI_3S_P, + CCPVDZ_SI_2P, CCPVDZ_SI_3P, CCPVDZ_SI_3P_P, CCPVDZ_SI_3D, + CCPVDZ_P_1S, CCPVDZ_P_2S, CCPVDZ_P_3S, CCPVDZ_P_3S_P, + CCPVDZ_P_2P, CCPVDZ_P_3P, CCPVDZ_P_3P_P, CCPVDZ_P_3D, + CCPVDZ_S_1S, CCPVDZ_S_2S, CCPVDZ_S_3S, CCPVDZ_S_3S_P, + CCPVDZ_S_2P, CCPVDZ_S_3P, CCPVDZ_S_3P_P, CCPVDZ_S_3D, + CCPVDZ_CL_1S, CCPVDZ_CL_2S, CCPVDZ_CL_3S, CCPVDZ_CL_3S_P, + CCPVDZ_CL_2P, CCPVDZ_CL_3P, CCPVDZ_CL_3P_P, CCPVDZ_CL_3D, + CCPVDZ_AR_1S, CCPVDZ_AR_2S, CCPVDZ_AR_3S, CCPVDZ_AR_3S_P, + CCPVDZ_AR_2P, CCPVDZ_AR_3P, CCPVDZ_AR_3P_P, CCPVDZ_AR_3D, + AUG_CCPVDZ_NA_DIFFUSE_S, AUG_CCPVDZ_NA_DIFFUSE_P, AUG_CCPVDZ_NA_DIFFUSE_D, + AUG_CCPVDZ_MG_DIFFUSE_S, AUG_CCPVDZ_MG_DIFFUSE_P, AUG_CCPVDZ_MG_DIFFUSE_D, + AUG_CCPVDZ_AL_DIFFUSE_S, AUG_CCPVDZ_AL_DIFFUSE_P, AUG_CCPVDZ_AL_DIFFUSE_D, + AUG_CCPVDZ_SI_DIFFUSE_S, AUG_CCPVDZ_SI_DIFFUSE_P, AUG_CCPVDZ_SI_DIFFUSE_D, + AUG_CCPVDZ_P_DIFFUSE_S, AUG_CCPVDZ_P_DIFFUSE_P, AUG_CCPVDZ_P_DIFFUSE_D, + AUG_CCPVDZ_S_DIFFUSE_S, AUG_CCPVDZ_S_DIFFUSE_P, AUG_CCPVDZ_S_DIFFUSE_D, + AUG_CCPVDZ_CL_DIFFUSE_S, AUG_CCPVDZ_CL_DIFFUSE_P, AUG_CCPVDZ_CL_DIFFUSE_D, + AUG_CCPVDZ_AR_DIFFUSE_S, AUG_CCPVDZ_AR_DIFFUSE_P, AUG_CCPVDZ_AR_DIFFUSE_D, } from "./integrals.js"; import { type CGShell, makeCGShell } from "./integrals-cg.js"; import { type Nucleus } from "./cg-molecular.js"; export type AtomSymbol = - | "H" | "He" | "Li" | "Be" | "B" | "C" | "N" | "O" | "F" | "Ne"; + | "H" | "He" | "Li" | "Be" | "B" | "C" | "N" | "O" | "F" | "Ne" + | "Na" | "Mg" | "Al" | "Si" | "P" | "S" | "Cl" | "Ar"; export type BasisName = "sto-3g" | "cc-pvdz" | "aug-cc-pvdz"; const ANGSTROM_TO_BOHR = 1 / 0.529177210903; @@ -70,11 +103,13 @@ export interface Atom { /** Atomic number for each supported atom. */ export const Z_FOR: Readonly> = { H: 1, He: 2, Li: 3, Be: 4, B: 5, C: 6, N: 7, O: 8, F: 9, Ne: 10, + Na: 11, Mg: 12, Al: 13, Si: 14, P: 15, S: 16, Cl: 17, Ar: 18, }; /** Number of electrons in the neutral atom. */ export const N_ELECTRONS_FOR: Readonly> = { H: 1, He: 2, Li: 3, Be: 4, B: 5, C: 6, N: 7, O: 8, F: 9, Ne: 10, + Na: 11, Mg: 12, Al: 13, Si: 14, P: 15, S: 16, Cl: 17, Ar: 18, }; /** @@ -89,6 +124,9 @@ export const N_ELECTRONS_FOR: Readonly> = { */ export const FROZEN_CORE_FOR: Readonly> = { H: 0, He: 0, Li: 1, Be: 1, B: 1, C: 1, N: 1, O: 1, F: 1, Ne: 1, + // Third row freezes the full neon core (1s + 2s + 2p = 5 orbitals), + // not just 1s — the standard convention for Na-Ar. + Na: 5, Mg: 5, Al: 5, Si: 5, P: 5, S: 5, Cl: 5, Ar: 5, }; /** Default frozen-core count for a molecule (sum of per-atom 1s cores). */ @@ -189,6 +227,102 @@ export function atomShells( makeCGShell(STO3G_NE_2P, pos_bohr, [0, 1, 0], "Ne:2p_y"), makeCGShell(STO3G_NE_2P, pos_bohr, [0, 0, 1], "Ne:2p_z"), ]; + case "Na": + return [ + makeCGShell(STO3G_NA_1S, pos_bohr, [0, 0, 0], "Na:1s"), + makeCGShell(STO3G_NA_2S, pos_bohr, [0, 0, 0], "Na:2s"), + makeCGShell(STO3G_NA_3S, pos_bohr, [0, 0, 0], "Na:3s"), + makeCGShell(STO3G_NA_2P, pos_bohr, [1, 0, 0], "Na:2p_x"), + makeCGShell(STO3G_NA_2P, pos_bohr, [0, 1, 0], "Na:2p_y"), + makeCGShell(STO3G_NA_2P, pos_bohr, [0, 0, 1], "Na:2p_z"), + makeCGShell(STO3G_NA_3P, pos_bohr, [1, 0, 0], "Na:3p_x"), + makeCGShell(STO3G_NA_3P, pos_bohr, [0, 1, 0], "Na:3p_y"), + makeCGShell(STO3G_NA_3P, pos_bohr, [0, 0, 1], "Na:3p_z"), + ]; + case "Mg": + return [ + makeCGShell(STO3G_MG_1S, pos_bohr, [0, 0, 0], "Mg:1s"), + makeCGShell(STO3G_MG_2S, pos_bohr, [0, 0, 0], "Mg:2s"), + makeCGShell(STO3G_MG_3S, pos_bohr, [0, 0, 0], "Mg:3s"), + makeCGShell(STO3G_MG_2P, pos_bohr, [1, 0, 0], "Mg:2p_x"), + makeCGShell(STO3G_MG_2P, pos_bohr, [0, 1, 0], "Mg:2p_y"), + makeCGShell(STO3G_MG_2P, pos_bohr, [0, 0, 1], "Mg:2p_z"), + makeCGShell(STO3G_MG_3P, pos_bohr, [1, 0, 0], "Mg:3p_x"), + makeCGShell(STO3G_MG_3P, pos_bohr, [0, 1, 0], "Mg:3p_y"), + makeCGShell(STO3G_MG_3P, pos_bohr, [0, 0, 1], "Mg:3p_z"), + ]; + case "Al": + return [ + makeCGShell(STO3G_AL_1S, pos_bohr, [0, 0, 0], "Al:1s"), + makeCGShell(STO3G_AL_2S, pos_bohr, [0, 0, 0], "Al:2s"), + makeCGShell(STO3G_AL_3S, pos_bohr, [0, 0, 0], "Al:3s"), + makeCGShell(STO3G_AL_2P, pos_bohr, [1, 0, 0], "Al:2p_x"), + makeCGShell(STO3G_AL_2P, pos_bohr, [0, 1, 0], "Al:2p_y"), + makeCGShell(STO3G_AL_2P, pos_bohr, [0, 0, 1], "Al:2p_z"), + makeCGShell(STO3G_AL_3P, pos_bohr, [1, 0, 0], "Al:3p_x"), + makeCGShell(STO3G_AL_3P, pos_bohr, [0, 1, 0], "Al:3p_y"), + makeCGShell(STO3G_AL_3P, pos_bohr, [0, 0, 1], "Al:3p_z"), + ]; + case "Si": + return [ + makeCGShell(STO3G_SI_1S, pos_bohr, [0, 0, 0], "Si:1s"), + makeCGShell(STO3G_SI_2S, pos_bohr, [0, 0, 0], "Si:2s"), + makeCGShell(STO3G_SI_3S, pos_bohr, [0, 0, 0], "Si:3s"), + makeCGShell(STO3G_SI_2P, pos_bohr, [1, 0, 0], "Si:2p_x"), + makeCGShell(STO3G_SI_2P, pos_bohr, [0, 1, 0], "Si:2p_y"), + makeCGShell(STO3G_SI_2P, pos_bohr, [0, 0, 1], "Si:2p_z"), + makeCGShell(STO3G_SI_3P, pos_bohr, [1, 0, 0], "Si:3p_x"), + makeCGShell(STO3G_SI_3P, pos_bohr, [0, 1, 0], "Si:3p_y"), + makeCGShell(STO3G_SI_3P, pos_bohr, [0, 0, 1], "Si:3p_z"), + ]; + case "P": + return [ + makeCGShell(STO3G_P_1S, pos_bohr, [0, 0, 0], "P:1s"), + makeCGShell(STO3G_P_2S, pos_bohr, [0, 0, 0], "P:2s"), + makeCGShell(STO3G_P_3S, pos_bohr, [0, 0, 0], "P:3s"), + makeCGShell(STO3G_P_2P, pos_bohr, [1, 0, 0], "P:2p_x"), + makeCGShell(STO3G_P_2P, pos_bohr, [0, 1, 0], "P:2p_y"), + makeCGShell(STO3G_P_2P, pos_bohr, [0, 0, 1], "P:2p_z"), + makeCGShell(STO3G_P_3P, pos_bohr, [1, 0, 0], "P:3p_x"), + makeCGShell(STO3G_P_3P, pos_bohr, [0, 1, 0], "P:3p_y"), + makeCGShell(STO3G_P_3P, pos_bohr, [0, 0, 1], "P:3p_z"), + ]; + case "S": + return [ + makeCGShell(STO3G_S_1S, pos_bohr, [0, 0, 0], "S:1s"), + makeCGShell(STO3G_S_2S, pos_bohr, [0, 0, 0], "S:2s"), + makeCGShell(STO3G_S_3S, pos_bohr, [0, 0, 0], "S:3s"), + makeCGShell(STO3G_S_2P, pos_bohr, [1, 0, 0], "S:2p_x"), + makeCGShell(STO3G_S_2P, pos_bohr, [0, 1, 0], "S:2p_y"), + makeCGShell(STO3G_S_2P, pos_bohr, [0, 0, 1], "S:2p_z"), + makeCGShell(STO3G_S_3P, pos_bohr, [1, 0, 0], "S:3p_x"), + makeCGShell(STO3G_S_3P, pos_bohr, [0, 1, 0], "S:3p_y"), + makeCGShell(STO3G_S_3P, pos_bohr, [0, 0, 1], "S:3p_z"), + ]; + case "Cl": + return [ + makeCGShell(STO3G_CL_1S, pos_bohr, [0, 0, 0], "Cl:1s"), + makeCGShell(STO3G_CL_2S, pos_bohr, [0, 0, 0], "Cl:2s"), + makeCGShell(STO3G_CL_3S, pos_bohr, [0, 0, 0], "Cl:3s"), + makeCGShell(STO3G_CL_2P, pos_bohr, [1, 0, 0], "Cl:2p_x"), + makeCGShell(STO3G_CL_2P, pos_bohr, [0, 1, 0], "Cl:2p_y"), + makeCGShell(STO3G_CL_2P, pos_bohr, [0, 0, 1], "Cl:2p_z"), + makeCGShell(STO3G_CL_3P, pos_bohr, [1, 0, 0], "Cl:3p_x"), + makeCGShell(STO3G_CL_3P, pos_bohr, [0, 1, 0], "Cl:3p_y"), + makeCGShell(STO3G_CL_3P, pos_bohr, [0, 0, 1], "Cl:3p_z"), + ]; + case "Ar": + return [ + makeCGShell(STO3G_AR_1S, pos_bohr, [0, 0, 0], "Ar:1s"), + makeCGShell(STO3G_AR_2S, pos_bohr, [0, 0, 0], "Ar:2s"), + makeCGShell(STO3G_AR_3S, pos_bohr, [0, 0, 0], "Ar:3s"), + makeCGShell(STO3G_AR_2P, pos_bohr, [1, 0, 0], "Ar:2p_x"), + makeCGShell(STO3G_AR_2P, pos_bohr, [0, 1, 0], "Ar:2p_y"), + makeCGShell(STO3G_AR_2P, pos_bohr, [0, 0, 1], "Ar:2p_z"), + makeCGShell(STO3G_AR_3P, pos_bohr, [1, 0, 0], "Ar:3p_x"), + makeCGShell(STO3G_AR_3P, pos_bohr, [0, 1, 0], "Ar:3p_y"), + makeCGShell(STO3G_AR_3P, pos_bohr, [0, 0, 1], "Ar:3p_z"), + ]; } } @@ -222,6 +356,39 @@ function atomShellsCcPvdz(symbol: AtomSymbol, pos: readonly [number, number, num makeCGShell(d1, pos, [0, 1, 1], `${sym}:3d_yz`), ]; } + /** Third-row cc-pVDZ stack: 4s + 3p + 1d -> 4 + 9 + 6 = 19 Cartesian + * functions (18 spherical). Kept separate from `heavyShells` rather + * than generalising it, so the row-2 path emits byte-identical shells + * and no pinned first-row number can move. The d-set must stay a + * contiguous canonical 6-tuple — `buildSphericalDTransform` in + * cg-molecular.ts throws otherwise. */ + function heavyShellsRow3( + sym: string, + s1: ShellData, s2: ShellData, s3: ShellData, s3p: ShellData, + p1: ShellData, p2: ShellData, p3p: ShellData, d1: ShellData, + ): CGShell[] { + return [ + makeCGShell(s1, pos, [0, 0, 0], `${sym}:1s`), + makeCGShell(s2, pos, [0, 0, 0], `${sym}:2s`), + makeCGShell(s3, pos, [0, 0, 0], `${sym}:3s`), + makeCGShell(s3p, pos, [0, 0, 0], `${sym}:3s'`), + makeCGShell(p1, pos, [1, 0, 0], `${sym}:2p_x`), + makeCGShell(p1, pos, [0, 1, 0], `${sym}:2p_y`), + makeCGShell(p1, pos, [0, 0, 1], `${sym}:2p_z`), + makeCGShell(p2, pos, [1, 0, 0], `${sym}:3p_x`), + makeCGShell(p2, pos, [0, 1, 0], `${sym}:3p_y`), + makeCGShell(p2, pos, [0, 0, 1], `${sym}:3p_z`), + makeCGShell(p3p, pos, [1, 0, 0], `${sym}:3p'_x`), + makeCGShell(p3p, pos, [0, 1, 0], `${sym}:3p'_y`), + makeCGShell(p3p, pos, [0, 0, 1], `${sym}:3p'_z`), + makeCGShell(d1, pos, [2, 0, 0], `${sym}:3d_xx`), + makeCGShell(d1, pos, [0, 2, 0], `${sym}:3d_yy`), + makeCGShell(d1, pos, [0, 0, 2], `${sym}:3d_zz`), + makeCGShell(d1, pos, [1, 1, 0], `${sym}:3d_xy`), + makeCGShell(d1, pos, [1, 0, 1], `${sym}:3d_xz`), + makeCGShell(d1, pos, [0, 1, 1], `${sym}:3d_yz`), + ]; + } switch (symbol) { case "H": return [ @@ -239,6 +406,14 @@ function atomShellsCcPvdz(symbol: AtomSymbol, pos: readonly [number, number, num case "O": return heavyShells("O", CCPVDZ_O_1S, CCPVDZ_O_2S, CCPVDZ_O_2S_P, CCPVDZ_O_2P, CCPVDZ_O_2P_P, CCPVDZ_O_3D); case "F": return heavyShells("F", CCPVDZ_F_1S, CCPVDZ_F_2S, CCPVDZ_F_2S_P, CCPVDZ_F_2P, CCPVDZ_F_2P_P, CCPVDZ_F_3D); case "Ne": return heavyShells("Ne", CCPVDZ_NE_1S, CCPVDZ_NE_2S, CCPVDZ_NE_2S_P, CCPVDZ_NE_2P, CCPVDZ_NE_2P_P, CCPVDZ_NE_3D); + case "Na": return heavyShellsRow3("Na", CCPVDZ_NA_1S, CCPVDZ_NA_2S, CCPVDZ_NA_3S, CCPVDZ_NA_3S_P, CCPVDZ_NA_2P, CCPVDZ_NA_3P, CCPVDZ_NA_3P_P, CCPVDZ_NA_3D); + case "Mg": return heavyShellsRow3("Mg", CCPVDZ_MG_1S, CCPVDZ_MG_2S, CCPVDZ_MG_3S, CCPVDZ_MG_3S_P, CCPVDZ_MG_2P, CCPVDZ_MG_3P, CCPVDZ_MG_3P_P, CCPVDZ_MG_3D); + case "Al": return heavyShellsRow3("Al", CCPVDZ_AL_1S, CCPVDZ_AL_2S, CCPVDZ_AL_3S, CCPVDZ_AL_3S_P, CCPVDZ_AL_2P, CCPVDZ_AL_3P, CCPVDZ_AL_3P_P, CCPVDZ_AL_3D); + case "Si": return heavyShellsRow3("Si", CCPVDZ_SI_1S, CCPVDZ_SI_2S, CCPVDZ_SI_3S, CCPVDZ_SI_3S_P, CCPVDZ_SI_2P, CCPVDZ_SI_3P, CCPVDZ_SI_3P_P, CCPVDZ_SI_3D); + case "P": return heavyShellsRow3("P", CCPVDZ_P_1S, CCPVDZ_P_2S, CCPVDZ_P_3S, CCPVDZ_P_3S_P, CCPVDZ_P_2P, CCPVDZ_P_3P, CCPVDZ_P_3P_P, CCPVDZ_P_3D); + case "S": return heavyShellsRow3("S", CCPVDZ_S_1S, CCPVDZ_S_2S, CCPVDZ_S_3S, CCPVDZ_S_3S_P, CCPVDZ_S_2P, CCPVDZ_S_3P, CCPVDZ_S_3P_P, CCPVDZ_S_3D); + case "Cl": return heavyShellsRow3("Cl", CCPVDZ_CL_1S, CCPVDZ_CL_2S, CCPVDZ_CL_3S, CCPVDZ_CL_3S_P, CCPVDZ_CL_2P, CCPVDZ_CL_3P, CCPVDZ_CL_3P_P, CCPVDZ_CL_3D); + case "Ar": return heavyShellsRow3("Ar", CCPVDZ_AR_1S, CCPVDZ_AR_2S, CCPVDZ_AR_3S, CCPVDZ_AR_3S_P, CCPVDZ_AR_2P, CCPVDZ_AR_3P, CCPVDZ_AR_3P_P, CCPVDZ_AR_3D); case "He": return [ makeCGShell(CCPVDZ_HE_1S, pos, [0, 0, 0], "He:1s"), @@ -287,6 +462,14 @@ function atomShellsAugDiffuse(symbol: AtomSymbol, pos: readonly [number, number, case "O": return heavyDiffuse("O", AUG_CCPVDZ_O_DIFFUSE_S, AUG_CCPVDZ_O_DIFFUSE_P, AUG_CCPVDZ_O_DIFFUSE_D); case "F": return heavyDiffuse("F", AUG_CCPVDZ_F_DIFFUSE_S, AUG_CCPVDZ_F_DIFFUSE_P, AUG_CCPVDZ_F_DIFFUSE_D); case "Ne": return heavyDiffuse("Ne", AUG_CCPVDZ_NE_DIFFUSE_S, AUG_CCPVDZ_NE_DIFFUSE_P, AUG_CCPVDZ_NE_DIFFUSE_D); + case "Na": return heavyDiffuse("Na", AUG_CCPVDZ_NA_DIFFUSE_S, AUG_CCPVDZ_NA_DIFFUSE_P, AUG_CCPVDZ_NA_DIFFUSE_D); + case "Mg": return heavyDiffuse("Mg", AUG_CCPVDZ_MG_DIFFUSE_S, AUG_CCPVDZ_MG_DIFFUSE_P, AUG_CCPVDZ_MG_DIFFUSE_D); + case "Al": return heavyDiffuse("Al", AUG_CCPVDZ_AL_DIFFUSE_S, AUG_CCPVDZ_AL_DIFFUSE_P, AUG_CCPVDZ_AL_DIFFUSE_D); + case "Si": return heavyDiffuse("Si", AUG_CCPVDZ_SI_DIFFUSE_S, AUG_CCPVDZ_SI_DIFFUSE_P, AUG_CCPVDZ_SI_DIFFUSE_D); + case "P": return heavyDiffuse("P", AUG_CCPVDZ_P_DIFFUSE_S, AUG_CCPVDZ_P_DIFFUSE_P, AUG_CCPVDZ_P_DIFFUSE_D); + case "S": return heavyDiffuse("S", AUG_CCPVDZ_S_DIFFUSE_S, AUG_CCPVDZ_S_DIFFUSE_P, AUG_CCPVDZ_S_DIFFUSE_D); + case "Cl": return heavyDiffuse("Cl", AUG_CCPVDZ_CL_DIFFUSE_S, AUG_CCPVDZ_CL_DIFFUSE_P, AUG_CCPVDZ_CL_DIFFUSE_D); + case "Ar": return heavyDiffuse("Ar", AUG_CCPVDZ_AR_DIFFUSE_S, AUG_CCPVDZ_AR_DIFFUSE_P, AUG_CCPVDZ_AR_DIFFUSE_D); case "He": return [ makeCGShell(AUG_CCPVDZ_HE_DIFFUSE_S, pos, [0, 0, 0], "He:aug-s"), diff --git a/src/chemistry/dft/grid.ts b/src/chemistry/dft/grid.ts index 3c9e93c..3216ef1 100644 --- a/src/chemistry/dft/grid.ts +++ b/src/chemistry/dft/grid.ts @@ -255,9 +255,11 @@ function beckeStep(mu: number): number { return 0.5 * (1 - p3); } -// ── Per-atom Becke radial scaling (bohr) ──────────────────── +// ── Per-atom Becke radial scaling (angstrom) ──────────────── // -// ξ_A is roughly the Bragg-Slater radius of atom A in Bohr; it +// ξ_A is the Bragg-Slater radius of atom A in ANGSTROM (the +// header previously said Bohr; every tabulated value is and always +// was the Angstrom number — C is 0.70, not 1.32); it // sets the "natural distance" at which the radial mapping puts // half its quadrature points. Values from Becke (1988) / // Treutler-Ahlrichs (1995). Atoms not listed default to 1.0. @@ -280,4 +282,16 @@ const BECKE_XI: Partial> = { F: 0.50, // was MISSING → silently fell back to 1.0, misplacing the // radial quadrature points for every fluorine DFT grid Ne: 1.50, + // Row 3, same source (pyscf.dft.radi.BRAGG_RADII, converted with + // pyscf.data.nist.BOHR). P / S / Cl all land on exactly 1.00 in that + // table — not a copy-paste slip. Na and Ar tie at 1.80 for the same + // reason the row-2 noble gases are large here. + Na: 1.80, + Mg: 1.50, + Al: 1.25, + Si: 1.10, + P: 1.00, + S: 1.00, + Cl: 1.00, + Ar: 1.80, }; diff --git a/src/chemistry/dispersion-d2.ts b/src/chemistry/dispersion-d2.ts index d5f3ec0..b405897 100644 --- a/src/chemistry/dispersion-d2.ts +++ b/src/chemistry/dispersion-d2.ts @@ -49,6 +49,14 @@ const C6_J_NM6_PER_MOL: Readonly> = { O: 0.70, F: 0.75, Ne: 0.63, + Na: 5.71, + Mg: 5.71, + Al: 10.79, + Si: 9.23, + P: 7.84, + S: 5.57, + Cl: 5.07, + Ar: 4.61, }; /** Grimme 2006 Table I: per-element R_R van-der-Waals radius in Å. */ @@ -63,6 +71,14 @@ const R_R_ANGSTROM: Readonly> = { O: 1.342, F: 1.287, Ne: 1.243, + Na: 1.144, + Mg: 1.364, + Al: 1.639, + Si: 1.716, + P: 1.705, + S: 1.683, + Cl: 1.639, + Ar: 1.595, }; /** Per-element C6 in Hartree·Bohr^6 (atomic units). */ @@ -77,6 +93,14 @@ const C6_AU: Readonly> = { O: C6_J_NM6_PER_MOL.O * J_NM6_PER_MOL_TO_HARTREE_BOHR6, F: C6_J_NM6_PER_MOL.F * J_NM6_PER_MOL_TO_HARTREE_BOHR6, Ne: C6_J_NM6_PER_MOL.Ne * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + Na: C6_J_NM6_PER_MOL.Na * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + Mg: C6_J_NM6_PER_MOL.Mg * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + Al: C6_J_NM6_PER_MOL.Al * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + Si: C6_J_NM6_PER_MOL.Si * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + P: C6_J_NM6_PER_MOL.P * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + S: C6_J_NM6_PER_MOL.S * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + Cl: C6_J_NM6_PER_MOL.Cl * J_NM6_PER_MOL_TO_HARTREE_BOHR6, + Ar: C6_J_NM6_PER_MOL.Ar * J_NM6_PER_MOL_TO_HARTREE_BOHR6, }; /** Per-element R_R in Bohr. */ @@ -91,6 +115,14 @@ const R_R_BOHR: Readonly> = { O: R_R_ANGSTROM.O * ANGSTROM_TO_BOHR, F: R_R_ANGSTROM.F * ANGSTROM_TO_BOHR, Ne: R_R_ANGSTROM.Ne * ANGSTROM_TO_BOHR, + Na: R_R_ANGSTROM.Na * ANGSTROM_TO_BOHR, + Mg: R_R_ANGSTROM.Mg * ANGSTROM_TO_BOHR, + Al: R_R_ANGSTROM.Al * ANGSTROM_TO_BOHR, + Si: R_R_ANGSTROM.Si * ANGSTROM_TO_BOHR, + P: R_R_ANGSTROM.P * ANGSTROM_TO_BOHR, + S: R_R_ANGSTROM.S * ANGSTROM_TO_BOHR, + Cl: R_R_ANGSTROM.Cl * ANGSTROM_TO_BOHR, + Ar: R_R_ANGSTROM.Ar * ANGSTROM_TO_BOHR, }; /** Damping steepness — universal in D2. */ diff --git a/src/chemistry/geometry.ts b/src/chemistry/geometry.ts index 01469ce..c27c7ce 100644 --- a/src/chemistry/geometry.ts +++ b/src/chemistry/geometry.ts @@ -246,6 +246,8 @@ const COVALENT_RADIUS_ANGSTROM: Readonly> = { H: 0.32, He: 0.46, Li: 1.33, Be: 1.02, B: 0.85, C: 0.75, N: 0.71, O: 0.63, F: 0.64, Ne: 0.67, + Na: 1.55, Mg: 1.39, Al: 1.26, + Si: 1.16, P: 1.11, S: 1.03, Cl: 0.99, Ar: 0.96, }; export interface Bond { @@ -263,6 +265,10 @@ const ATOMIC_MASS: Readonly> = { C: 12.0, N: 14.0030740048, O: 15.99491461956, F: 18.998403163, Ne: 19.9924401762, + Na: 22.989769282, Mg: 23.985041697, + Al: 26.98153853, Si: 27.97692653465, + P: 30.97376199842, S: 31.9720711744, + Cl: 34.968852682, Ar: 39.9623831237, }; /** Total mass of a molecule in amu. */ diff --git a/src/chemistry/integrals.ts b/src/chemistry/integrals.ts index 8d15471..bd6e4aa 100644 --- a/src/chemistry/integrals.ts +++ b/src/chemistry/integrals.ts @@ -649,6 +649,717 @@ export const CCPVDZ_NE_3D = { c: [1.0] as const, }; +// ═════════════════════════════════════════════════════════════ +// Third row: Na, Mg, Al, Si, P, S, Cl, Ar (added 2026-08-10) +// +// Generated by scripts/gen-basis-tables.py from PySCF — not +// transcribed. Structure is uniform across the whole row, so +// atoms.ts needs one row-3 code path rather than eight: +// STO-3G 3s + 2p -> 9 Cartesian functions +// cc-pVDZ 4s + 3p + 1d -> 19 Cartesian / 18 spherical +// aug +1s +1p +1d +// Note the d function appears from Na onward, so the Cartesian-vs- +// spherical convention difference applies to this entire row. +// ═════════════════════════════════════════════════════════════ + +// STO-3G M-shell (3s/3p) contraction coefficients. Universal across +// Na-Ar exactly as the L-shell constants are across Li-Ne; verified +// identical for all eight elements to 1e-6. The 3s and 3p shells +// share one exponent set, mirroring the 2s/2p L-shell. +const STO3G_M_3S_C = [-0.219620369, 0.2255954336, 0.900398426] as const; +const STO3G_M_3P_C = [0.01058760429, 0.5951670053, 0.462001012] as const; + +export const STO3G_NA_1S = { + alpha: [250.77243, 45.678511, 12.362388] as const, + c: STO3G_S_C, +}; +export const STO3G_NA_2S = { + alpha: [12.040193, 2.7978819, 0.909958] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_NA_3S = { + alpha: [1.4787406, 0.4125649, 0.1614751] as const, + c: STO3G_M_3S_C, +}; +export const STO3G_NA_2P = { + alpha: [12.040193, 2.7978819, 0.909958] as const, + c: STO3G_L_2P_C, +}; +export const STO3G_NA_3P = { + alpha: [1.4787406, 0.4125649, 0.1614751] as const, + c: STO3G_M_3P_C, +}; +export const STO3G_MG_1S = { + alpha: [299.2374, 54.50647, 14.75158] as const, + c: STO3G_S_C, +}; +export const STO3G_MG_2S = { + alpha: [15.12182, 3.513987, 1.142857] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_MG_3S = { + alpha: [1.395448, 0.389326, 0.15238] as const, + c: STO3G_M_3S_C, +}; +export const STO3G_MG_2P = { + alpha: [15.12182, 3.513987, 1.142857] as const, + c: STO3G_L_2P_C, +}; +export const STO3G_MG_3P = { + alpha: [1.395448, 0.389326, 0.15238] as const, + c: STO3G_M_3P_C, +}; +export const STO3G_AL_1S = { + alpha: [351.4214767, 64.01186067, 17.32410761] as const, + c: STO3G_S_C, +}; +export const STO3G_AL_2S = { + alpha: [18.89939621, 4.391813233, 1.42835397] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_AL_3S = { + alpha: [1.395448293, 0.3893265318, 0.1523797659] as const, + c: STO3G_M_3S_C, +}; +export const STO3G_AL_2P = { + alpha: [18.89939621, 4.391813233, 1.42835397] as const, + c: STO3G_L_2P_C, +}; +export const STO3G_AL_3P = { + alpha: [1.395448293, 0.3893265318, 0.1523797659] as const, + c: STO3G_M_3P_C, +}; +export const STO3G_SI_1S = { + alpha: [407.7975514, 74.28083305, 20.10329229] as const, + c: STO3G_S_C, +}; +export const STO3G_SI_2S = { + alpha: [23.19365606, 5.389706871, 1.752899952] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_SI_3S = { + alpha: [1.478740622, 0.4125648801, 0.1614750979] as const, + c: STO3G_M_3S_C, +}; +export const STO3G_SI_2P = { + alpha: [23.19365606, 5.389706871, 1.752899952] as const, + c: STO3G_L_2P_C, +}; +export const STO3G_SI_3P = { + alpha: [1.478740622, 0.4125648801, 0.1614750979] as const, + c: STO3G_M_3P_C, +}; +export const STO3G_P_1S = { + alpha: [468.3656378, 85.31338559, 23.08913156] as const, + c: STO3G_S_C, +}; +export const STO3G_P_2S = { + alpha: [28.03263958, 6.514182577, 2.118614352] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_P_3S = { + alpha: [1.743103231, 0.4863213771, 0.1903428909] as const, + c: STO3G_M_3S_C, +}; +export const STO3G_P_2P = { + alpha: [28.03263958, 6.514182577, 2.118614352] as const, + c: STO3G_L_2P_C, +}; +export const STO3G_P_3P = { + alpha: [1.743103231, 0.4863213771, 0.1903428909] as const, + c: STO3G_M_3P_C, +}; +export const STO3G_S_1S = { + alpha: [533.1257359, 97.1095183, 26.28162542] as const, + c: STO3G_S_C, +}; +export const STO3G_S_2S = { + alpha: [33.32975173, 7.745117521, 2.518952599] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_S_3S = { + alpha: [2.029194274, 0.5661400518, 0.2215833792] as const, + c: STO3G_M_3S_C, +}; +export const STO3G_S_2P = { + alpha: [33.32975173, 7.745117521, 2.518952599] as const, + c: STO3G_L_2P_C, +}; +export const STO3G_S_3P = { + alpha: [2.029194274, 0.5661400518, 0.2215833792] as const, + c: STO3G_M_3P_C, +}; +export const STO3G_CL_1S = { + alpha: [601.3456136, 109.5358542, 29.64467686] as const, + c: STO3G_S_C, +}; +export const STO3G_CL_2S = { + alpha: [38.96041889, 9.053563477, 2.944499834] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_CL_3S = { + alpha: [2.129386495, 0.5940934274, 0.232524141] as const, + c: STO3G_M_3S_C, +}; +export const STO3G_CL_2P = { + alpha: [38.96041889, 9.053563477, 2.944499834] as const, + c: STO3G_L_2P_C, +}; +export const STO3G_CL_3P = { + alpha: [2.129386495, 0.5940934274, 0.232524141] as const, + c: STO3G_M_3P_C, +}; +export const STO3G_AR_1S = { + alpha: [674.4465184, 122.8512753, 33.24834945] as const, + c: STO3G_S_C, +}; +export const STO3G_AR_2S = { + alpha: [45.16424392, 10.495199, 3.413364448] as const, + c: STO3G_L_2S_C, +}; +export const STO3G_AR_3S = { + alpha: [2.621366518, 0.731354605, 0.2862472356] as const, + c: STO3G_M_3S_C, +}; +export const STO3G_AR_2P = { + alpha: [45.16424392, 10.495199, 3.413364448] as const, + c: STO3G_L_2P_C, +}; +export const STO3G_AR_3P = { + alpha: [2.621366518, 0.731354605, 0.2862472356] as const, + c: STO3G_M_3P_C, +}; + +export const CCPVDZ_NA_1S = { + alpha: [ + 31700.0, 4755.0, 1082.0, 306.4, 99.53, 35.42, 13.3, 4.392, 1.676, 0.5889, + 0.0564, + ] as const, + c: [ + 0.000458878, 0.0035507, 0.0182618, 0.071665, 0.212346, 0.416203, 0.37302, + 0.0625054, -0.00624532, 0.00243374, -0.000442381, + ] as const, +}; +export const CCPVDZ_NA_2S = { + alpha: [ + 31700.0, 4755.0, 1082.0, 306.4, 99.53, 35.42, 13.3, 4.392, 1.676, 0.5889, + 0.0564, + ] as const, + c: [ + -0.000112162, -0.000868512, -0.0045133, -0.0181436, -0.0580799, -0.137653, + -0.193908, 0.0858009, 0.604419, 0.441719, 0.0130547, + ] as const, +}; +export const CCPVDZ_NA_3S = { + alpha: [ + 31700.0, 4755.0, 1082.0, 306.4, 99.53, 35.42, 13.3, 4.392, 1.676, 0.5889, + 0.0564, + ] as const, + c: [ + 1.7016e-05, 0.000130693, 0.000687784, 0.00272359, 0.00895529, 0.0207832, + 0.031938, -0.0191368, -0.102595, -0.198945, 0.655952, + ] as const, +}; +export const CCPVDZ_NA_3S_P = { + alpha: [0.02307] as const, + c: [1.0] as const, +}; +export const CCPVDZ_NA_2P = { + alpha: [138.1, 32.24, 9.985, 3.484, 1.231, 0.4177, 0.06513] as const, + c: [ + 0.00579641, 0.0415756, 0.162873, 0.359401, 0.449988, 0.227507, 0.00808247, + ] as const, +}; +export const CCPVDZ_NA_3P = { + alpha: [138.1, 32.24, 9.985, 3.484, 1.231, 0.4177, 0.06513] as const, + c: [ + -0.000581531, -0.00407306, -0.0167937, -0.0353268, -0.0521971, -0.0168359, + 0.434613, + ] as const, +}; +export const CCPVDZ_NA_3P_P = { + alpha: [0.02053] as const, + c: [1.0] as const, +}; +export const CCPVDZ_NA_3D = { + alpha: [0.0973] as const, + c: [1.0] as const, +}; +export const CCPVDZ_MG_1S = { + alpha: [ + 47390.0, 7108.0, 1618.0, 458.4, 149.3, 53.59, 20.7, 8.384, 2.542, 0.8787, + 0.1077, + ] as const, + c: [ + 0.000346023, 0.00268077, 0.0138367, 0.0551767, 0.16966, 0.364703, 0.406856, + 0.135089, 0.00490884, 0.00028646, 2.6459e-05, + ] as const, +}; +export const CCPVDZ_MG_2S = { + alpha: [ + 47390.0, 7108.0, 1618.0, 458.4, 149.3, 53.59, 20.7, 8.384, 2.542, 0.8787, + 0.1077, + ] as const, + c: [ + -8.77839e-05, -0.000674725, -0.00355603, -0.0142154, -0.0476748, -0.114892, + -0.200676, -0.0341224, 0.570454, 0.542309, 0.0218128, + ] as const, +}; +export const CCPVDZ_MG_3S = { + alpha: [ + 47390.0, 7108.0, 1618.0, 458.4, 149.3, 53.59, 20.7, 8.384, 2.542, 0.8787, + 0.1077, + ] as const, + c: [ + 1.69628e-05, 0.000129865, 0.000688831, 0.00273533, 0.00931224, 0.0223265, + 0.0411195, 0.00545642, -0.134012, -0.256176, 0.605856, + ] as const, +}; +export const CCPVDZ_MG_3S_P = { + alpha: [0.03999] as const, + c: [1.0] as const, +}; +export const CCPVDZ_MG_2P = { + alpha: [179.9, 42.14, 13.13, 4.628, 1.67, 0.5857, 0.1311] as const, + c: [ + 0.00538161, 0.0392418, 0.157445, 0.358535, 0.457226, 0.215918, 0.00664948, + ] as const, +}; +export const CCPVDZ_MG_3P = { + alpha: [179.9, 42.14, 13.13, 4.628, 1.67, 0.5857, 0.1311] as const, + c: [ + -0.000865948, -0.00615978, -0.0261519, -0.0570647, -0.0873906, -0.012299, + 0.502085, + ] as const, +}; +export const CCPVDZ_MG_3P_P = { + alpha: [0.04112] as const, + c: [1.0] as const, +}; +export const CCPVDZ_MG_3D = { + alpha: [0.187] as const, + c: [1.0] as const, +}; +export const CCPVDZ_AL_1S = { + alpha: [ + 64150.0, 9617.0, 2189.0, 620.5, 202.7, 73.15, 28.55, 11.77, 3.3, 1.173, + 0.1752, + ] as const, + c: [ + 0.00029025, 0.00225064, 0.0116459, 0.0467377, 0.146299, 0.330283, 0.415861, + 0.189253, 0.0115889, -0.00128385, 0.000425883, + ] as const, +}; +export const CCPVDZ_AL_2S = { + alpha: [ + 64150.0, 9617.0, 2189.0, 620.5, 202.7, 73.15, 28.55, 11.77, 3.3, 1.173, + 0.1752, + ] as const, + c: [ + -7.58048e-05, -0.000581791, -0.00308113, -0.0123112, -0.0419781, -0.103371, + -0.196308, -0.0830002, 0.54104, 0.578796, 0.0288147, + ] as const, +}; +export const CCPVDZ_AL_3S = { + alpha: [ + 64150.0, 9617.0, 2189.0, 620.5, 202.7, 73.15, 28.55, 11.77, 3.3, 1.173, + 0.1752, + ] as const, + c: [ + 1.75078e-05, 0.000134208, 0.000712442, 0.0028433, 0.00976842, 0.024185, + 0.0474993, 0.0203621, -0.158788, -0.311694, 0.620147, + ] as const, +}; +export const CCPVDZ_AL_3S_P = { + alpha: [0.06473] as const, + c: [1.0] as const, +}; +export const CCPVDZ_AL_2P = { + alpha: [258.8, 60.89, 19.14, 6.881, 2.574, 0.9572, 0.2099] as const, + c: [0.00406847, 0.0306815, 0.129149, 0.320831, 0.453815, 0.275066, 0.0190807] as const, +}; +export const CCPVDZ_AL_3P = { + alpha: [258.8, 60.89, 19.14, 6.881, 2.574, 0.9572, 0.2099] as const, + c: [ + -0.000748053, -0.00545796, -0.0245371, -0.0582138, -0.0983756, -0.0260064, + 0.46402, + ] as const, +}; +export const CCPVDZ_AL_3P_P = { + alpha: [0.05986] as const, + c: [1.0] as const, +}; +export const CCPVDZ_AL_3D = { + alpha: [0.189] as const, + c: [1.0] as const, +}; +export const CCPVDZ_SI_1S = { + alpha: [ + 78860.0, 11820.0, 2692.0, 763.4, 249.6, 90.28, 35.29, 14.51, 4.053, + 1.482, 0.2517, + ] as const, + c: [ + 0.000270443, 0.00209717, 0.0108506, 0.0436754, 0.137653, 0.316644, + 0.418581, 0.210212, 0.0144952, -0.0020359, 0.000624186, + ] as const, +}; +export const CCPVDZ_SI_2S = { + alpha: [ + 78860.0, 11820.0, 2692.0, 763.4, 249.6, 90.28, 35.29, 14.51, 4.053, + 1.482, 0.2517, + ] as const, + c: [ + -7.23177e-05, -0.000555116, -0.00293805, -0.0117687, -0.0402907, -0.100609, + -0.196528, -0.102382, 0.52719, 0.593251, 0.0332652, + ] as const, +}; +export const CCPVDZ_SI_3S = { + alpha: [ + 78860.0, 11820.0, 2692.0, 763.4, 249.6, 90.28, 35.29, 14.51, 4.053, + 1.482, 0.2517, + ] as const, + c: [ + 1.85113e-05, 0.000142236, 0.000752185, 0.00302279, 0.0103677, 0.0262563, + 0.0523989, 0.0290959, -0.178003, -0.346874, 0.62302, + ] as const, +}; +export const CCPVDZ_SI_3S_P = { + alpha: [0.09243] as const, + c: [1.0] as const, +}; +export const CCPVDZ_SI_2P = { + alpha: [315.9, 74.42, 23.48, 8.488, 3.217, 1.229, 0.2964] as const, + c: [0.00392656, 0.0298811, 0.127212, 0.320943, 0.455429, 0.268563, 0.0188336] as const, +}; +export const CCPVDZ_SI_3P = { + alpha: [315.9, 74.42, 23.48, 8.488, 3.217, 1.229, 0.2964] as const, + c: [ + -0.000858302, -0.00630328, -0.0288255, -0.069456, -0.119493, -0.0199581, + 0.510268, + ] as const, +}; +export const CCPVDZ_SI_3P_P = { + alpha: [0.08768] as const, + c: [1.0] as const, +}; +export const CCPVDZ_SI_3D = { + alpha: [0.275] as const, + c: [1.0] as const, +}; +export const CCPVDZ_P_1S = { + alpha: [ + 94840.0, 14220.0, 3236.0, 917.1, 299.5, 108.1, 42.18, 17.28, 4.858, + 1.818, 0.3372, + ] as const, + c: [ + 0.000255509, 0.00198193, 0.010276, 0.0414823, 0.131984, 0.308662, 0.420647, + 0.222878, 0.0164035, -0.00254255, 0.00074805, + ] as const, +}; +export const CCPVDZ_P_2S = { + alpha: [ + 94840.0, 14220.0, 3236.0, 917.1, 299.5, 108.1, 42.18, 17.28, 4.858, + 1.818, 0.3372, + ] as const, + c: [ + -6.96939e-05, -0.000535266, -0.00283709, -0.0113983, -0.0392929, -0.0996364, + -0.197983, -0.11486, 0.518595, 0.601847, 0.0368612, + ] as const, +}; +export const CCPVDZ_P_3S = { + alpha: [ + 94840.0, 14220.0, 3236.0, 917.1, 299.5, 108.1, 42.18, 17.28, 4.858, + 1.818, 0.3372, + ] as const, + c: [ + 1.91199e-05, 0.000147223, 0.000777912, 0.00314546, 0.01082, 0.0279957, + 0.0563978, 0.035819, -0.193387, -0.372097, 0.624246, + ] as const, +}; +export const CCPVDZ_P_3S_P = { + alpha: [0.1232] as const, + c: [1.0] as const, +}; +export const CCPVDZ_P_2P = { + alpha: [370.5, 87.33, 27.59, 10.0, 3.825, 1.494, 0.3921] as const, + c: [0.00395005, 0.0302492, 0.129554, 0.327594, 0.456992, 0.253086, 0.0168798] as const, +}; +export const CCPVDZ_P_3P = { + alpha: [370.5, 87.33, 27.59, 10.0, 3.825, 1.494, 0.3921] as const, + c: [ + -0.000959832, -0.00711177, -0.0327122, -0.0795784, -0.135016, -0.00910585, + 0.537802, + ] as const, +}; +export const CCPVDZ_P_3P_P = { + alpha: [0.1186] as const, + c: [1.0] as const, +}; +export const CCPVDZ_P_3D = { + alpha: [0.373] as const, + c: [1.0] as const, +}; +export const CCPVDZ_S_1S = { + alpha: [ + 110800.0, 16610.0, 3781.0, 1071.0, 349.8, 126.3, 49.26, 20.16, 5.72, + 2.182, 0.4327, + ] as const, + c: [ + 0.000247635, 0.00192026, 0.00996192, 0.0402975, 0.128604, 0.30348, + 0.421432, 0.230781, 0.0178971, -0.00297516, 0.000849522, + ] as const, +}; +export const CCPVDZ_S_2S = { + alpha: [ + 110800.0, 16610.0, 3781.0, 1071.0, 349.8, 126.3, 49.26, 20.16, 5.72, + 2.182, 0.4327, + ] as const, + c: [ + -6.87039e-05, -0.000527681, -0.00279671, -0.0112651, -0.0388834, -0.0995025, + -0.19974, -0.12336, 0.513194, 0.60712, 0.0396753, + ] as const, +}; +export const CCPVDZ_S_3S = { + alpha: [ + 110800.0, 16610.0, 3781.0, 1071.0, 349.8, 126.3, 49.26, 20.16, 5.72, + 2.182, 0.4327, + ] as const, + c: [ + 1.99077e-05, 0.000153483, 0.000809503, 0.00328974, 0.0112967, 0.0296385, + 0.0599851, 0.0413248, -0.207474, -0.392889, 0.63284, + ] as const, +}; +export const CCPVDZ_S_3S_P = { + alpha: [0.157] as const, + c: [1.0] as const, +}; +export const CCPVDZ_S_2P = { + alpha: [399.7, 94.19, 29.75, 10.77, 4.119, 1.625, 0.4726] as const, + c: [0.00447541, 0.0341708, 0.14425, 0.353928, 0.459085, 0.206383, 0.0102141] as const, +}; +export const CCPVDZ_S_3P = { + alpha: [399.7, 94.19, 29.75, 10.77, 4.119, 1.625, 0.4726] as const, + c: [ + -0.00116251, -0.00865664, -0.0390886, -0.0934625, -0.147994, 0.0301904, + 0.561573, + ] as const, +}; +export const CCPVDZ_S_3P_P = { + alpha: [0.1407] as const, + c: [1.0] as const, +}; +export const CCPVDZ_S_3D = { + alpha: [0.479] as const, + c: [1.0] as const, +}; +export const CCPVDZ_CL_1S = { + alpha: [ + 127900.0, 19170.0, 4363.0, 1236.0, 403.6, 145.7, 56.81, 23.23, 6.644, + 2.575, 0.5371, + ] as const, + c: [ + 0.000241153, 0.00187095, 0.00970827, 0.0393153, 0.125932, 0.299341, + 0.421886, 0.237201, 0.0191531, -0.00334792, 0.000929883, + ] as const, +}; +export const CCPVDZ_CL_2S = { + alpha: [ + 127900.0, 19170.0, 4363.0, 1236.0, 403.6, 145.7, 56.81, 23.23, 6.644, + 2.575, 0.5371, + ] as const, + c: [ + -6.78922e-05, -0.000521836, -0.00276513, -0.0111537, -0.0385919, -0.0994848, + -0.201392, -0.130313, 0.509443, 0.610725, 0.0421549, + ] as const, +}; +export const CCPVDZ_CL_3S = { + alpha: [ + 127900.0, 19170.0, 4363.0, 1236.0, 403.6, 145.7, 56.81, 23.23, 6.644, + 2.575, 0.5371, + ] as const, + c: [ + 2.04986e-05, 0.000158298, 0.000833639, 0.0033988, 0.0116738, 0.0309622, + 0.0629533, 0.0460257, -0.219312, -0.408773, 0.638465, + ] as const, +}; +export const CCPVDZ_CL_3S_P = { + alpha: [0.1938] as const, + c: [1.0] as const, +}; +export const CCPVDZ_CL_2P = { + alpha: [417.6, 98.33, 31.04, 11.19, 4.249, 1.624, 0.5322] as const, + c: [ + 0.00525982, 0.0398332, 0.164655, 0.387322, 0.457072, 0.151636, 0.00181615, + ] as const, +}; +export const CCPVDZ_CL_3P = { + alpha: [417.6, 98.33, 31.04, 11.19, 4.249, 1.624, 0.5322] as const, + c: [ + -0.0014357, -0.0107796, -0.0470075, -0.11103, -0.153275, 0.0894609, + 0.579444, + ] as const, +}; +export const CCPVDZ_CL_3P_P = { + alpha: [0.162] as const, + c: [1.0] as const, +}; +export const CCPVDZ_CL_3D = { + alpha: [0.6] as const, + c: [1.0] as const, +}; +export const CCPVDZ_AR_1S = { + alpha: [ + 145700.0, 21840.0, 4972.0, 1408.0, 459.7, 165.9, 64.69, 26.44, 7.628, + 2.996, 0.6504, + ] as const, + c: [ + 0.0002367, 0.00183523, 0.0095286, 0.0386283, 0.124081, 0.296471, 0.422068, + 0.241711, 0.0200509, -0.00361, 0.000975607, + ] as const, +}; +export const CCPVDZ_AR_2S = { + alpha: [ + 145700.0, 21840.0, 4972.0, 1408.0, 459.7, 165.9, 64.69, 26.44, 7.628, + 2.996, 0.6504, + ] as const, + c: [ + -6.7491e-05, -0.000518522, -0.00274825, -0.0111007, -0.038482, -0.0997599, + -0.203088, -0.135608, 0.507195, 0.612898, 0.0442968, + ] as const, +}; +export const CCPVDZ_AR_3S = { + alpha: [ + 145700.0, 21840.0, 4972.0, 1408.0, 459.7, 165.9, 64.69, 26.44, 7.628, + 2.996, 0.6504, + ] as const, + c: [ + 2.10457e-05, 0.000162565, 0.000855463, 0.00349745, 0.0120156, 0.0321368, + 0.0655279, 0.049937, -0.229769, -0.421006, 0.642331, + ] as const, +}; +export const CCPVDZ_AR_3S_P = { + alpha: [0.2337] as const, + c: [1.0] as const, +}; +export const CCPVDZ_AR_2P = { + alpha: [453.7, 106.8, 33.73, 12.13, 4.594, 1.678, 0.5909] as const, + c: [ + 0.00570555, 0.043046, 0.176591, 0.406863, 0.452549, 0.122801, -0.00445996, + ] as const, +}; +export const CCPVDZ_AR_3P = { + alpha: [453.7, 106.8, 33.73, 12.13, 4.594, 1.678, 0.5909] as const, + c: [ + -0.00160655, -0.0121714, -0.0520789, -0.123737, -0.151619, 0.142425, + 0.584501, + ] as const, +}; +export const CCPVDZ_AR_3P_P = { + alpha: [0.1852] as const, + c: [1.0] as const, +}; +export const CCPVDZ_AR_3D = { + alpha: [0.738] as const, + c: [1.0] as const, +}; + +export const AUG_CCPVDZ_NA_DIFFUSE_S = { + alpha: [0.00725] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_NA_DIFFUSE_P = { + alpha: [0.00633] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_NA_DIFFUSE_D = { + alpha: [0.0468] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_MG_DIFFUSE_S = { + alpha: [0.01488] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_MG_DIFFUSE_P = { + alpha: [0.00935] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_MG_DIFFUSE_D = { + alpha: [0.0595] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_AL_DIFFUSE_S = { + alpha: [0.0231] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_AL_DIFFUSE_P = { + alpha: [0.0153] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_AL_DIFFUSE_D = { + alpha: [0.0535] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_SI_DIFFUSE_S = { + alpha: [0.0332] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_SI_DIFFUSE_P = { + alpha: [0.025] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_SI_DIFFUSE_D = { + alpha: [0.0823] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_P_DIFFUSE_S = { + alpha: [0.0417] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_P_DIFFUSE_P = { + alpha: [0.0343] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_P_DIFFUSE_D = { + alpha: [0.113] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_S_DIFFUSE_S = { + alpha: [0.0507] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_S_DIFFUSE_P = { + alpha: [0.0399] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_S_DIFFUSE_D = { + alpha: [0.152] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_CL_DIFFUSE_S = { + alpha: [0.0608] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_CL_DIFFUSE_P = { + alpha: [0.0466] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_CL_DIFFUSE_D = { + alpha: [0.196] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_AR_DIFFUSE_S = { + alpha: [0.0709] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_AR_DIFFUSE_P = { + alpha: [0.0533] as const, + c: [1.0] as const, +}; +export const AUG_CCPVDZ_AR_DIFFUSE_D = { + alpha: [0.24] as const, + c: [1.0] as const, +}; + // ── Generic shell type + multi-shell integrals ─────────────── // // A Shell is a contracted s-Gaussian: a center plus a list of diff --git a/src/chemistry/periodic-table.ts b/src/chemistry/periodic-table.ts index 1a146eb..2f153c9 100644 --- a/src/chemistry/periodic-table.ts +++ b/src/chemistry/periodic-table.ts @@ -1,6 +1,7 @@ // ───────────────────────────────────────────────────────────── // periodic-table.ts — atomic property database for the elements -// supported by webgpu-q (H, He, Li, Be, B, C, N, O, F, Ne). +// supported by webgpu-q (H, He, Li, Be, B, C, N, O, F, Ne, +// Na, Mg, Al, Si, P, S, Cl, Ar) — periods 1-3 complete. // // Data sources: // - Atomic number, isotope mass: NIST 2024. @@ -94,6 +95,54 @@ export const PERIODIC_TABLE: Readonly> = { covalentRadius: 0.67, vdwRadius: 1.54, period: 2, group: 18, family: "noble", }, + Na: { + symbol: "Na", name: "Sodium", atomicNumber: 11, + mass: 22.989769282, electronegativity: 0.93, + covalentRadius: 1.55, vdwRadius: 2.27, + period: 3, group: 1, family: "alkali", + }, + Mg: { + symbol: "Mg", name: "Magnesium", atomicNumber: 12, + mass: 23.985041697, electronegativity: 1.31, + covalentRadius: 1.39, vdwRadius: 1.73, + period: 3, group: 2, family: "alkaline-earth", + }, + Al: { + symbol: "Al", name: "Aluminium", atomicNumber: 13, + mass: 26.98153853, electronegativity: 1.61, + covalentRadius: 1.26, vdwRadius: 1.84, + period: 3, group: 13, family: "main", + }, + Si: { + symbol: "Si", name: "Silicon", atomicNumber: 14, + mass: 27.97692653465, electronegativity: 1.90, + covalentRadius: 1.16, vdwRadius: 2.10, + period: 3, group: 14, family: "main", + }, + P: { + symbol: "P", name: "Phosphorus", atomicNumber: 15, + mass: 30.97376199842, electronegativity: 2.19, + covalentRadius: 1.11, vdwRadius: 1.80, + period: 3, group: 15, family: "main", + }, + S: { + symbol: "S", name: "Sulfur", atomicNumber: 16, + mass: 31.9720711744, electronegativity: 2.58, + covalentRadius: 1.03, vdwRadius: 1.80, + period: 3, group: 16, family: "main", + }, + Cl: { + symbol: "Cl", name: "Chlorine", atomicNumber: 17, + mass: 34.968852682, electronegativity: 3.16, + covalentRadius: 0.99, vdwRadius: 1.75, + period: 3, group: 17, family: "halogen", + }, + Ar: { + symbol: "Ar", name: "Argon", atomicNumber: 18, + mass: 39.9623831237, electronegativity: null, + covalentRadius: 0.96, vdwRadius: 1.88, + period: 3, group: 18, family: "noble", + }, }; /** Look up element data. Throws on unknown symbol. */ diff --git a/src/chemistry/raman.ts b/src/chemistry/raman.ts index 61a4855..9741158 100644 --- a/src/chemistry/raman.ts +++ b/src/chemistry/raman.ts @@ -62,6 +62,14 @@ const ATOMIC_MASS: Readonly> = { O: 15.99491461956, F: 18.998403163, Ne: 19.9924401762, + Na: 22.989769282, + Mg: 23.985041697, + Al: 26.98153853, + Si: 27.97692653465, + P: 30.97376199842, + S: 31.9720711744, + Cl: 34.968852682, + Ar: 39.9623831237, }; export type EnergyMethod = "hf" | FunctionalKind; diff --git a/src/chemistry/thermochemistry.ts b/src/chemistry/thermochemistry.ts index 02b806c..cb17c86 100644 --- a/src/chemistry/thermochemistry.ts +++ b/src/chemistry/thermochemistry.ts @@ -83,6 +83,14 @@ const ATOMIC_MASS: Readonly> = { O: 15.99491461956, F: 18.998403163, Ne: 19.9924401762, + Na: 22.989769282, + Mg: 23.985041697, + Al: 26.98153853, + Si: 27.97692653465, + P: 30.97376199842, + S: 31.9720711744, + Cl: 34.968852682, + Ar: 39.9623831237, }; export interface ThermochemistryOpts { diff --git a/src/chemistry/vibrations.ts b/src/chemistry/vibrations.ts index ca4b0e5..8fb6f9b 100644 --- a/src/chemistry/vibrations.ts +++ b/src/chemistry/vibrations.ts @@ -47,7 +47,8 @@ const HA_BOHR2_AMU_TO_CM_INV = 5140.4865; const E2_PER_AMU_TO_KM_PER_MOL = 974.864; /** Atomic masses for the supported elements. Isotope-pure values - * (^1H, ^7Li, ^9Be, ^11B, ^12C, ^14N, ^16O, ^20Ne) — the standard + * (^1H, ^7Li, ^9Be, ^11B, ^12C, ^14N, ^16O, ^20Ne, ^23Na, ^24Mg, + * ^27Al, ^28Si, ^31P, ^32S, ^35Cl, ^40Ar) — the standard * choice for theoretical reference frequencies. */ const ATOMIC_MASS: Readonly> = { H: 1.00782503207, @@ -60,6 +61,14 @@ const ATOMIC_MASS: Readonly> = { O: 15.99491461956, F: 18.998403163, Ne: 19.9924401762, + Na: 22.989769282, + Mg: 23.985041697, + Al: 26.98153853, + Si: 27.97692653465, + P: 30.97376199842, + S: 31.9720711744, + Cl: 34.968852682, + Ar: 39.9623831237, }; export type EnergyMethod = "hf" | FunctionalKind; diff --git a/src/chemistry/xyz.ts b/src/chemistry/xyz.ts index 1ebf4ca..f4856b3 100644 --- a/src/chemistry/xyz.ts +++ b/src/chemistry/xyz.ts @@ -20,8 +20,9 @@ // Limitations: // - Single-frame only. Multi-frame trajectory XYZ (concatenated) // would need explicit frame splitting; queued. -// - Element set restricted to the AtomSymbol union (H, He, Li, Be, C, N, -// O, F). Throws for unsupported elements with a clear message. +// - Element set restricted to the AtomSymbol union: H through Ar +// (Z = 1-18), complete. Throws for anything heavier with a clear +// message naming the supported range. // - Comment line preserved as opt return; no parsing of charge / // multiplicity / properties (those vary too much across codes). // ───────────────────────────────────────────────────────────── @@ -34,9 +35,12 @@ import { type Atom, type AtomSymbol } from "./atoms.js"; const SYMBOL_BY_Z: Record = { 1: "H", 2: "He", 3: "Li", 4: "Be", 5: "B", 6: "C", 7: "N", 8: "O", 9: "F", 10: "Ne", + 11: "Na", 12: "Mg", 13: "Al", 14: "Si", + 15: "P", 16: "S", 17: "Cl", 18: "Ar", }; const SUPPORTED_SYMBOLS = new Set([ "H", "He", "Li", "Be", "B", "C", "N", "O", "F", "Ne", + "Na", "Mg", "Al", "Si", "P", "S", "Cl", "Ar", ]); export interface XYZParseResult { @@ -171,14 +175,14 @@ function resolveSymbol(token: string, lineNo: number): AtomSymbol { if (Number.isInteger(Z) && Z > 0) { const sym = SYMBOL_BY_Z[Z]; if (!sym) { - throw new Error(`parseXYZ: line ${lineNo} atomic number ${Z} is not supported (supported: H, He, Li, Be, B, C, N, O, F, Ne)`); + throw new Error(`parseXYZ: line ${lineNo} atomic number ${Z} is not supported (supported: H through Ar, Z = 1-18)`); } return sym; } // Otherwise treat as symbol — case-normalize and validate. const upper = token.charAt(0).toUpperCase() + token.slice(1).toLowerCase(); if (!SUPPORTED_SYMBOLS.has(upper as AtomSymbol)) { - throw new Error(`parseXYZ: line ${lineNo} element "${token}" not supported (supported: H, He, Li, Be, B, C, N, O, F, Ne)`); + throw new Error(`parseXYZ: line ${lineNo} element "${token}" not supported (supported: H through Ar, Z = 1-18)`); } return upper as AtomSymbol; } diff --git a/src/molecule/import-formats.ts b/src/molecule/import-formats.ts index 892ca54..b57f439 100644 --- a/src/molecule/import-formats.ts +++ b/src/molecule/import-formats.ts @@ -4,11 +4,12 @@ import type { Atom, AtomSymbol } from "../chemistry/atoms.js"; -// The chemistry stack currently supports periods 1-2 complete (see +// The chemistry stack currently supports periods 1-3 complete (see // src/chemistry/atoms.ts → AtomSymbol). Imports outside this set are // rejected with a clear error. const SUPPORTED_SYMBOLS = new Set([ "H", "He", "Li", "Be", "B", "C", "N", "O", "F", "Ne", + "Na", "Mg", "Al", "Si", "P", "S", "Cl", "Ar", ]); export interface ParseResult { @@ -83,17 +84,22 @@ function parsePDB(text: string): ParseResult { if (line.length < 54) continue; // Element column 77-78 (PDB v3.30); 1-indexed 77,78 = 0-indexed 76,77. let symRaw = line.slice(76, 78).trim(); + // The element column is authoritative; the atom-name column is a + // guess ("CA" = alpha-carbon, not calcium) and is the only place the + // one-character truncation in canonicalSymbol may be applied. + let fromAtomName = false; if (!symRaw) { // Fall back to first 1-2 alpha chars of atom name (cols 13-16). const name = line.slice(12, 16).trim(); symRaw = name.replace(/[^A-Za-z]/g, "").slice(0, 2); + fromAtomName = true; } if (!symRaw) continue; const x = parseFloat(line.slice(30, 38)); const y = parseFloat(line.slice(38, 46)); const z = parseFloat(line.slice(46, 54)); if (!Number.isFinite(x) || !Number.isFinite(y) || !Number.isFinite(z)) continue; - atoms.push({ symbol: canonicalSymbol(symRaw), pos: [x, y, z] }); + atoms.push({ symbol: canonicalSymbol(symRaw, fromAtomName), pos: [x, y, z] }); } if (atoms.length === 0) { throw new Error("PDB: no ATOM / HETATM records found"); @@ -135,15 +141,30 @@ function parseMOL(text: string): ParseResult { } // ── Symbol normalization ─────────────────────────────────────── -function canonicalSymbol(raw: string): AtomSymbol { +/** + * Normalize a raw token to a supported element symbol. + * + * `allowNameFallback` truncates a two-character token to one character + * when the pair isn't a supported symbol. That is ONLY correct for PDB + * *atom names* (cols 13-16), where "CA" means alpha-carbon, not calcium. + * + * It is actively dangerous anywhere else: the truncation is silent and + * yields a real-but-wrong element. "Fe" → "F" imports iron as fluorine. + * The hazard scales with the supported single-letter set, which this + * repo has grown from {H,B,C,N,O,F} to {H,B,C,N,O,F,P,S} — so "Pd"/"Pt" + * → P and "Sc"/"Se"/"Sn"/"Sr" → S all became newly wrong. ("Na" → N was + * already wrong before that.) Hence: opt-in, not default. + */ +function canonicalSymbol(raw: string, allowNameFallback = false): AtomSymbol { // "C" → "C", "ca" → "Ca", "12C" → "C", "Cα" → "C". const clean = raw.replace(/[^A-Za-z]/g, ""); if (!clean) throw new Error(`Unrecognized atom symbol "${raw}"`); const norm = clean.slice(0, 1).toUpperCase() + clean.slice(1, 2).toLowerCase(); - // Try the 2-char form first; if not in the supported set, fall back to 1-char. if (SUPPORTED_SYMBOLS.has(norm)) return norm as AtomSymbol; - const single = clean.slice(0, 1).toUpperCase(); - if (SUPPORTED_SYMBOLS.has(single)) return single as AtomSymbol; + if (allowNameFallback) { + const single = clean.slice(0, 1).toUpperCase(); + if (SUPPORTED_SYMBOLS.has(single)) return single as AtomSymbol; + } throw new Error(`Unsupported atom symbol "${raw}" (parsed "${norm}"). ` + `Supported: ${Array.from(SUPPORTED_SYMBOLS).join(", ")}`); } diff --git a/src/viz/molecular-graph-2d.ts b/src/viz/molecular-graph-2d.ts index bedf3c3..70f81ff 100644 --- a/src/viz/molecular-graph-2d.ts +++ b/src/viz/molecular-graph-2d.ts @@ -52,12 +52,13 @@ const ATOM_COLOR: Record = { P: "#ffb56e", S: "#fbbf24", Cl: "#34d399", + Ar: "#a78bfa", // noble — same hue as Ne, matching the group-inherits-colour pattern }; const ATOM_RADIUS: Record = { H: 6, He: 6, Li: 10, Be: 9, B: 9, C: 9, N: 9, O: 9, F: 9, Ne: 9, - Na: 12, Mg: 12, Al: 12, Si: 12, P: 12, S: 12, Cl: 11, + Na: 12, Mg: 12, Al: 12, Si: 12, P: 12, S: 12, Cl: 11, Ar: 11, }; export function moleculeGraph2D( diff --git a/tests/chemistry/elements/pyscf-reference.json b/tests/chemistry/elements/pyscf-reference.json index d321675..55d56a1 100644 --- a/tests/chemistry/elements/pyscf-reference.json +++ b/tests/chemistry/elements/pyscf-reference.json @@ -1956,6 +1956,1686 @@ "ok": true, "E_HF": -100.03442833006022, "nao": 34 + }, + { + "element": "Na", + "molecule": "NaH", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "Na", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.8874 + ] + } + ], + "ok": true, + "E_HF": -160.3021762604238, + "nao": 10 + }, + { + "element": "Na", + "molecule": "NaH", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Na", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.8874 + ] + } + ], + "ok": true, + "E_HF": -160.3021762604238, + "nao": 10 + }, + { + "element": "Na", + "molecule": "NaH", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Na", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.8874 + ] + } + ], + "ok": true, + "E_HF": -162.38394919200044, + "nao": 23 + }, + { + "element": "Na", + "molecule": "NaH", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Na", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.8874 + ] + } + ], + "ok": true, + "E_HF": -162.38410508781988, + "nao": 24 + }, + { + "element": "Na", + "molecule": "NaH", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Na", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.8874 + ] + } + ], + "ok": true, + "E_HF": -162.38455935886105, + "nao": 36 + }, + { + "element": "Na", + "molecule": "NaH", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Na", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.8874 + ] + } + ], + "ok": true, + "E_HF": -162.38473196292642, + "nao": 38 + }, + { + "element": "Mg", + "molecule": "MgH2", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "Mg", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.7297 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.7297 + ] + } + ], + "ok": true, + "E_HF": -198.19476027964825, + "nao": 11 + }, + { + "element": "Mg", + "molecule": "MgH2", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Mg", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.7297 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.7297 + ] + } + ], + "ok": true, + "E_HF": -198.19476027964825, + "nao": 11 + }, + { + "element": "Mg", + "molecule": "MgH2", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Mg", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.7297 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.7297 + ] + } + ], + "ok": true, + "E_HF": -200.72930433623605, + "nao": 28 + }, + { + "element": "Mg", + "molecule": "MgH2", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Mg", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.7297 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.7297 + ] + } + ], + "ok": true, + "E_HF": -200.7293760924933, + "nao": 29 + }, + { + "element": "Mg", + "molecule": "MgH2", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Mg", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.7297 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.7297 + ] + } + ], + "ok": true, + "E_HF": -200.73000650867198, + "nao": 45 + }, + { + "element": "Mg", + "molecule": "MgH2", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Mg", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.7297 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + -1.7297 + ] + } + ], + "ok": true, + "E_HF": -200.73010720998403, + "nao": 47 + }, + { + "element": "Al", + "molecule": "AlH3", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "Al", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 9.699202649247038e-17, + 1.584, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.3717842395945508, + -0.7920000000000003, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.3717842395945505, + -0.7920000000000007, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -240.71343276972328, + "nao": 12 + }, + { + "element": "Al", + "molecule": "AlH3", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Al", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 9.699202649247038e-17, + 1.584, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.3717842395945508, + -0.7920000000000003, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.3717842395945505, + -0.7920000000000007, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -240.71343276972328, + "nao": 12 + }, + { + "element": "Al", + "molecule": "AlH3", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Al", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 9.699202649247038e-17, + 1.584, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.3717842395945508, + -0.7920000000000003, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.3717842395945505, + -0.7920000000000007, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -243.6317557129841, + "nao": 33 + }, + { + "element": "Al", + "molecule": "AlH3", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Al", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 9.699202649247038e-17, + 1.584, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.3717842395945508, + -0.7920000000000003, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.3717842395945505, + -0.7920000000000007, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -243.63196936345872, + "nao": 34 + }, + { + "element": "Al", + "molecule": "AlH3", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Al", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 9.699202649247038e-17, + 1.584, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.3717842395945508, + -0.7920000000000003, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.3717842395945505, + -0.7920000000000007, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -243.63255345764384, + "nao": 54 + }, + { + "element": "Al", + "molecule": "AlH3", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Al", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 9.699202649247038e-17, + 1.584, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + -1.3717842395945508, + -0.7920000000000003, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 1.3717842395945505, + -0.7920000000000007, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -243.63264673280918, + "nao": 56 + }, + { + "element": "Si", + "molecule": "SiH4", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "Si", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + 0.8543629283468083, + 0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + -0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + 0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + -0.8543629283468083, + 0.8543629283468083 + ] + } + ], + "ok": true, + "E_HF": -287.9102282043203, + "nao": 13 + }, + { + "element": "Si", + "molecule": "SiH4", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Si", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + 0.8543629283468083, + 0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + -0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + 0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + -0.8543629283468083, + 0.8543629283468083 + ] + } + ], + "ok": true, + "E_HF": -287.9102282043203, + "nao": 13 + }, + { + "element": "Si", + "molecule": "SiH4", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Si", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + 0.8543629283468083, + 0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + -0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + 0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + -0.8543629283468083, + 0.8543629283468083 + ] + } + ], + "ok": true, + "E_HF": -291.24284859977655, + "nao": 38 + }, + { + "element": "Si", + "molecule": "SiH4", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Si", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + 0.8543629283468083, + 0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + -0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + 0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + -0.8543629283468083, + 0.8543629283468083 + ] + } + ], + "ok": true, + "E_HF": -291.2430035069773, + "nao": 39 + }, + { + "element": "Si", + "molecule": "SiH4", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Si", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + 0.8543629283468083, + 0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + -0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + 0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + -0.8543629283468083, + 0.8543629283468083 + ] + } + ], + "ok": true, + "E_HF": -291.2431798099662, + "nao": 63 + }, + { + "element": "Si", + "molecule": "SiH4", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Si", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + 0.8543629283468083, + 0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + 0.8543629283468083, + -0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + 0.8543629283468083, + -0.8543629283468083 + ] + }, + { + "symbol": "H", + "pos": [ + -0.8543629283468083, + -0.8543629283468083, + 0.8543629283468083 + ] + } + ], + "ok": true, + "E_HF": -291.2433039442788, + "nao": 65 + }, + { + "element": "P", + "molecule": "PH3", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "P", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.312926304492749e-17, + 1.1942914985094957, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0342867772330082, + -0.597145749254748, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + 1.034286777233008, + -0.5971457492547484, + -0.7681587183570484 + ] + } + ], + "ok": true, + "E_HF": -338.63361401512805, + "nao": 12 + }, + { + "element": "P", + "molecule": "PH3", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "P", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.312926304492749e-17, + 1.1942914985094957, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0342867772330082, + -0.597145749254748, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + 1.034286777233008, + -0.5971457492547484, + -0.7681587183570484 + ] + } + ], + "ok": true, + "E_HF": -338.63361401512805, + "nao": 12 + }, + { + "element": "P", + "molecule": "PH3", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "P", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.312926304492749e-17, + 1.1942914985094957, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0342867772330082, + -0.597145749254748, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + 1.034286777233008, + -0.5971457492547484, + -0.7681587183570484 + ] + } + ], + "ok": true, + "E_HF": -342.47042849794116, + "nao": 33 + }, + { + "element": "P", + "molecule": "PH3", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "P", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.312926304492749e-17, + 1.1942914985094957, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0342867772330082, + -0.597145749254748, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + 1.034286777233008, + -0.5971457492547484, + -0.7681587183570484 + ] + } + ], + "ok": true, + "E_HF": -342.47054747932475, + "nao": 34 + }, + { + "element": "P", + "molecule": "PH3", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "P", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.312926304492749e-17, + 1.1942914985094957, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0342867772330082, + -0.597145749254748, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + 1.034286777233008, + -0.5971457492547484, + -0.7681587183570484 + ] + } + ], + "ok": true, + "E_HF": -342.4720623161136, + "nao": 54 + }, + { + "element": "P", + "molecule": "PH3", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "P", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 7.312926304492749e-17, + 1.1942914985094957, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + -1.0342867772330082, + -0.597145749254748, + -0.7681587183570484 + ] + }, + { + "symbol": "H", + "pos": [ + 1.034286777233008, + -0.5971457492547484, + -0.7681587183570484 + ] + } + ], + "ok": true, + "E_HF": -342.4721254986151, + "nao": 56 + }, + { + "element": "S", + "molecule": "H2S", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "S", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + }, + { + "symbol": "H", + "pos": [ + -0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + } + ], + "ok": true, + "E_HF": -394.3115557633474, + "nao": 11 + }, + { + "element": "S", + "molecule": "H2S", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "S", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + }, + { + "symbol": "H", + "pos": [ + -0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + } + ], + "ok": true, + "E_HF": -394.3115557633474, + "nao": 11 + }, + { + "element": "S", + "molecule": "H2S", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "S", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + }, + { + "symbol": "H", + "pos": [ + -0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + } + ], + "ok": true, + "E_HF": -398.6945380060537, + "nao": 28 + }, + { + "element": "S", + "molecule": "H2S", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "S", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + }, + { + "symbol": "H", + "pos": [ + -0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + } + ], + "ok": true, + "E_HF": -398.69474102043057, + "nao": 29 + }, + { + "element": "S", + "molecule": "H2S", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "S", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + }, + { + "symbol": "H", + "pos": [ + -0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + } + ], + "ok": true, + "E_HF": -398.69803596335436, + "nao": 45 + }, + { + "element": "S", + "molecule": "H2S", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "S", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + }, + { + "symbol": "H", + "pos": [ + -0.9616404049011908, + 0.0, + 0.9268630382432312 + ] + } + ], + "ok": true, + "E_HF": -398.69820133904415, + "nao": 47 + }, + { + "element": "Cl", + "molecule": "HCl", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "Cl", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.2746 + ] + } + ], + "ok": true, + "E_HF": -455.1348097021576, + "nao": 10 + }, + { + "element": "Cl", + "molecule": "HCl", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Cl", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.2746 + ] + } + ], + "ok": true, + "E_HF": -455.1348097021576, + "nao": 10 + }, + { + "element": "Cl", + "molecule": "HCl", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Cl", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.2746 + ] + } + ], + "ok": true, + "E_HF": -460.08944519174275, + "nao": 23 + }, + { + "element": "Cl", + "molecule": "HCl", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Cl", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.2746 + ] + } + ], + "ok": true, + "E_HF": -460.089733039711, + "nao": 24 + }, + { + "element": "Cl", + "molecule": "HCl", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Cl", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.2746 + ] + } + ], + "ok": true, + "E_HF": -460.0926154947685, + "nao": 36 + }, + { + "element": "Cl", + "molecule": "HCl", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Cl", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + }, + { + "symbol": "H", + "pos": [ + 0.0, + 0.0, + 1.2746 + ] + } + ], + "ok": true, + "E_HF": -460.0929435882411, + "nao": 38 + }, + { + "element": "Ar", + "molecule": "Ar", + "basis": "sto-3g", + "convention": "spherical", + "atoms": [ + { + "symbol": "Ar", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -521.2228808035303, + "nao": 9 + }, + { + "element": "Ar", + "molecule": "Ar", + "basis": "sto-3g", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Ar", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -521.2228808035303, + "nao": 9 + }, + { + "element": "Ar", + "molecule": "Ar", + "basis": "cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Ar", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -526.7998653097462, + "nao": 18 + }, + { + "element": "Ar", + "molecule": "Ar", + "basis": "cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Ar", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -526.800080054873, + "nao": 19 + }, + { + "element": "Ar", + "molecule": "Ar", + "basis": "aug-cc-pvdz", + "convention": "spherical", + "atoms": [ + { + "symbol": "Ar", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -526.8009724025385, + "nao": 27 + }, + { + "element": "Ar", + "molecule": "Ar", + "basis": "aug-cc-pvdz", + "convention": "cartesian", + "atoms": [ + { + "symbol": "Ar", + "pos": [ + 0.0, + 0.0, + 0.0 + ] + } + ], + "ok": true, + "E_HF": -526.8013706297443, + "nao": 29 } ] } \ No newline at end of file diff --git a/tests/chemistry/periodic-table.test.ts b/tests/chemistry/periodic-table.test.ts index 5a5c039..b5f8d39 100644 --- a/tests/chemistry/periodic-table.test.ts +++ b/tests/chemistry/periodic-table.test.ts @@ -12,8 +12,10 @@ describe("Periodic table", () => { expect(PERIODIC_TABLE.F.atomicNumber).toBe(9); expect(PERIODIC_TABLE.B.atomicNumber).toBe(5); expect(PERIODIC_TABLE.Ne.atomicNumber).toBe(10); - // H, He, Li, Be, B, C, N, O, F, Ne — periods 1-2 complete. - expect(Object.keys(PERIODIC_TABLE).length).toBe(10); + expect(PERIODIC_TABLE.Na.atomicNumber).toBe(11); + expect(PERIODIC_TABLE.Ar.atomicNumber).toBe(18); + // H…Ne plus Na, Mg, Al, Si, P, S, Cl, Ar — periods 1-3 complete. + expect(Object.keys(PERIODIC_TABLE).length).toBe(18); }); test("Pauling electronegativity ordering: H < C < N < O < F", () => { diff --git a/tests/molecule/import-formats.test.ts b/tests/molecule/import-formats.test.ts index 0bf7e5e..9cecdb8 100644 --- a/tests/molecule/import-formats.test.ts +++ b/tests/molecule/import-formats.test.ts @@ -32,10 +32,41 @@ H -0.757160 0.586260 0.000000`; }); test("XYZ — rejects unsupported atom", () => { - const xyz = `1\nsi atom\nSi 0 0 0`; + // Titanium: period 4, far outside anything the basis tables cover. This + // test previously used Si, which became supported when period 3 landed. + // The negative example must ALSO not begin with a supported one-letter + // symbol — canonicalSymbol falls back to the 1-char form, so "Fe" would + // silently parse as F and "Sc"/"Pd" as S/P. "T" is not a symbol, so Ti + // reaches the throw. (tests/chemistry/xyz.test.ts uses Kr for the same + // role; a distinct element keeps the two failures distinguishable.) + const xyz = `1\nti atom\nTi 0 0 0`; expect(() => parseGeometry(xyz, "xyz")).toThrow(/Unsupported atom symbol/); }); + test("XYZ — an unsupported 2-char symbol never truncates to a supported 1-char one", () => { + // Regression guard. canonicalSymbol used to fall back from the 2-char + // form to the 1-char form unconditionally, so "Fe" imported as + // FLUORINE — a real, wrong, silent element substitution. The hazard + // grew when period 3 made P and S supported single letters. The + // fallback is now opt-in and reserved for PDB atom names. + for (const sym of ["Fe", "Pd", "Pt", "Sc", "Se", "Sn", "Sr", "Ni"]) { + const xyz = `1\ntest\n${sym} 0 0 0`; + expect(() => parseGeometry(xyz, "xyz"), `${sym} must throw, not truncate`) + .toThrow(/Unsupported atom symbol/); + } + }); + + test("PDB — atom-name fallback still resolves CA to carbon", () => { + // The one place the 1-char truncation IS correct: PDB cols 13-16 carry + // atom names, where "CA" means alpha-carbon rather than calcium. The + // element column (77-78) is left blank here to force that path. + const pdb = +`ATOM 1 CA ALA A 1 1.000 2.000 3.000 1.00 0.00`; + const r = parseGeometry(pdb, "pdb"); + expect(r.atoms).toHaveLength(1); + expect(r.atoms[0]!.symbol).toBe("C"); + }); + test("PDB — H₂O from cols", () => { // Standard PDB atom-record format. Element field at cols 77-78. const pdb = From 14d0cee4553da553a73bb72702d217cbbbb83c24 Mon Sep 17 00:00:00 2001 From: Ahmet Gunaydin Date: Mon, 10 Aug 2026 10:24:46 +0700 Subject: [PATCH 07/13] docs: run report + correct LIMITATIONS coverage; test the curriculum unlock MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Phase 3 of docs/RUN-PLAN-24H-ELEMENTS.md. The report leads with what could NOT be verified, per the plan, because on a run this size that section is worth more than the list of what shipped. Unverified, in short: Level A of the two-level scheme was never built (only Level B exists, so wrong-data and wrong-integrals are not separable); its 1e-9 Ha bar is unreachable anyway because boys0 uses an Abramowitz-Stegun erf fit with a measured 2.66e-7 relative floor; all 108 reference cells are RHF only, so no correlated method is validated for any new element; the row-3 frozen-core path (5 orbitals) is implemented but untested and interacts with a known throw at rhf-auto.ts:452; gradients were not hardened; and reference geometries cannot be validated by a test that feeds identical coordinates to both engines. LIMITATIONS.md basis coverage rewritten. It claimed 6-31G* was "wired, spot-checked" for H/C/N/O — the string "6-31g" appears nowhere in src/ and BasisName has only three members. That claim was never true, and it was sitting in the document whose entire job is honest disclosure. The section also cited ccpvdz-firstrow.test.ts as its precision evidence at "~10 mHa"; that test's real slack is 116-170 mHa and it failed to catch Li/Be being wrong by 1.29 mHa. Precision claims now rest on the element agreement test instead. curriculum-molecules.test.ts checks what the run was actually for: CH3Cl (the SN2 substrate), CH3SH, H2S and PH3 all reproduce PySCF to better than 0.1 mHa at STO-3G, in ~1.2 s combined. Element coverage is only interesting if it unlocks chemistry people teach. Co-Authored-By: Claude Opus 5 (1M context) --- LIMITATIONS.md | 48 ++++--- docs/RUN-REPORT-2026-08-10.md | 135 ++++++++++++++++++ .../elements/curriculum-molecules.test.ts | 101 +++++++++++++ 3 files changed, 268 insertions(+), 16 deletions(-) create mode 100644 docs/RUN-REPORT-2026-08-10.md create mode 100644 tests/chemistry/elements/curriculum-molecules.test.ts diff --git a/LIMITATIONS.md b/LIMITATIONS.md index aa8b3a8..d9587d9 100644 --- a/LIMITATIONS.md +++ b/LIMITATIONS.md @@ -43,24 +43,40 @@ reviewer or chemist would discover anyway. ### Basis-set atom coverage +**Periods 1–3 complete as of 2026-08-10: H through Ar (Z = 1–18).** + | basis | atoms wired | |---|---| -| STO-3G | H, Li, Be, C, N, O, F | -| 6-31G* | H, C, N, O (spot-checked) | -| **cc-pVDZ** | **H, Li, Be, C, N, O, F** (full first-row coverage — Tier 3 shipped 2026-05) | -| aug-cc-pVDZ | **H, Li, Be, C, N, O, F** (diffuse tables for full first row wired 2026-05) | - -LiH / BeH₂ / CH₄ / NH₃ / HF now first-class cc-pVDZ targets. -Verified by `tests/chemistry/ccpvdz-firstrow.test.ts` — each -molecule's cc-pVDZ HF energy converges and lies variationally -below its STO-3G counterpart. Tolerance is ~10 mHa vs PySCF 2.13.0 -reference values (loose, because the test is verifying basis -wiring not SCF precision; the existing H₂O cc-pVDZ tests cover -the precision case to 35 µHa). - -Aug-cc-pVDZ diffuse functions for Li, Be, C, N, F → ~30 minutes -each from the same EMSL source, queued as a follow-up if needed -for anions or excited-state work on those systems. +| **STO-3G** | **H–Ar** (all 18) | +| **cc-pVDZ** | **H–Ar** (all 18) | +| **aug-cc-pVDZ** | **H–Ar** (all 18) | + +Nothing beyond Z = 18. No transition metals, no fourth row. +**6-31G\* is NOT implemented** — `BasisName` is exactly +`"sto-3g" | "cc-pvdz" | "aug-cc-pvdz"` and the string "6-31g" appears +nowhere in `src/`. (An earlier revision of this table claimed 6-31G* +was "wired, spot-checked" for H/C/N/O. That was never true.) + +Every (element × basis) cell is validated against PySCF by +`tests/chemistry/elements/reference-agreement.test.ts`: 18 elements × +3 bases × both d-conventions, one real molecule each (hydrides +wherever one exists, so two-center integrals are actually exercised), +**bar 0.1 mHa**. `tests/chemistry/elements/pyscf-reference.json` holds +the committed reference table; regenerate it with +`scripts/run-element-reference.py`. + +The d-convention matters from Na onward as well as across Li–Ne: our +Cartesian path must be compared against PySCF with `mol.cart = True` +and our spherical path against PySCF's default, or a ~0.34 mHa phantom +error appears on every element carrying d functions. STO-3G has no d +functions below Z = 19, so it is convention-free throughout. + +Superseded caveat, kept for provenance: this section previously +described `ccpvdz-firstrow.test.ts` as the verification, at "~10 mHa" +tolerance. That test's real slack was 116–170 mHa, and it did not +detect that the Li and Be cc-pVDZ tables were wrong by up to 1.29 mHa +(fixed 2026-08-10). Precision claims here now rest on the element +agreement test, not on it. --- diff --git a/docs/RUN-REPORT-2026-08-10.md b/docs/RUN-REPORT-2026-08-10.md new file mode 100644 index 0000000..0c5952f --- /dev/null +++ b/docs/RUN-REPORT-2026-08-10.md @@ -0,0 +1,135 @@ +# Run report — periodic-table expansion, 2026-08-10 + +Branch: `feat/periodic-table-expansion`. Nothing pushed, nothing deployed. + +--- + +## 1. What could NOT be verified + +Listed first deliberately — this is the most useful section. + +### 1.1 The Level-A bar in the plan is unreachable, and was never run + +`docs/RUN-PLAN-24H-ELEMENTS.md` Gate 0.3 specified a two-level scheme: +Level A (our primitives fed to PySCF via `gto.basis.parse`, bar 1e-9 Ha) +to isolate engine bugs, Level B (our tables vs PySCF's) to isolate data +errors. + +**Only Level B was built.** The committed harness +(`tests/chemistry/elements/reference-agreement.test.ts`) compares our +engine against PySCF's own basis at 0.1 mHa. That catches wrong data and +wrong integrals together but cannot separate them. + +Level A's 1e-9 Ha bar is also unreachable as specified. `boys0` in +`integrals-cg.ts` uses the Abramowitz–Stegun 7.1.26 rational fit for +erf, with a measured worst-case relative error of **2.66e-7** (the code +admits ~1.5e-7 at `integrals-cg.ts:195`). Row 3 makes it worse: Ar's +STO-3G 1s exponent is 674.45 against carbon's 71.62, pushing more +core-core pairs onto the approximate branch. A real Level A needs either +an f64-accurate `erf` or a bar set from measurement. + +### 1.2 Correlated methods are unvalidated for every new element + +All 108 committed reference cells are **RHF only**. MP2, CCSD, CCSD(T), +DFT and EOM have no per-element reference for B, Ne, or any of Na–Ar. +Wrong basis data would show up in HF, so this is not blind — but +"HF agrees" is not "CCSD(T) agrees". + +### 1.3 The row-3 frozen-core path is untested + +`FROZEN_CORE_FOR` now returns 5 for Na–Ar (the full neon core) and +`defaultFrozenCore` sums it. No test exercises it. +`tests/chemistry/frozen-core-audit.test.ts` hardcodes `nFrozenCore: 1` +and only uses H2O/CH4, where 1 is correct — so it passes and will keep +passing while the row-3 path stays uncovered. + +Known interaction: `rhf-auto.ts:452` throws for `nFrozenCore > 0` on the +exact-ERI UMP2 path. Any caller that inherits 5 from `defaultFrozenCore` +on a small row-3 molecule will hit that throw. Not triggered by anything +currently in the suite. + +### 1.4 Gradients were not hardened (Phase 2 not started) + +No analytic-vs-finite-difference check was run for any new element. + +### 1.5 No geometry validation + +Reference geometries for the new hydrides (NaH 1.8874, MgH2 1.7297, +AlH3 1.5840, SiH4 1.4798, PH3 1.42/93.5, H2S 1.3356/92.11, HCl 1.2746, +BH3 1.19 A) were taken as experimental-ish values and used identically +on both sides. Since PySCF and webgpu-q consume the same coordinates, +a wrong geometry cannot cause a test failure — it just means the cell +validates a molecule that isn't quite the real one. + +--- + +### 1.6 The 6-31G* claim was false, and had been for a long time + +`LIMITATIONS.md` listed 6-31G* as "wired, spot-checked" for H/C/N/O. +The string "6-31g" appears nowhere in `src/`; `BasisName` is exactly +`"sto-3g" | "cc-pvdz" | "aug-cc-pvdz"`. Corrected in this run. Worth +noting *where* it was found: the honest-limitations document was itself +carrying a false capability claim. + +--- + +## 2. What landed + +18 elements, H through Ar, complete. 109/109 agreement cells green +(18 elements x 3 bases x 2 d-conventions, plus a convention sanity +check), bar 0.1 mHa. Full suite 140 files / 1065 passed / 1 skipped. + +The point of the exercise was curriculum reach, so that is checked +directly in `tests/chemistry/elements/curriculum-molecules.test.ts`: + +| molecule | why it matters | vs PySCF | +|---|---|---| +| CH3Cl | the SN2 substrate — Cl was the blocking element | < 0.1 mHa | +| CH3SH | real organosulfur, not a toy hydride | < 0.1 mHa | +| H2S | the "why is H2S a gas and H2O a liquid" comparison | < 0.1 mHa | +| PH3 | phosphorus; entry point to phosphates and DNA | < 0.1 mHa | + +All four converge in ~1.2 s combined at STO-3G. + +--- + +## 3. Defects found in existing code + +See commit messages for full detail. Summary: + +| defect | worst measured impact | +|---|---| +| Li cc-pVDZ 2p coefficients + diffuse exponent | 1.293 mHa on LiH | +| Be cc-pVDZ 1s/2s/2p coefficients | 0.367 mHa on BeH2 | +| aug-cc-pVDZ diffuse: only H and C were correct of 8 | 0.171 mHa on HF | +| CCPVDZ_H exponents non-canonical | ~5 uHa, every H molecule | +| CCPVDZ_N_1S coefficients | ~0.03 uHa | +| BECKE_XI missing F; He/C/O off-by-one | 0.64 uHa | +| import-formats: "Fe" silently imported as fluorine | wrong element, silent | + +--- + +## 4. Corrections to the plan itself + +- Na and Mg carry d functions in cc-pVDZ, so the Cartesian/spherical + convention difference applies from **Na**, not from Al as written. +- STO-3G row 3 is **9** Cartesian functions (3s + 2p), not 13. +- The plan said "fix the generator — never the existing constants". + Inverted: the constants were wrong and PySCF was the authority. + +--- + +## 5. What the next run should start with + +1. Phase 2 — gradient hardening (analytic vs FD) for the 10 new elements. +2. Correlated-method references (MP2/CCSD/CCSD(T)) per element. +3. A row-3 frozen-core test, plus resolving the `rhf-auto.ts:452` throw. +4. The loose-assertion cleanup, starting with + `e2e/swarm-hf-anthracene-ccpvdz.spec.ts`, whose `(-1500, -100)` window + certifies a documented-wrong -880 Ha result as green. +5. Re-source the ~15 hardcoded literals from + `experiments/results/2026-07-06/level-6/E34-pyscf.json`, which holds + 40 full-precision references that no test currently uses. +6. Wire `scripts/check-basis-vs-pyscf.py` into CI — it is the only + artifact that catches bad basis digits at the source, and it runs + from no workflow today. diff --git a/tests/chemistry/elements/curriculum-molecules.test.ts b/tests/chemistry/elements/curriculum-molecules.test.ts new file mode 100644 index 0000000..d455451 --- /dev/null +++ b/tests/chemistry/elements/curriculum-molecules.test.ts @@ -0,0 +1,101 @@ +// The molecules period 3 actually buys. +// +// Element coverage is only interesting if it unlocks chemistry people +// teach. These four are the payoff for adding Cl, S and P: +// +// CH3Cl — the SN2 substrate. Chloromethane is the canonical example +// in every introductory organic course; Cl was the single +// element blocking it. +// CH3SH — methanethiol. Real organosulfur, not a toy hydride. +// H2S — the sulfur analogue of water; the classic "why is H2S a gas +// and H2O a liquid" comparison. +// PH3 — phosphine; the entry point to phosphates and the DNA +// backbone. +// +// Each is checked against PySCF at STO-3G, Cartesian convention, on +// identical geometries. Bar is the same 0.1 mHa used by the element +// agreement test. +import { describe, expect, test } from "vitest"; +import { computeMolecularIntegrals } from "../../../src/chemistry/cg-molecular.js"; +import { moleculeToShellsNuclei, type Atom } from "../../../src/chemistry/atoms.js"; +import { runRHFSCF } from "../../../src/chemistry/hf-scf.js"; + +const D = Math.PI / 180; +const tet = (b: number) => { + const a = b / Math.sqrt(3); + return [[a, a, a], [a, -a, -a], [-a, a, -a], [-a, -a, a]] as const; +}; + +function ch3cl(): Atom[] { + const th = 108.0 * D, ch = 1.087; + const out: Atom[] = [ + { symbol: "C", pos: [0, 0, 0] }, + { symbol: "Cl", pos: [0, 0, 1.785] }, + ]; + for (let k = 0; k < 3; k++) { + const phi = 120 * k * D; + const r = ch * Math.sin(Math.PI - th); + const z = ch * Math.cos(Math.PI - th); + out.push({ symbol: "H", pos: [r * Math.cos(phi), r * Math.sin(phi), z] }); + } + return out; +} + +function ch3sh(): Atom[] { + const out: Atom[] = [ + { symbol: "S", pos: [0, 0, 0] }, + { symbol: "C", pos: [0, 0, 1.819] }, + { symbol: "H", pos: [1.336 * Math.sin(96.5 * D), 0, -1.336 * Math.cos(96.5 * D)] }, + ]; + for (let k = 0; k < 3; k++) { + const phi = 120 * k * D; + out.push({ symbol: "H", pos: [0.63 * Math.cos(phi), 0.63 * Math.sin(phi), 2.18] }); + } + return out; +} + +function bent(sym: Atom["symbol"], bond: number, angleDeg: number): Atom[] { + const h = (angleDeg / 2) * D; + const x = bond * Math.sin(h), z = bond * Math.cos(h); + return [ + { symbol: sym, pos: [0, 0, 0] }, + { symbol: "H", pos: [x, 0, z] }, + { symbol: "H", pos: [-x, 0, z] }, + ]; +} + +function pyramidal(sym: Atom["symbol"], bond: number, angleDeg: number): Atom[] { + const ha = (angleDeg / 2) * D; + const r = bond * Math.sin(ha) / Math.sin(60 * D); + const h = Math.sqrt(Math.max(bond * bond - r * r, 0)); + const out: Atom[] = [{ symbol: sym, pos: [0, 0, 0] }]; + for (let k = 0; k < 3; k++) { + const a = (90 + 120 * k) * D; + out.push({ symbol: "H", pos: [r * Math.cos(a), r * Math.sin(a), -h] }); + } + return out; +} + +// PySCF RHF/STO-3G, mol.cart = True, conv_tol 1e-12. +const CASES: { name: string; atoms: Atom[]; nao: number; ref: number }[] = [ + { name: "CH3Cl (SN2 substrate)", atoms: ch3cl(), nao: 17, ref: -493.509047901 }, + { name: "CH3SH (methanethiol)", atoms: ch3sh(), nao: 18, ref: -432.209032619 }, + { name: "H2S", atoms: bent("S", 1.3356, 92.11), nao: 11, ref: -394.311555763 }, + { name: "PH3 (phosphine)", atoms: pyramidal("P", 1.42, 93.5), nao: 12, ref: -338.633614015 }, +]; + +describe("Curriculum molecules unlocked by period 3", () => { + for (const c of CASES) { + test(`${c.name}: RHF/STO-3G within 0.1 mHa of PySCF`, () => { + const { shells, nuclei, nElectrons } = moleculeToShellsNuclei(c.atoms, "sto-3g"); + const ints = computeMolecularIntegrals(shells, nuclei); + expect(ints.n).toBe(c.nao); + const hf = runRHFSCF(ints, nElectrons, { + useDIIS: true, maxIter: 300, energyTol: 1e-11, densityTol: 1e-9, + }); + // A non-converged SCF is a failure, never a datum. + expect(hf.converged).toBe(true); + expect(Math.abs((hf.energy - c.ref) * 1000)).toBeLessThan(0.1); + }, 300_000); + } +}); From 7de9b74d4ca5bb03927cccb9bbb6526203b070dd Mon Sep 17 00:00:00 2001 From: Ahmet Gunaydin Date: Mon, 10 Aug 2026 10:56:11 +0700 Subject: [PATCH 08/13] =?UTF-8?q?feat:=20labs.html=20=E2=80=94=20three=20c?= =?UTF-8?q?omp-chem=20teaching=20labs,=20computed=20live=20in=20the=20tab?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit The education deliverable. Three standard undergraduate computational chemistry labs, each running real SCF / MP2 / CCSD / CCSD(T) in a Web Worker on the reader's machine. Nothing is precomputed and no number is hard-coded into the page; the interpretation text is generated FROM the results, so a student reads their own run rather than a promise about one. Lab 1 Bond dissociation — RHF vs exact FCI on H2 across 14 bond lengths. The teaching moment is c_g / c_u going 0.994 / -0.113 at equilibrium to 0.707 / -0.707 at 5 A: two configurations contributing equally is static correlation, and it is why single-reference theory fails there and why CASSCF exists. Lab 2 Basis-set convergence — H2O at STO-3G / cc-pVDZ / aug-cc-pVDZ. Minimal to double-zeta buys 1064 mHa; adding diffuse functions buys 14.8 mHa for 18 more basis functions. Diminishing returns, measured rather than asserted. Lab 3 The correlation ladder — HF to MP2 to CCSD to CCSD(T) on H2O. MP2 recovers 71.7% of the correlation CCSD(T) finds, CCSD 99.7%. STO-3G runs in seconds; cc-pVDZ is selectable and is the calculation a methods course actually assigns. Verified in a real browser against the built bundle, not just in unit tests. Lab 3's STO-3G numbers (MP2 -74.998420, CCSD -75.012287) match the repo's own full-precision PySCF artifact (experiments/results/2026-07-06/level-6/E34-pyscf.json: -74.99842043, -75.01228678) exactly. Lab 2's aug-cc-pVDZ value is -76.041983 against PySCF's -76.0419832 — that agreement is new this session; before the oxygen diffuse-function repair it read -76.041842. The worker treats a non-converged SCF or CCSD amplitude set as an ERROR and the page says so explicitly. That is deliberate: an earlier attempt at Lab 1 printed unconverged CCSD numbers as data, and a teaching page that does this teaches the opposite of what it should. Honest scoping, stated on the page: the exact curve in Lab 1 is a 16 determinant full CI specific to H2 in STO-3G. There is no general FCI solver here, so this lab does not generalise to N2. Co-Authored-By: Claude Opus 5 (1M context) --- labs.html | 166 +++++++++++++++ src/labs/lab-worker.ts | 146 +++++++++++++ src/labs/main.ts | 200 ++++++++++++++++++ .../elements/curriculum-molecules.test.ts | 4 - vite.config.ts | 1 + 5 files changed, 513 insertions(+), 4 deletions(-) create mode 100644 labs.html create mode 100644 src/labs/lab-worker.ts create mode 100644 src/labs/main.ts diff --git a/labs.html b/labs.html new file mode 100644 index 0000000..5dac853 --- /dev/null +++ b/labs.html @@ -0,0 +1,166 @@ + + + + + + + + Three comp-chem labs, run in your tab · webgpu-q + + + + + + + + + + + + + + + + + + +
+ + +

Three labs, computed in this tab

+

+ The standard opening sequence of any computational chemistry course. + Nothing on this page is precomputed — pressing Run starts a real SCF, + MP2, CCSD and CCSD(T) in a Web Worker on your machine. +

+ +
+ Why this is not a simulation of a lab. These are the same + algorithms, at the same precision, that a cluster job would run — + validated against PySCF to better than 0.1 mHa. What is missing is + throughput, not correctness: a laptop tab does small molecules. + For a teaching lab, small molecules are the entire syllabus. +
+ +
+

Lab 1 · Why single-reference methods break

+

Pull an H₂ molecule apart. Where does Hartree–Fock stop working, and why exactly there?

+

+ Two curves: restricted Hartree–Fock, and the exact answer from a full + CI in the same basis. They agree near equilibrium. Watch the gap open, + and watch the two configuration coefficients c_g and c_u while it does. +

+
+ + 14 bond lengths, ~1 s +
+
+
+
+ +
+

Lab 2 · What a bigger basis actually buys

+

Same molecule, same method, three basis sets. How fast do the returns diminish?

+

+ Water at fixed experimental geometry, Hartree–Fock throughout. Only + the basis changes: minimal (STO-3G), double-zeta (cc-pVDZ), then + double-zeta plus diffuse functions (aug-cc-pVDZ). +

+
+ + 3 SCF calculations, ~5 s +
+
+
+
+ +
+

Lab 3 · The correlation ladder

+

Hartree–Fock is missing something. How much, and what does each rung of theory buy back?

+

+ Water again, climbing HF → MP2 → CCSD → CCSD(T). STO-3G runs in + seconds. cc-pVDZ is the real thing and takes a couple of minutes — + it is the calculation a methods course actually assigns. +

+
+ + + + ready +
+
+
+
+ +
+ Real Hartree–Fock, MP2, CCSD and CCSD(T) — in your browser, no install, no backend. + Powered by webgpu-q. +
+
+ + + diff --git a/src/labs/lab-worker.ts b/src/labs/lab-worker.ts new file mode 100644 index 0000000..51d8956 --- /dev/null +++ b/src/labs/lab-worker.ts @@ -0,0 +1,146 @@ +// ───────────────────────────────────────────────────────────── +// labs/lab-worker.ts — off-main-thread compute for the three +// teaching labs. Each lab streams partial results back as they +// land so the page fills in row by row rather than freezing and +// then dumping everything at the end. +// +// A non-converged calculation is reported as an ERROR, never as a +// number. This is the whole pedagogical point of Lab 1: CCSD stops +// converging once the HOMO–LUMO gap closes, and a lab that silently +// printed those amplitudes as data would teach the opposite of what +// it should. +// ───────────────────────────────────────────────────────────── + +import type { Atom, BasisName } from "../chemistry/atoms.js"; +import { moleculeToShellsNuclei } from "../chemistry/atoms.js"; +import { computeMolecularIntegrals } from "../chemistry/cg-molecular.js"; +import { runRHFSCF } from "../chemistry/hf-scf.js"; +import { runMP2 } from "../chemistry/mp2.js"; +import { runCCSD } from "../chemistry/ccsd.js"; +import { runCCSDT } from "../chemistry/ccsd-t.js"; +import { fciState } from "../viz/h2-fci-state.js"; + +export type LabRequest = + | { kind: "dissociation"; seq: number } + | { kind: "basis"; seq: number } + | { kind: "ladder"; basis: BasisName; seq: number }; + +export type LabReply = + | { kind: "row"; seq: number; lab: string; row: Record } + | { kind: "done"; seq: number; lab: string; seconds: number } + | { kind: "error"; seq: number; lab: string; message: string }; + +const post = (m: LabReply) => (self as unknown as Worker).postMessage(m); + +const SCF = { useDIIS: true, maxIter: 300, energyTol: 1e-11, densityTol: 1e-9 } as const; + +function water(): Atom[] { + const h = (104.52 / 2) * Math.PI / 180; + const x = 0.9572 * Math.sin(h), z = 0.9572 * Math.cos(h); + return [ + { symbol: "O", pos: [0, 0, 0] }, + { symbol: "H", pos: [x, 0, z] }, + { symbol: "H", pos: [-x, 0, z] }, + ]; +} + +function rhf(atoms: Atom[], basis: BasisName) { + const { shells, nuclei, nElectrons } = moleculeToShellsNuclei(atoms, basis); + const integrals = computeMolecularIntegrals(shells, nuclei); + const hf = runRHFSCF(integrals, nElectrons, SCF); + if (!hf.converged) throw new Error(`RHF did not converge (${basis})`); + return { hf, integrals, n: integrals.n }; +} + +// ── Lab 1: H₂ dissociation, RHF vs exact ───────────────────── +// The exact curve comes from `fciState`, a 16-dimensional full CI in +// the H₂ STO-3G minimal basis. That is H₂-SPECIFIC (it builds a dense +// 4-qubit Hamiltonian); it is not a general FCI solver and this lab +// does not generalise to N₂ or anything else. +function labDissociation(seq: number) { + const t0 = performance.now(); + const Rs = [0.4, 0.5, 0.6, 0.7414, 0.9, 1.1, 1.3, 1.6, 2.0, 2.5, 3.0, 3.5, 4.0, 5.0]; + for (const R of Rs) { + const exact = fciState(R); + let eRHF = Number.NaN; + let rhfOk = false; + try { + const { hf } = rhf([ + { symbol: "H", pos: [0, 0, 0] }, + { symbol: "H", pos: [0, 0, R] }, + ], "sto-3g"); + eRHF = hf.energy; + rhfOk = true; + } catch { /* recorded as non-converged, not as a number */ } + post({ + kind: "row", seq, lab: "dissociation", + row: { + R, eExact: exact.energy, eRHF, rhfOk, + cG: exact.cG, cU: exact.cU, + errMHa: rhfOk ? (eRHF - exact.energy) * 1000 : Number.NaN, + }, + }); + } + post({ kind: "done", seq, lab: "dissociation", seconds: (performance.now() - t0) / 1000 }); +} + +// ── Lab 2: basis-set convergence on H₂O ────────────────────── +function labBasis(seq: number) { + const t0 = performance.now(); + const atoms = water(); + let prev: number | null = null; + for (const basis of ["sto-3g", "cc-pvdz", "aug-cc-pvdz"] as const) { + const { hf, n } = rhf(atoms, basis); + post({ + kind: "row", seq, lab: "basis", + row: { + basis, n, energy: hf.energy, + gainMHa: prev === null ? Number.NaN : (hf.energy - prev) * 1000, + }, + }); + prev = hf.energy; + } + post({ kind: "done", seq, lab: "basis", seconds: (performance.now() - t0) / 1000 }); +} + +// ── Lab 3: the correlation ladder on H₂O ───────────────────── +function labLadder(basis: BasisName, seq: number) { + const t0 = performance.now(); + const atoms = water(); + const { hf, integrals, n } = rhf(atoms, basis); + post({ kind: "row", seq, lab: "ladder", row: { method: "RHF", energy: hf.energy, corrMHa: 0, pct: 0, n } }); + + const mp2 = runMP2(hf, integrals); + const ccsd = runCCSD(hf, integrals); + if (!ccsd.converged) throw new Error("CCSD amplitudes did not converge"); + const t = runCCSDT(ccsd, hf, integrals); + + // CCSD(T) is the reference for "100% of recovered correlation" here — + // it is the best number this ladder produces, not the exact answer. + const total = t.totalEnergy - hf.energy; + const rows: [string, number][] = [ + ["MP2", mp2.totalEnergy], ["CCSD", ccsd.totalEnergy], ["CCSD(T)", t.totalEnergy], + ]; + for (const [method, energy] of rows) { + const corr = energy - hf.energy; + post({ + kind: "row", seq, lab: "ladder", + row: { method, energy, corrMHa: corr * 1000, pct: (corr / total) * 100, n }, + }); + } + post({ kind: "done", seq, lab: "ladder", seconds: (performance.now() - t0) / 1000 }); +} + +self.addEventListener("message", (ev: MessageEvent) => { + const req = ev.data; + try { + if (req.kind === "dissociation") labDissociation(req.seq); + else if (req.kind === "basis") labBasis(req.seq); + else labLadder(req.basis, req.seq); + } catch (e) { + post({ + kind: "error", seq: req.seq, lab: req.kind, + message: e instanceof Error ? e.message : String(e), + }); + } +}); diff --git a/src/labs/main.ts b/src/labs/main.ts new file mode 100644 index 0000000..fa686ff --- /dev/null +++ b/src/labs/main.ts @@ -0,0 +1,200 @@ +// ───────────────────────────────────────────────────────────── +// labs/main.ts — three standard undergraduate computational +// chemistry labs, each computed live in the tab. +// +// Nothing here is precomputed or hard-coded: pressing Run starts a +// real SCF / MP2 / CCSD / CCSD(T) in a worker and the table fills in +// as results arrive. The interpretation text is revealed only after +// the numbers exist, so a student reads their own result rather than +// a promise about it. +// ───────────────────────────────────────────────────────────── + +import type { LabReply, LabRequest } from "./lab-worker.js"; +import type { BasisName } from "../chemistry/atoms.js"; + +type Row = Record; + +const worker = new Worker(new URL("./lab-worker.ts", import.meta.url), { type: "module" }); +let seq = 0; +const pending = new Map void; + onDone: (s: number) => void; + onError: (m: string) => void; +}>(); + +worker.addEventListener("message", (ev: MessageEvent) => { + const m = ev.data; + const h = pending.get(m.seq); + if (!h) return; + if (m.kind === "row") h.onRow(m.row); + else if (m.kind === "done") { h.onDone(m.seconds); pending.delete(m.seq); } + else { h.onError(m.message); pending.delete(m.seq); } +}); + +function run( + req: Omit, + onRow: (r: Row) => void, + onDone: (s: number) => void, + onError: (m: string) => void, +) { + const s = ++seq; + pending.set(s, { onRow, onDone, onError }); + worker.postMessage({ ...req, seq: s } as LabRequest); +} + +const $ = (id: string) => document.getElementById(id)!; +const num = (v: unknown, d = 6) => + typeof v === "number" && Number.isFinite(v) ? v.toFixed(d) : "—"; + +function wire( + labId: string, + /** Built at click time, so controls (e.g. the basis dropdown) are read fresh. */ + reqFor: () => Omit, + header: string[], + toCells: (r: Row) => string[], + interpret: (rows: Row[]) => string, +) { + const btn = $(`${labId}-run`) as HTMLButtonElement; + const status = $(`${labId}-status`); + const tbody = $(`${labId}-body`); + const thead = $(`${labId}-head`); + const note = $(`${labId}-note`); + + btn.addEventListener("click", () => { + btn.disabled = true; + btn.textContent = "Running…"; + status.textContent = "computing — this is a real calculation, not a lookup"; + status.className = "status busy"; + tbody.innerHTML = ""; + note.innerHTML = ""; + thead.innerHTML = `${header.map((h) => `${h}`).join("")}`; + const rows: Row[] = []; + + run( + reqFor(), + (r) => { + rows.push(r); + const tr = document.createElement("tr"); + tr.innerHTML = toCells(r).map((c) => `${c}`).join(""); + tbody.appendChild(tr); + }, + (s) => { + btn.disabled = false; + btn.textContent = "Run again"; + status.textContent = `done in ${s.toFixed(1)} s`; + status.className = "status ok"; + note.innerHTML = interpret(rows); + }, + (m) => { + btn.disabled = false; + btn.textContent = "Run"; + status.textContent = `stopped: ${m}`; + status.className = "status err"; + note.innerHTML = + `

The calculation refused to return a number. That is the + library behaving correctly — an unconverged SCF or CCSD amplitude set + is not a result, and reporting it as one is how wrong numbers end up + in papers.

`; + }, + ); + }); +} + +// ── Lab 1 ──────────────────────────────────────────────────── +wire( + "diss", + () => ({ kind: "dissociation" }), + ["R (Å)", "E_exact (Ha)", "E_RHF (Ha)", "error (mHa)", "c_g", "c_u"], + (r) => [ + num(r["R"], 4), + num(r["eExact"]), + r["rhfOk"] ? num(r["eRHF"]) : 'no convergence', + r["rhfOk"] ? num(r["errMHa"], 2) : "—", + num(r["cG"], 3), + num(r["cU"], 3), + ], + (rows) => { + const far = rows[rows.length - 1]; + const eq = rows.find((r) => Math.abs(Number(r["R"]) - 0.7414) < 1e-6); + if (!far || !eq) return ""; + return ` +

Read the last row. At R = ${num(far["R"], 1)} Å the two + coefficients have gone to c_g ≈ ${num(far["cG"], 3)} and + c_u ≈ ${num(far["cU"], 3)} — equal magnitude. Near equilibrium + (${num(eq["R"], 4)} Å) they were ${num(eq["cG"], 3)} and + ${num(eq["cU"], 3)}.

+

That crossover is static correlation. A single Slater + determinant assumes one configuration dominates. Once two contribute + equally, no amount of single-reference theory can fix it — which is why + RHF's error grows to ${num(far["errMHa"], 0)} mHa at dissociation while + being nearly exact at equilibrium. It is also the reason CASSCF exists.

+

The exact curve here is a full CI in the H₂ STO-3G + minimal basis — 16 determinants, solved directly. That is specific to + H₂; this page does not contain a general FCI solver, so the same lab + cannot be run on N₂.

`; + }, +); + +// ── Lab 2 ──────────────────────────────────────────────────── +wire( + "basis", + () => ({ kind: "basis" }), + ["basis", "functions", "E_HF (Ha)", "gain (mHa)"], + (r) => [ + String(r["basis"]), + String(r["n"]), + num(r["energy"]), + Number.isFinite(Number(r["gainMHa"])) ? num(r["gainMHa"], 2) : "—", + ], + (rows) => { + if (rows.length < 3) return ""; + const [a, b, c] = rows as [Row, Row, Row]; + return ` +

Diminishing returns, quantified. Going minimal → double-zeta + bought ${num(b["gainMHa"], 0)} mHa. Adding diffuse functions on top + bought only ${num(c["gainMHa"], 1)} mHa — roughly + ${Math.round(Number(b["gainMHa"]) / Number(c["gainMHa"]))}× less, for + ${Number(c["n"]) - Number(b["n"])} more basis functions.

+

Every energy is variational, so each row must lie below the one above + it. It does. But notice none of them is the exact answer: the basis-set + limit and the correlation problem are two different errors, and this lab + only attacks the first. Lab 3 attacks the second.

+

${a["basis"]} at n = ${a["n"]} is a minimal basis — + one function per occupied atomic orbital. It is useful for teaching and + almost never for publishing.

`; + }, +); + +// ── Lab 3 ──────────────────────────────────────────────────── +const ladderBasis = () => + ((document.getElementById("ladder-basis") as HTMLSelectElement).value) as BasisName; + +wire( + "ladder", + () => ({ kind: "ladder", basis: ladderBasis() }), + ["method", "E (Ha)", "E_corr (mHa)", "% of CCSD(T)"], + (r) => [ + String(r["method"]), + num(r["energy"]), + Number(r["corrMHa"]) === 0 ? "0" : num(r["corrMHa"], 2), + Number(r["corrMHa"]) === 0 ? "—" : `${num(r["pct"], 1)}%`, + ], + (rows) => { + const mp2 = rows.find((r) => r["method"] === "MP2"); + const ccsd = rows.find((r) => r["method"] === "CCSD"); + const t = rows.find((r) => r["method"] === "CCSD(T)"); + if (!mp2 || !ccsd || !t) return ""; + return ` +

The ladder, in one table. MP2 — the cheapest correlated method + — already recovers ${num(mp2["pct"], 1)}% of what CCSD(T) finds. CCSD + gets to ${num(ccsd["pct"], 1)}%. The perturbative triples supply the + last ${num(100 - Number(ccsd["pct"]), 1)}%.

+

Total correlation energy here is ${num(t["corrMHa"], 1)} mHa. For + scale, chemical accuracy is 1.594 mHa — so the correlation you are + recovering is over a hundred times larger than the accuracy you need, + which is why you cannot simply ignore it.

+

CCSD(T) is the reference for the percentages, not the + exact answer. It is the best number this ladder produces; the true + correlation energy in this basis is slightly larger still.

`; + }, +); diff --git a/tests/chemistry/elements/curriculum-molecules.test.ts b/tests/chemistry/elements/curriculum-molecules.test.ts index d455451..c8033a9 100644 --- a/tests/chemistry/elements/curriculum-molecules.test.ts +++ b/tests/chemistry/elements/curriculum-molecules.test.ts @@ -21,10 +21,6 @@ import { moleculeToShellsNuclei, type Atom } from "../../../src/chemistry/atoms. import { runRHFSCF } from "../../../src/chemistry/hf-scf.js"; const D = Math.PI / 180; -const tet = (b: number) => { - const a = b / Math.sqrt(3); - return [[a, a, a], [a, -a, -a], [-a, a, -a], [-a, -a, a]] as const; -}; function ch3cl(): Atom[] { const th = 108.0 * D, ch = 1.087; diff --git a/vite.config.ts b/vite.config.ts index 05c5203..0fb0eb8 100644 --- a/vite.config.ts +++ b/vite.config.ts @@ -38,6 +38,7 @@ export default defineConfig({ viz: resolve(__dirname, "viz.html"), molecule: resolve(__dirname, "molecule.html"), learn: resolve(__dirname, "learn.html"), + labs: resolve(__dirname, "labs.html"), swarm: resolve(__dirname, "swarm.html"), screening: resolve(__dirname, "screening.html"), gpuMps: resolve(__dirname, "experiments/gpu-mps/index.html"), From 2e05ff8635fa47fe793294140469dcf6d3152894 Mon Sep 17 00:00:00 2001 From: Ahmet Gunaydin Date: Mon, 10 Aug 2026 10:56:57 +0700 Subject: [PATCH 09/13] feat: link labs.html from the site navigation A page nothing links to does not exist. Added to both nav bars in index.html and to the learn and viz navs. Co-Authored-By: Claude Opus 5 (1M context) --- index.html | 2 ++ learn.html | 1 + viz.html | 1 + 3 files changed, 4 insertions(+) diff --git a/index.html b/index.html index c07cf45..4723614 100644 --- a/index.html +++ b/index.html @@ -650,6 +650,7 @@ hyperscope molecule learn + labs swarm screen experiments @@ -1149,6 +1150,7 @@

Open the URL.
Run the science.

hyperscope molecule learn + labs swarm screen experiments diff --git a/learn.html b/learn.html index 2bc6cd0..8176852 100644 --- a/learn.html +++ b/learn.html @@ -102,6 +102,7 @@