Metadata-Version: 2.4
Name: nebcert
Version: 1.2.0
Summary: Quality checks for NEB reaction paths, transition states and tunnelling-corrected rate constants (ORCA NEB/NEB-TS outputs)
Author-email: Andrés Monreal-Hernández <andres.monreal@ues.mx>
License: MIT
Project-URL: Homepage, https://github.com/sircalch/nebcert
Project-URL: Documentation, https://github.com/sircalch/nebcert#readme
Project-URL: Repository, https://github.com/sircalch/nebcert.git
Project-URL: Bug Tracker, https://github.com/sircalch/nebcert/issues
Project-URL: Changelog, https://github.com/sircalch/nebcert/blob/main/CHANGELOG.md
Keywords: chemical-kinetics,transition-state-theory,nudged-elastic-band,neb,ci-neb,reaction-pathway,minimum-energy-path,quantum-tunneling,eckart-barrier,wigner-correction,irc,vasp,orca,gaussian,computational-chemistry,reproducibility
Classifier: Development Status :: 5 - Production/Stable
Classifier: Intended Audience :: Science/Research
Classifier: Topic :: Scientific/Engineering :: Chemistry
Classifier: Topic :: Scientific/Engineering :: Physics
Classifier: License :: OSI Approved :: MIT License
Classifier: Programming Language :: Python :: 3
Classifier: Programming Language :: Python :: 3.9
Classifier: Programming Language :: Python :: 3.10
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Operating System :: OS Independent
Requires-Python: >=3.9
Description-Content-Type: text/markdown
License-File: LICENSE
Requires-Dist: numpy>=1.20.0
Requires-Dist: scipy>=1.7.0
Requires-Dist: pandas>=1.3.0
Requires-Dist: matplotlib>=3.5.0
Requires-Dist: jinja2>=3.0.0
Provides-Extra: dev
Requires-Dist: pytest>=7.0.0; extra == "dev"
Requires-Dist: pytest-cov>=4.0.0; extra == "dev"
Requires-Dist: black>=23.0.0; extra == "dev"
Requires-Dist: flake8>=6.0.0; extra == "dev"
Requires-Dist: build>=0.10.0; extra == "dev"
Requires-Dist: twine>=4.0.0; extra == "dev"
Dynamic: license-file

# NEBCert

[![CI](https://github.com/sircalch/nebcert/actions/workflows/test.yml/badge.svg)](https://github.com/sircalch/nebcert/actions)
[![License: MIT](https://img.shields.io/badge/License-MIT-yellow.svg)](https://opensource.org/licenses/MIT)
[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.22217586.svg)](https://doi.org/10.5281/zenodo.22217586)

Quality checks for nudged-elastic-band (NEB) reaction paths, transition states and tunnelling-corrected
rate constants. NEBCert reads an ORCA NEB / NEB-CI / NEB-TS output (or a table of image energies), checks the
band and the transition state, and computes Eyring rate constants with Wigner and Eckart tunnelling factors.

A PASS means that the checks that were run found nothing wrong. It does not show that the path is the
minimum-energy path, that the transition state connects the intended minima, or that the barrier is accurate
at the chosen level of theory.

## What it checks

| Check | PASS | WARNING | FAIL |
|---|---|---|---|
| Band (NEB) | converged band, no interior minimum deeper than 0.5 kcal/mol | force > 0.05 eV/Å at the image nearest the barrier (skipped after a converged TS optimisation); interior minimum (more than one step) | band or TS optimisation not converged; end points identical; flat profile |
| Transition state | exactly one imaginary mode with \|ν\| ≥ 50 cm⁻¹ and no other imaginary mode | extra small imaginary modes (\|ν\| < 50 cm⁻¹); only small ones; IRC reported as failed | no imaginary mode; two or more with \|ν\| ≥ 50 cm⁻¹ |
| Eyring rate | barrier given as a Gibbs energy of activation | — | — |
| Tunnelling | T = 298.15 K above the crossover temperature T_c = hc\|ν\|/(2πk_B) | 298.15 K below T_c: Wigner invalid, one-dimensional Eckart only an estimate | barrier height ≤ 0 |

With an electronic barrier (no zero-point, thermal or entropic terms, as from a band) the Eyring number is
reported as NOT_APPLICABLE: it is not a TST rate constant, and it does not count towards the overall status,
which is the worst of the individual statuses. NEBCert does not run or read IRC calculations;
`--irc` records the user's statement that an IRC connected the intended minima.

For ORCA outputs NEBCert reads the last `PATH SUMMARY` table, the climbing image, the NEB and TS convergence
banners, the `PATH SUMMARY FOR NEB-TS` transition state and the final frequencies. When an optimised NEB-TS
transition state is present its energy is the barrier; otherwise the energy of a converged climbing image,
and otherwise the maximum of a natural cubic spline through the images. A table printed before the images were computed (E = 0, as left by a job whose image calculations
failed) is rejected.

The Eckart factor is κ(T) = exp(V₁/kT)/kT ∫ P(E) exp(−E/kT) dE with the transmission probability of the
asymmetric Eckart barrier (Eckart 1930; Johnston and Heicklen 1962), the barrier width fixed by the imaginary
frequency, and the integral evaluated by adaptive quadrature (scipy `quad`).

## Validation (`validation/`)

All scripts and results are in `validation/`; every number below is regenerated by them.

- **Tunnelling and rate constants** (`compare_overreact.py`, `eckart_highprecision.py`): on 720 combinations of
  imaginary frequency (200–3000 cm⁻¹), forward and reverse barriers (2–40 and 2–60 kcal/mol) and temperature
  (200–1000 K), κ_Eckart agrees with a 40-digit mpmath reference to 1.7 × 10⁻⁹ (relative). overreact deviates
  from the same reference by up to 1.3 × 10⁻², above 10⁻³ at 73 points, all with κ < 1.9. κ_Wigner and the
  Eyring rate agree with overreact to 3 × 10⁻¹¹ and 1.3 × 10⁻⁹.
- **Version 1.0.0** (the only release published before 1.2.0; `legacy100_grid.py`): its κ_Eckart differs from the
  reference by more than a factor of 2 at 529 of the 720 points (median ratio 6.8, range 3 × 10⁻⁸ to 2 × 10⁴³).
- **Real ORCA 6.1.1 calculations** (`make_neb_inputs.py`, `benchmark_neb.py`; B3LYP-D3(BJ)/def2-SVP; results in
  `validation/results/neb_benchmark.csv`, figure `validation/figures/fig3_profiles.png`). Eight reactions with
  NEB-TS and frequencies (H + HCl, OH + H₂, CH₃ + H₂, Cl⁻ + CH₃Cl, HCOH → H₂CO, H₂CO → H₂ + CO, HCN → HNC,
  NH₃ inversion), three degraded runs and one failed job.
  - Values read from the `.out` file agree with ORCA's full-precision files: band energies within 0.006 kcal/mol
    (`.final.interp`; the output table has five decimals), TS energies within 0.003 kcal/mol
    (`_NEB-TS_converged.xyz`), frequencies within 0.005 cm⁻¹ (`.hess`). κ_Eckart for the eight transition states
    agrees with the mpmath reference to 3 × 10⁻⁹.
  - Verdicts: PASS for six reactions; WARNING for HCOH → H₂CO and H₂CO → H₂ + CO, whose imaginary frequencies
    (2060 and 1880 cm⁻¹) put 298.15 K below the crossover temperature (κ_Eckart = 1.8 × 10⁵ and 7.6 × 10⁶).
    FAIL for a band stopped after three iterations (its highest image, 45.8 kcal/mol, against a converged barrier
    of 2.3 kcal/mol) and for a band whose end points collapsed into the same structure; the failed job is
    rejected. A 3-image NEB-CI passes, with the climbing-image barrier equal to the 8-image NEB-TS barrier
    (84.36 kcal/mol; the spline maximum is 85.62).
  - Version 1.1.0 warned on every reaction (the `--irc` bug), passed the failed job and crashed on the collapsed
    band. Version 1.0.0 passed the failed job, crashed on the collapsed band, and its κ_Eckart was 1.6 to
    6 × 10⁴⁶ times the reference for the eight transition states.

### Not validated
- The VASP `neb.dat` and Gaussian IRC parsers have not been tested on real output; the Gaussian IRC parser is
  not used by the command-line interface.
- Paths with several steps are flagged, not split.

## Installation

```bash
pip install "git+https://github.com/sircalch/nebcert@v1.2.0"
```

From source:

```bash
git clone https://github.com/sircalch/nebcert.git
cd nebcert
pip install -e .[dev]
```

## Command line

```bash
# ORCA NEB / NEB-CI / NEB-TS output (frequencies are read from a NEB-TS ... Freq run)
nebcert assess -i reaction.out --reaction "H + HCl -> H2 + Cl" --functional "B3LYP-D3(BJ)/def2-SVP" -o report/

# Gibbs energy of activation for the Eyring rate constant
nebcert assess -i reaction.out --gibbs-barrier-kcal 12.4 -o report/

# Table of image energies plus frequencies
nebcert assess -i neb.dat --frequencies "-1250,150,300,800" -o report/

# Synthetic demonstration (not a real calculation)
nebcert demo -o demo/
```

The report directory contains `report.html`, figures, a summary table (CSV and LaTeX), a draft methods
paragraph that states the results of the checks, and `citation.bib`.

## Python

```python
from nebcert import (calculate_neb_profile_analysis, verify_ts_frequency_and_irc,
                     calculate_eyring_tst_rates, calculate_quantum_tunneling_corrections,
                     assess_reaction_pathway_quality)
from nebcert.parsers.orca_neb import parse_orca_neb_output

d = parse_orca_neb_output("reaction.out")
ts = d["ts"]
neb = calculate_neb_profile_analysis(d["energies_ev"], d["coordinates_s_ang"], d["tangent_forces"],
                                     band_converged=d["neb_converged"],
                                     ts_energy_rel_ev=ts["barrier_fwd_ev"] if ts else None,
                                     ts_optimisation_converged=ts["converged"] if ts else None)
freq = verify_ts_frequency_and_irc(d["frequencies"])
tst = calculate_eyring_tst_rates(neb.e_forward_barrier_ev, barrier_type="electronic")
tun = calculate_quantum_tunneling_corrections(abs(freq.imaginary_frequency_cm1),
                                              neb.e_forward_barrier_ev, neb.e_reverse_barrier_ev)
report = assess_reaction_pathway_quality({"reaction": "H + HCl -> H2 + Cl"}, neb, freq, tst, tun)
print(report.overall_status, report.validation_score)
for message in report.recommendations:
    print("-", message)
```

## Citation

```bash
nebcert cite
```

prints the APA and BibTeX entries for the installed version (concept DOI 10.5281/zenodo.22217586).

## License

MIT, see [LICENSE](LICENSE).
