From 1aac8a966f39f62041c664dbbb06b77a1058ec0c Mon Sep 17 00:00:00 2001 From: OutisLi_Bot <198997215+OutisLi-Bot@users.noreply.github.com> Date: Thu, 15 Oct 2026 16:51:49 +0800 Subject: [PATCH 1/8] feat(lammps): add Cartesian magnetic-force mode for dpa4spin Offer an explicit magnetic_force precession|cartesian selection on pair_style dpa4spin and dpa4spin/kk so Newtonian spin integrators can consume -dE/dS in atom->fm without |m|/hbar scaling. Default remains precession. Expose the choice via extract("fm_cartesian_force"). Fixes #6077 --- doc/model/dpa4c.md | 14 ++ examples/spin/dpa4c/lmp/README.md | 6 + examples/spin/dpa4c/lmp/in.lammps | 6 + source/lmp/pair_dpa4spin.cpp | 146 +++++++++++----- source/lmp/pair_dpa4spin.h | 14 +- source/lmp/pair_dpa4spin_kokkos.cpp | 26 +-- ...test_dpa4spin_magnetic_force_convention.py | 163 ++++++++++++++++++ 7 files changed, 320 insertions(+), 55 deletions(-) create mode 100644 source/lmp/tests/test_dpa4spin_magnetic_force_convention.py diff --git a/doc/model/dpa4c.md b/doc/model/dpa4c.md index e75f286b54..df66b340d0 100644 --- a/doc/model/dpa4c.md +++ b/doc/model/dpa4c.md @@ -398,6 +398,20 @@ already performs, and the magnetic force is reduced back onto owning atoms alongside the conservative force, so domain decomposition needs no additional exchange. +By default the pair styles store the precession magnetic force expected by the +LAMMPS SPIN package (the model energy gradient scaled by the moment magnitude +over `hbar`). A Newtonian Cartesian integrator instead needs the raw energy +gradient `-dE/dS`. Select it with `magnetic_force cartesian` (requires +`units metal`); the default is `magnetic_force precession`: + +```lammps +pair_style dpa4spin/kk compressed_model.pt2 magnetic_force cartesian +``` + +The selected convention is exposed as the scalar integer +`Pair::extract("fm_cartesian_force", dim)` (zero for precession, one for +Cartesian) so a consuming integrator can reject an incompatible configuration. + `min_style spin` reads the magnetic force from the pair style and relaxes the moment directions. Spin dynamics through `fix nve/spin` requires a LAMMPS build whose fix recognizes this pair style, because the stock fix accumulates the diff --git a/examples/spin/dpa4c/lmp/README.md b/examples/spin/dpa4c/lmp/README.md index aa4153f4f6..42924a6c0b 100644 --- a/examples/spin/dpa4c/lmp/README.md +++ b/examples/spin/dpa4c/lmp/README.md @@ -67,3 +67,9 @@ they stay fixed through that stage. Spin dynamics through `fix nve/spin` requires a LAMMPS build whose fix recognizes this pair style, because the stock fix accumulates the magnetic force only from pair styles matching its own name pattern. + +The default magnetic-force convention is precession (`|m| / hbar` scaling of +the model energy gradient). A separately maintained Newtonian Cartesian +integrator can request the unconverted gradient with +`magnetic_force cartesian` on the `pair_style` line (`units metal`) and can +confirm the choice through `Pair::extract("fm_cartesian_force", dim)`. diff --git a/examples/spin/dpa4c/lmp/in.lammps b/examples/spin/dpa4c/lmp/in.lammps index 050ce53903..7c0d283ba5 100644 --- a/examples/spin/dpa4c/lmp/in.lammps +++ b/examples/spin/dpa4c/lmp/in.lammps @@ -11,6 +11,12 @@ # graph and the moment in device memory: # # lmp -k on g 1 -sf kk -in in.lammps +# +# By default atom->fm holds the precession magnetic force. For a Newtonian +# Cartesian integrator that needs -dE/dS, append `magnetic_force cartesian` +# (units metal) and read Pair::extract("fm_cartesian_force", dim): +# +# pair_style dpa4spin/kk frozen_model.pt2 magnetic_force cartesian units metal dimension 3 diff --git a/source/lmp/pair_dpa4spin.cpp b/source/lmp/pair_dpa4spin.cpp index 0a8d6e1ba8..916f72323b 100644 --- a/source/lmp/pair_dpa4spin.cpp +++ b/source/lmp/pair_dpa4spin.cpp @@ -107,9 +107,12 @@ static const char cite_user_deepmd_package[] = "}\n\n"; namespace { -// Reduced Planck constant in eV.ps. The model reports the magnetic force as -// the energy gradient with respect to the magnetic moment, while LAMMPS stores -// the precession force; the two differ by the factor hbar / |m|. +// Reduced Planck constant in eV.ps (metal units). The model reports the +// magnetic force as the energy gradient with respect to the magnetic moment. +// Under magnetic_force precession (the default), that gradient is stored in +// atom->fm after multiplication by |m| / hbar. Under magnetic_force cartesian +// the gradient is stored with only the energy-per-model-spin unit conversion, +// and this constant is unused. constexpr double kHBar = 6.5821191e-04; // Positions of the LAMMPS global virial components (xx, yy, zz, xy, xz, yz) // within the nine-component tensor the model reports. @@ -120,7 +123,11 @@ constexpr int kCentroidVirialMap[9] = {0, 4, 8, 3, 6, 7, 1, 2, 5}; } // namespace PairDPA4Spin::PairDPA4Spin(LAMMPS* lmp) - : Pair(lmp), scale(nullptr), cutoff(0.0), commdata_(nullptr) { + : Pair(lmp), + scale(nullptr), + cutoff(0.0), + fm_cartesian_force(0), + commdata_(nullptr) { if (lmp->citeme) { lmp->citeme->add(cite_user_deepmd_package); } @@ -237,46 +244,85 @@ void PairDPA4Spin::settings(int narg, char** arg) { } cutoff = deep_spin.cutoff() * dist_unit_cvt_factor; - // How many values name a charge state is a property of the artifact, so the - // keyword can only be read once the model is loaded. + // Optional keywords follow the artifact path. How many values name a charge + // state is a property of the artifact, so charge_spin can only be read once + // the model is loaded. magnetic_force selects the convention stored in + // atom->fm and is independent of whether the model uses native spin. const int dim_chg_spin = deep_spin.dim_chg_spin(); std::vector charge_spin; + bool magnetic_force_seen = false; + fm_cartesian_force = 0; int iarg = 1; while (iarg < narg) { - if (std::string(arg[iarg]) != "charge_spin") { + const std::string key(arg[iarg]); + if (key == "charge_spin") { + // One charge state holds for the whole run, so a second occurrence names + // no state the style could serve. Values already read mark the keyword + // as seen: a successful read always contributes at least one value. + if (!charge_spin.empty()) { + error->all(FLERR, + "Illegal pair_style command: keyword charge_spin names the " + "single charge state of the whole run, so it may appear " + "only once."); + } + if (dim_chg_spin == 0) { + error->all(FLERR, + "Illegal pair_style command: the artifact served by pair " + "style " + + style + + " carries no charge/spin condition, so keyword " + "charge_spin names nothing it can serve."); + } + if (iarg + dim_chg_spin >= narg) { + error->all(FLERR, + "Illegal pair_style command: keyword charge_spin names a " + "charge state with " + + std::to_string(dim_chg_spin) + " value(s)."); + } + for (int ii = 0; ii < dim_chg_spin; ++ii) { + charge_spin.push_back( + utils::numeric(FLERR, arg[iarg + 1 + ii], false, lmp)); + } + iarg += 1 + dim_chg_spin; + } else if (key == "magnetic_force") { + if (magnetic_force_seen) { + error->all(FLERR, + "Illegal pair_style command: keyword magnetic_force may " + "appear only once."); + } + if (iarg + 1 >= narg) { + error->all(FLERR, + "Illegal pair_style command: keyword magnetic_force takes " + "one of precession or cartesian."); + } + const std::string value(arg[iarg + 1]); + if (value == "precession") { + fm_cartesian_force = 0; + } else if (value == "cartesian") { + fm_cartesian_force = 1; + } else { + error->all(FLERR, + "Illegal pair_style command: keyword magnetic_force takes " + "one of precession or cartesian, not '" + + value + "'."); + } + magnetic_force_seen = true; + iarg += 2; + } else { error->all(FLERR, "Illegal pair_style command: pair style " + style + " takes the artifact path followed by the optional " - "keyword charge_spin, not '" + - std::string(arg[iarg]) + "'."); - } - // One charge state holds for the whole run, so a second occurrence names - // no state the style could serve. Values already read mark the keyword - // as seen: a successful read always contributes at least one value. - if (!charge_spin.empty()) { - error->all(FLERR, - "Illegal pair_style command: keyword charge_spin names the " - "single charge state of the whole run, so it may appear " - "only once."); - } - if (dim_chg_spin == 0) { - error->all(FLERR, - "Illegal pair_style command: the artifact served by pair " - "style " + - style + - " carries no charge/spin condition, so keyword " - "charge_spin names nothing it can serve."); - } - if (iarg + dim_chg_spin >= narg) { - error->all(FLERR, - "Illegal pair_style command: keyword charge_spin names a " - "charge state with " + - std::to_string(dim_chg_spin) + " value(s)."); + "keywords charge_spin and magnetic_force, not '" + + key + "'."); } - for (int ii = 0; ii < dim_chg_spin; ++ii) { - charge_spin.push_back( - utils::numeric(FLERR, arg[iarg + 1 + ii], false, lmp)); - } - iarg += 1 + dim_chg_spin; + } + + // The Cartesian convention stores an energy-per-spin gradient. Its only + // validated unit system is metal, where the model eV and the LAMMPS energy + // unit coincide and the precession path's metal hbar remains consistent. + if (fm_cartesian_force && strcmp(update->unit_style, "metal") != 0) { + error->all(FLERR, + "pair style " + style + + " keyword magnetic_force cartesian requires units metal."); } // A charge/spin condition named on the pair_style line holds for the whole @@ -305,6 +351,8 @@ void PairDPA4Spin::settings(int narg, char** arg) { } utils::logmesg(lmp, " using charge_spin: {}\n", values); } + utils::logmesg(lmp, " magnetic_force: {}\n", + fm_cartesian_force ? "cartesian" : "precession"); } /* ---------------------------------------------------------------------- @@ -536,16 +584,24 @@ void PairDPA4Spin::compute(int eflag, int vflag) { // === Step 3. Accumulate the model outputs === // The model reports a force and a magnetic force for every atom the neighbor // list covers, ghosts included; the spin atom style folds the ghost rows - // onto their owners in its reverse communication. Dividing the magnetic - // force by hbar / |m| turns the energy gradient with respect to the moment - // into the precession force LAMMPS stores. Multiplication by |m| preserves - // that relation while making a zero moment well-defined. + // onto their owners in its reverse communication. In the default precession + // mode, dividing the magnetic force by hbar / |m| turns the energy gradient + // with respect to the moment into the precession force LAMMPS stores; + // multiplication by |m| preserves that relation while making a zero moment + // well-defined. In Cartesian mode the energy gradient is stored directly, + // with only the energy-per-model-spin unit conversion, including at zero + // moment. const double force_scale = scale[1][1] * force_unit_cvt_factor; - const double magnetic_force_scale = force_scale / kHBar; + const double magnetic_force_scale = + fm_cartesian_force ? scale[1][1] * ener_unit_cvt_factor + : force_scale / kHBar; for (int ii = 0; ii < nall; ++ii) { + const double moment_factor = + fm_cartesian_force ? 1.0 : sp[ii][3]; for (int dd = 0; dd < 3; ++dd) { f[ii][dd] += dforce[3 * ii + dd] * force_scale; - fm[ii][dd] += dforce_mag[3 * ii + dd] * sp[ii][3] * magnetic_force_scale; + fm[ii][dd] += + dforce_mag[3 * ii + dd] * moment_factor * magnetic_force_scale; } } @@ -579,5 +635,9 @@ void* PairDPA4Spin::extract(const char* str, int& dim) { dim = 2; return (void*)scale; } + if (strcmp(str, "fm_cartesian_force") == 0) { + dim = 0; + return (void*)&fm_cartesian_force; + } return nullptr; } diff --git a/source/lmp/pair_dpa4spin.h b/source/lmp/pair_dpa4spin.h index 7d0aba8977..63b3fb6e64 100644 --- a/source/lmp/pair_dpa4spin.h +++ b/source/lmp/pair_dpa4spin.h @@ -66,8 +66,10 @@ class PairDPA4Spin : public Pair { PairDPA4Spin(class LAMMPS*); ~PairDPA4Spin() override; - // Load the artifact named by ``pair_style dpa4spin [charge_spin - // ... ]`` and fix the charge state the keyword names. + // Load the artifact named by ``pair_style dpa4spin + // [charge_spin ... ] [magnetic_force precession|cartesian]`` and + // fix the optional charge state and magnetic-force convention the keywords + // name. ``magnetic_force`` defaults to precession. void settings(int, char**) override; // Resolve the LAMMPS atom types onto the element list of the model. void coeff(int, char**) override; @@ -98,6 +100,14 @@ class PairDPA4Spin : public Pair { std::vector type_idx_map; // LAMMPS unit system relative to the eV / Angstrom system of the model. double ener_unit_cvt_factor, dist_unit_cvt_factor, force_unit_cvt_factor; + // Magnetic-force convention stored in ``atom->fm``. Zero selects the + // precession force (energy gradient scaled by the moment magnitude over + // hbar), which is the default and what the stock LAMMPS SPIN integrator + // expects. One selects the Cartesian energy gradient ``-dE/dS`` with only + // the energy-per-model-spin unit conversion, which a Newtonian spin + // integrator consumes. Exposed to callers through + // ``extract("fm_cartesian_force", dim)``. + int fm_cartesian_force; private: CommBrickDPA4Spin* commdata_; diff --git a/source/lmp/pair_dpa4spin_kokkos.cpp b/source/lmp/pair_dpa4spin_kokkos.cpp index f479b27559..3b8571d9c7 100644 --- a/source/lmp/pair_dpa4spin_kokkos.cpp +++ b/source/lmp/pair_dpa4spin_kokkos.cpp @@ -18,9 +18,12 @@ using namespace LAMMPS_NS; namespace { -// Reduced Planck constant in eV.ps. The model reports the magnetic force as -// the energy gradient with respect to the magnetic moment, while LAMMPS stores -// the precession force; the two differ by the moment magnitude over hbar. +// Reduced Planck constant in eV.ps (metal units). The model reports the +// magnetic force as the energy gradient with respect to the magnetic moment. +// Under magnetic_force precession (the default), that gradient is stored in +// atom->fm after multiplication by |m| / hbar. Under magnetic_force cartesian +// the gradient is stored with only the energy-per-model-spin unit conversion, +// and this constant is unused. constexpr double kHBar = 6.5821191e-04; } // namespace @@ -333,7 +336,11 @@ void PairDPA4SpinKokkos::compute(int eflag, int vflag) { // communication stages. auto model2loc = compact_graph.d_model2loc; const double fscale = scale[1][1] * force_unit_cvt_factor; - const double fmscale = scale[1][1] * force_unit_cvt_factor / kHBar; + // Precession: energy gradient times |m| / hbar. Cartesian: energy gradient + // with only the energy-per-model-spin unit conversion (no |m|, no hbar). + const int cartesian = fm_cartesian_force; + const double fmscale = cartesian ? scale[1][1] * ener_unit_cvt_factor + : scale[1][1] * force_unit_cvt_factor / kHBar; reverse_virial = false; reverse_used_host = false; // The KOKKOS package runs 'newton off', so the integrator's force_clear only @@ -364,15 +371,14 @@ void PairDPA4SpinKokkos::compute(int eflag, int vflag) { "dpa4spin/kk:scatter_f", Kokkos::RangePolicy(0, nnode_m), KOKKOS_LAMBDA(const int m) { const int i = model2loc(m); - // The precession force carried by ``fm`` is the energy gradient with - // respect to the moment, rescaled by the moment magnitude over hbar. - const double moment = fmscale * sp(i, 3); + // Precession multiplies by the moment magnitude; Cartesian does not. + const double moment_factor = cartesian ? 1.0 : sp(i, 3); f(i, 0) += fscale * out_force(3 * m + 0); f(i, 1) += fscale * out_force(3 * m + 1); f(i, 2) += fscale * out_force(3 * m + 2); - fm(i, 0) += moment * out_force_mag(3 * m + 0); - fm(i, 1) += moment * out_force_mag(3 * m + 1); - fm(i, 2) += moment * out_force_mag(3 * m + 2); + fm(i, 0) += fmscale * moment_factor * out_force_mag(3 * m + 0); + fm(i, 1) += fmscale * moment_factor * out_force_mag(3 * m + 1); + fm(i, 2) += fmscale * moment_factor * out_force_mag(3 * m + 2); }); atomKK->modified(execution_space, F_MASK | FM_MASK); if (multi_rank) { diff --git a/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py new file mode 100644 index 0000000000..af3b4524dc --- /dev/null +++ b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py @@ -0,0 +1,163 @@ +# SPDX-License-Identifier: LGPL-3.0-or-later +"""Unit checks for the ``dpa4spin`` / ``dpa4spin/kk`` magnetic-force conventions. + +The LAMMPS adapter stores the model energy gradient with respect to the +magnetic moment in ``atom->fm`` under one of two explicit conventions, +selected by ``pair_style ... magnetic_force precession|cartesian`` and +reported by ``Pair::extract("fm_cartesian_force", dim)``: + +* precession (default, extract value 0): ``fm = scale * force_cvt * |m| / hbar + * force_mag``, matching the stock LAMMPS SPIN integrator; +* cartesian (extract value 1): ``fm = scale * ener_cvt * force_mag``, the raw + ``-dE/dS`` gradient with only the energy-per-model-spin unit conversion, + including at zero moment. + +These tests encode that conversion contract without a LAMMPS build. Full +adapter coverage (``dpa4spin`` / ``dpa4spin/kk``, single-rank and MPI, +built-in and plugin, finite-difference comparison of Cartesian ``fm`` to +energy) belongs in CI with a Kokkos-enabled LAMMPS and a native-spin +archive; this module does not substitute for that. +""" + +from __future__ import annotations + +from pathlib import Path + +import numpy as np +import pytest + +# Metal-units reduced Planck constant used by pair_dpa4spin{,_kokkos}.cpp. +HBAR_METAL = 6.5821191e-04 + + +def store_magnetic_force( + force_mag: np.ndarray, + spin_norm: np.ndarray, + *, + cartesian: bool, + scale: float = 1.0, + ener_unit_cvt_factor: float = 1.0, + force_unit_cvt_factor: float = 1.0, +) -> np.ndarray: + """Apply the adapter's magnetic-force conversion (metal path by default). + + Parameters + ---------- + force_mag : np.ndarray + Model energy gradient with respect to the moment, shape (natoms, 3). + spin_norm : np.ndarray + Per-atom moment magnitude ``sp[:, 3]``, shape (natoms,). + cartesian : bool + True for ``magnetic_force cartesian``, False for precession. + scale : float + Pair ``scale[1][1]`` prefactor. + ener_unit_cvt_factor, force_unit_cvt_factor : float + LAMMPS unit conversions relative to the model eV / Angstrom system. + + Returns + ------- + np.ndarray + Values that the adapter accumulates into ``atom->fm``, shape (natoms, 3). + """ + force_mag = np.asarray(force_mag, dtype=np.float64) + spin_norm = np.asarray(spin_norm, dtype=np.float64) + if cartesian: + return scale * ener_unit_cvt_factor * force_mag + return (scale * force_unit_cvt_factor / HBAR_METAL) * spin_norm[:, None] * force_mag + + +@pytest.mark.parametrize( + "spin_norm", + [ + np.array([1.2737, 0.05, 0.0, 0.0]), + np.array([0.0, 0.0]), + np.array([2.5]), + ], +) +def test_cartesian_preserves_model_gradient_including_zero_moment( + spin_norm: np.ndarray, +) -> None: + """Cartesian mode must not multiply by |m| or divide by hbar.""" + rng = np.random.default_rng(0) + force_mag = rng.normal(size=(spin_norm.shape[0], 3)) + stored = store_magnetic_force(force_mag, spin_norm, cartesian=True) + np.testing.assert_allclose(stored, force_mag) + # Zero-moment atoms keep a non-zero Cartesian force when the model does. + zero = spin_norm == 0.0 + if np.any(zero): + assert np.any(np.abs(stored[zero]) > 0.0) + + +@pytest.mark.parametrize( + "spin_norm", + [ + np.array([1.2737, 0.05, 0.0, 0.0]), + np.array([0.0, 0.0]), + ], +) +def test_precession_scales_by_moment_over_hbar(spin_norm: np.ndarray) -> None: + """Precession mode multiplies by |m| / hbar; zero moment stays zero.""" + rng = np.random.default_rng(1) + force_mag = rng.normal(size=(spin_norm.shape[0], 3)) + stored = store_magnetic_force(force_mag, spin_norm, cartesian=False) + expected = force_mag * (spin_norm / HBAR_METAL)[:, None] + np.testing.assert_allclose(stored, expected) + zero = spin_norm == 0.0 + if np.any(zero): + np.testing.assert_array_equal(stored[zero], 0.0) + + +def test_precession_matches_cartesian_times_moment_over_hbar() -> None: + """On metal units the two conventions differ only by |m| / hbar.""" + spin_norm = np.array([1.1, 0.02, 0.0]) + force_mag = np.array( + [ + [0.3, -0.1, 0.2], + [1.0, 0.0, -0.5], + [4.0, 5.0, 6.0], + ], + dtype=np.float64, + ) + cartesian = store_magnetic_force(force_mag, spin_norm, cartesian=True) + precession = store_magnetic_force(force_mag, spin_norm, cartesian=False) + np.testing.assert_allclose( + precession, cartesian * (spin_norm / HBAR_METAL)[:, None] + ) + + +def test_adapter_declares_extract_key_and_keyword() -> None: + """The C++ adapter must keep the discoverable selection and extract API.""" + root = Path(__file__).resolve().parents[1] + pair_cpp = (root / "pair_dpa4spin.cpp").read_text(encoding="utf-8") + pair_h = (root / "pair_dpa4spin.h").read_text(encoding="utf-8") + assert 'strcmp(str, "fm_cartesian_force")' in pair_cpp + assert "fm_cartesian_force" in pair_h + assert 'key == "magnetic_force"' in pair_cpp + assert 'value == "cartesian"' in pair_cpp + assert 'value == "precession"' in pair_cpp + assert "dim = 0" in pair_cpp + assert 'strcmp(update->unit_style, "metal")' in pair_cpp + + +def test_non_unit_energy_conversion_applies_only_to_cartesian() -> None: + """Cartesian uses ener_unit_cvt_factor; precession uses force_unit_cvt_factor.""" + spin_norm = np.array([1.0]) + force_mag = np.array([[2.0, 0.0, 0.0]]) + cartesian = store_magnetic_force( + force_mag, + spin_norm, + cartesian=True, + ener_unit_cvt_factor=23.060549, + force_unit_cvt_factor=99.0, + ) + precession = store_magnetic_force( + force_mag, + spin_norm, + cartesian=False, + ener_unit_cvt_factor=23.060549, + force_unit_cvt_factor=2.0, + ) + np.testing.assert_allclose(cartesian, [[46.121098, 0.0, 0.0]]) + np.testing.assert_allclose( + precession, force_mag * (2.0 / HBAR_METAL) + ) From 9993cc2299f63727f582aaef3f0f40538085ee4f Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 15 Oct 2026 08:53:00 +0000 Subject: [PATCH 2/8] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- source/lmp/pair_dpa4spin.cpp | 9 ++++----- source/lmp/pair_dpa4spin_kokkos.cpp | 5 +++-- .../tests/test_dpa4spin_magnetic_force_convention.py | 12 +++++++----- 3 files changed, 14 insertions(+), 12 deletions(-) diff --git a/source/lmp/pair_dpa4spin.cpp b/source/lmp/pair_dpa4spin.cpp index 916f72323b..442a5f620b 100644 --- a/source/lmp/pair_dpa4spin.cpp +++ b/source/lmp/pair_dpa4spin.cpp @@ -592,12 +592,11 @@ void PairDPA4Spin::compute(int eflag, int vflag) { // with only the energy-per-model-spin unit conversion, including at zero // moment. const double force_scale = scale[1][1] * force_unit_cvt_factor; - const double magnetic_force_scale = - fm_cartesian_force ? scale[1][1] * ener_unit_cvt_factor - : force_scale / kHBar; + const double magnetic_force_scale = fm_cartesian_force + ? scale[1][1] * ener_unit_cvt_factor + : force_scale / kHBar; for (int ii = 0; ii < nall; ++ii) { - const double moment_factor = - fm_cartesian_force ? 1.0 : sp[ii][3]; + const double moment_factor = fm_cartesian_force ? 1.0 : sp[ii][3]; for (int dd = 0; dd < 3; ++dd) { f[ii][dd] += dforce[3 * ii + dd] * force_scale; fm[ii][dd] += diff --git a/source/lmp/pair_dpa4spin_kokkos.cpp b/source/lmp/pair_dpa4spin_kokkos.cpp index 3b8571d9c7..efedbc5af4 100644 --- a/source/lmp/pair_dpa4spin_kokkos.cpp +++ b/source/lmp/pair_dpa4spin_kokkos.cpp @@ -339,8 +339,9 @@ void PairDPA4SpinKokkos::compute(int eflag, int vflag) { // Precession: energy gradient times |m| / hbar. Cartesian: energy gradient // with only the energy-per-model-spin unit conversion (no |m|, no hbar). const int cartesian = fm_cartesian_force; - const double fmscale = cartesian ? scale[1][1] * ener_unit_cvt_factor - : scale[1][1] * force_unit_cvt_factor / kHBar; + const double fmscale = cartesian + ? scale[1][1] * ener_unit_cvt_factor + : scale[1][1] * force_unit_cvt_factor / kHBar; reverse_virial = false; reverse_used_host = false; // The KOKKOS package runs 'newton off', so the integrator's force_clear only diff --git a/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py index af3b4524dc..a474545aa5 100644 --- a/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py +++ b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py @@ -19,9 +19,13 @@ archive; this module does not substitute for that. """ -from __future__ import annotations +from __future__ import ( + annotations, +) -from pathlib import Path +from pathlib import ( + Path, +) import numpy as np import pytest @@ -158,6 +162,4 @@ def test_non_unit_energy_conversion_applies_only_to_cartesian() -> None: force_unit_cvt_factor=2.0, ) np.testing.assert_allclose(cartesian, [[46.121098, 0.0, 0.0]]) - np.testing.assert_allclose( - precession, force_mag * (2.0 / HBAR_METAL) - ) + np.testing.assert_allclose(precession, force_mag * (2.0 / HBAR_METAL)) From 7e994f857018b31230acb1ae60e4bfb412873929 Mon Sep 17 00:00:00 2001 From: OutisLi_Bot <198997215+OutisLi-Bot@users.noreply.github.com> Date: Thu, 15 Oct 2026 16:57:50 +0800 Subject: [PATCH 3/8] review(lammps): tighten Cartesian magnetic-force dual-path contract Align Kokkos precession scaling with the host fscale path, require the Kokkos scatter source in the convention unit test, and correct the example README wording that called Cartesian fm unconverted. --- examples/spin/dpa4c/lmp/README.md | 7 ++++--- source/lmp/pair_dpa4spin_kokkos.cpp | 5 ++--- .../lmp/tests/test_dpa4spin_magnetic_force_convention.py | 5 +++++ 3 files changed, 11 insertions(+), 6 deletions(-) diff --git a/examples/spin/dpa4c/lmp/README.md b/examples/spin/dpa4c/lmp/README.md index 42924a6c0b..8b3b9b85e6 100644 --- a/examples/spin/dpa4c/lmp/README.md +++ b/examples/spin/dpa4c/lmp/README.md @@ -70,6 +70,7 @@ force only from pair styles matching its own name pattern. The default magnetic-force convention is precession (`|m| / hbar` scaling of the model energy gradient). A separately maintained Newtonian Cartesian -integrator can request the unconverted gradient with -`magnetic_force cartesian` on the `pair_style` line (`units metal`) and can -confirm the choice through `Pair::extract("fm_cartesian_force", dim)`. +integrator can request `-dE/dS` with only the energy-per-model-spin unit +conversion via `magnetic_force cartesian` on the `pair_style` line (`units +metal`) and can confirm the choice through +`Pair::extract("fm_cartesian_force", dim)`. diff --git a/source/lmp/pair_dpa4spin_kokkos.cpp b/source/lmp/pair_dpa4spin_kokkos.cpp index efedbc5af4..50f1360eec 100644 --- a/source/lmp/pair_dpa4spin_kokkos.cpp +++ b/source/lmp/pair_dpa4spin_kokkos.cpp @@ -339,9 +339,8 @@ void PairDPA4SpinKokkos::compute(int eflag, int vflag) { // Precession: energy gradient times |m| / hbar. Cartesian: energy gradient // with only the energy-per-model-spin unit conversion (no |m|, no hbar). const int cartesian = fm_cartesian_force; - const double fmscale = cartesian - ? scale[1][1] * ener_unit_cvt_factor - : scale[1][1] * force_unit_cvt_factor / kHBar; + const double fmscale = + cartesian ? scale[1][1] * ener_unit_cvt_factor : fscale / kHBar; reverse_virial = false; reverse_used_host = false; // The KOKKOS package runs 'newton off', so the integrator's force_clear only diff --git a/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py index a474545aa5..6ad75d4da9 100644 --- a/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py +++ b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py @@ -134,6 +134,7 @@ def test_adapter_declares_extract_key_and_keyword() -> None: root = Path(__file__).resolve().parents[1] pair_cpp = (root / "pair_dpa4spin.cpp").read_text(encoding="utf-8") pair_h = (root / "pair_dpa4spin.h").read_text(encoding="utf-8") + pair_kk = (root / "pair_dpa4spin_kokkos.cpp").read_text(encoding="utf-8") assert 'strcmp(str, "fm_cartesian_force")' in pair_cpp assert "fm_cartesian_force" in pair_h assert 'key == "magnetic_force"' in pair_cpp @@ -141,6 +142,10 @@ def test_adapter_declares_extract_key_and_keyword() -> None: assert 'value == "precession"' in pair_cpp assert "dim = 0" in pair_cpp assert 'strcmp(update->unit_style, "metal")' in pair_cpp + # Primary issue target is dpa4spin/kk: scatter must keep both conventions. + assert "fm_cartesian_force" in pair_kk + assert "ener_unit_cvt_factor" in pair_kk + assert "cartesian ? 1.0 : sp(i, 3)" in pair_kk def test_non_unit_energy_conversion_applies_only_to_cartesian() -> None: From b70be2c577e28b3baacc808fc6fd69f8c010604f Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 15 Oct 2026 08:58:51 +0000 Subject: [PATCH 4/8] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- examples/spin/dpa4c/lmp/README.md | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/examples/spin/dpa4c/lmp/README.md b/examples/spin/dpa4c/lmp/README.md index 8b3b9b85e6..4215cf5464 100644 --- a/examples/spin/dpa4c/lmp/README.md +++ b/examples/spin/dpa4c/lmp/README.md @@ -71,6 +71,5 @@ force only from pair styles matching its own name pattern. The default magnetic-force convention is precession (`|m| / hbar` scaling of the model energy gradient). A separately maintained Newtonian Cartesian integrator can request `-dE/dS` with only the energy-per-model-spin unit -conversion via `magnetic_force cartesian` on the `pair_style` line (`units -metal`) and can confirm the choice through +conversion via `magnetic_force cartesian` on the `pair_style` line (`units metal`) and can confirm the choice through `Pair::extract("fm_cartesian_force", dim)`. From 2c4baaecc2695438e446738bddd0595c59591e34 Mon Sep 17 00:00:00 2001 From: OutisLi-Bot <44916828+OutisLi-Bot@users.noreply.github.com> Date: Thu, 15 Oct 2026 17:43:37 +0800 Subject: [PATCH 5/8] docs(lammps): use compressed archive in dpa4spin/kk cartesian example dpa4spin/kk requires the compact canonical graph lower; point the magnetic_force cartesian comment at compressed_model.pt2 and the README compress step instead of the graph-lower freeze output. --- examples/spin/dpa4c/lmp/README.md | 7 +++++-- examples/spin/dpa4c/lmp/in.lammps | 7 +++++-- 2 files changed, 10 insertions(+), 4 deletions(-) diff --git a/examples/spin/dpa4c/lmp/README.md b/examples/spin/dpa4c/lmp/README.md index 4215cf5464..d13be52950 100644 --- a/examples/spin/dpa4c/lmp/README.md +++ b/examples/spin/dpa4c/lmp/README.md @@ -71,5 +71,8 @@ force only from pair styles matching its own name pattern. The default magnetic-force convention is precession (`|m| / hbar` scaling of the model energy gradient). A separately maintained Newtonian Cartesian integrator can request `-dE/dS` with only the energy-per-model-spin unit -conversion via `magnetic_force cartesian` on the `pair_style` line (`units metal`) and can confirm the choice through -`Pair::extract("fm_cartesian_force", dim)`. +conversion via `magnetic_force cartesian` on the `pair_style` line +(`units metal`) and can confirm the choice through +`Pair::extract("fm_cartesian_force", dim)`. On `dpa4spin/kk`, pair that keyword +with the compressed compact-canonical archive (not the graph-lower freeze +output); host `dpa4spin` may use either archive. diff --git a/examples/spin/dpa4c/lmp/in.lammps b/examples/spin/dpa4c/lmp/in.lammps index 7c0d283ba5..d55fde1cdf 100644 --- a/examples/spin/dpa4c/lmp/in.lammps +++ b/examples/spin/dpa4c/lmp/in.lammps @@ -14,9 +14,12 @@ # # By default atom->fm holds the precession magnetic force. For a Newtonian # Cartesian integrator that needs -dE/dS, append `magnetic_force cartesian` -# (units metal) and read Pair::extract("fm_cartesian_force", dim): +# (units metal) and read Pair::extract("fm_cartesian_force", dim). Host +# `dpa4spin` accepts the graph-lower archive from freeze; `dpa4spin/kk` +# requires the compact canonical archive from `dp --pt-expt compress` +# (see README.md): # -# pair_style dpa4spin/kk frozen_model.pt2 magnetic_force cartesian +# pair_style dpa4spin/kk compressed_model.pt2 magnetic_force cartesian units metal dimension 3 From 810abceedfd82094710efbe00f8a18b8f92da0e7 Mon Sep 17 00:00:00 2001 From: OutisLi-Bot <44916828+OutisLi-Bot@users.noreply.github.com> Date: Thu, 15 Oct 2026 18:50:20 +0800 Subject: [PATCH 6/8] test(lammps): exercise dpa4spin magnetic_force via real adapter OutisLi CHANGES_REQUESTED on #6078: NumPy/substring tests could not see a += to -= mutation in the host/Kokkos scatter. Add behavioral coverage that selects magnetic_force through pair_style dpa4spin, checks extract fm_cartesian_force, compares Cartesian fm to energy FD (including a spin-active zero-moment Ni), and verifies default precession scaling. Register host dpa4spin in the LAMMPS plugin so plugin-mode CI can run the tests; generate a no-fparam DPA4C fixture for dpa4spin via gen_dpa4c_spin_lmp.py. Kokkos path skips without a GPU backend. --- source/install/test_cc_local.sh | 6 + source/lmp/plugin/deepmdplugin.cpp | 11 + ...test_dpa4spin_magnetic_force_convention.py | 335 +++++++++++++++++- source/tests/infer/gen_dpa4c_spin_lmp.py | 91 +++++ 4 files changed, 437 insertions(+), 6 deletions(-) create mode 100644 source/tests/infer/gen_dpa4c_spin_lmp.py diff --git a/source/install/test_cc_local.sh b/source/install/test_cc_local.sh index 17b4f2f593..8fc4fcc4d9 100755 --- a/source/install/test_cc_local.sh +++ b/source/install/test_cc_local.sh @@ -134,6 +134,11 @@ else: # exports NO with-comm artifact (single-rank only by construction). env ${_GEN_ENV} python ${INFER_SCRIPT_PATH}/gen_dpa4_spin_zbl.py & PID14=$! + # Minimal native-spin DPA4C for LAMMPS pair_style dpa4spin (no + # fparam/aparam/chg_spin). CPU-safe; used by the magnetic-force + # convention behavioral tests. + env ${_GEN_ENV} python ${INFER_SCRIPT_PATH}/gen_dpa4c_spin_lmp.py & + PID14B=$! # Native-spin DPA4C: the archives above all ship a with-comm artifact # and are therefore served by DeepSpinPTExpt, which leaves # NativeSpinPTExpt -- the class the dispatch selects for a native-spin @@ -149,6 +154,7 @@ else: wait $PID12 wait $PID13 wait $PID14 + wait $PID14B if [ "${DP_VARIANT}" = "cuda" ]; then wait $PID15 fi diff --git a/source/lmp/plugin/deepmdplugin.cpp b/source/lmp/plugin/deepmdplugin.cpp index 9bc0b685d2..46922f4ff1 100644 --- a/source/lmp/plugin/deepmdplugin.cpp +++ b/source/lmp/plugin/deepmdplugin.cpp @@ -13,6 +13,7 @@ #endif #include "pair_deepmd.h" #include "pair_deepspin.h" +#include "pair_dpa4spin.h" #include "version.h" #if LAMMPS_VERSION_NUMBER >= 20220328 #include "pppm_dplr.h" @@ -22,6 +23,7 @@ using namespace LAMMPS_NS; static Pair* pairdeepmd(LAMMPS* lmp) { return new PairDeepMD(lmp); } static Pair* pairdeepspin(LAMMPS* lmp) { return new PairDeepSpin(lmp); } +static Pair* pairdpa4spin(LAMMPS* lmp) { return new PairDPA4Spin(lmp); } #ifdef LMP_KOKKOS // Runtime plugins do not consume the PairStyle declarations used by LAMMPS' @@ -109,6 +111,15 @@ extern "C" void lammpsplugin_init(void* lmp, void* handle, void* regfunc) { plugin.handle = handle; (*register_plugin)(&plugin, lmp); + plugin.version = LAMMPS_VERSION; + plugin.style = "pair"; + plugin.name = "dpa4spin"; + plugin.info = "dpa4spin pair style " STR_GIT_SUMM; + plugin.author = "Tiancheng Li"; + plugin.creator.v1 = (lammpsplugin_factory1*)&pairdpa4spin; + plugin.handle = handle; + (*register_plugin)(&plugin, lmp); + plugin.style = "compute"; plugin.name = "deeptensor/atom"; plugin.info = "compute deeptensor/atom " STR_GIT_SUMM; diff --git a/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py index 6ad75d4da9..93a9d2488b 100644 --- a/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py +++ b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py @@ -1,5 +1,5 @@ # SPDX-License-Identifier: LGPL-3.0-or-later -"""Unit checks for the ``dpa4spin`` / ``dpa4spin/kk`` magnetic-force conventions. +"""Checks for the ``dpa4spin`` / ``dpa4spin/kk`` magnetic-force conventions. The LAMMPS adapter stores the model energy gradient with respect to the magnetic moment in ``atom->fm`` under one of two explicit conventions, @@ -12,17 +12,21 @@ ``-dE/dS`` gradient with only the energy-per-model-spin unit conversion, including at zero moment. -These tests encode that conversion contract without a LAMMPS build. Full -adapter coverage (``dpa4spin`` / ``dpa4spin/kk``, single-rank and MPI, -built-in and plugin, finite-difference comparison of Cartesian ``fm`` to -energy) belongs in CI with a Kokkos-enabled LAMMPS and a native-spin -archive; this module does not substitute for that. +NumPy helpers below document the conversion algebra. Behavioral tests drive +the registered ``dpa4spin`` pair style against +``deeppot_dpa4c_spin_lmp.pt2`` (no fparam/aparam/chg_spin) when PyLammps and +the style are available; they skip otherwise so CI without a dpa4spin build +stays green. ``dpa4spin/kk`` requires a GPU Kokkos backend and a compressed +canonical archive, so the Kokkos case is an explicit skip when unavailable. """ from __future__ import ( annotations, ) +import ctypes +import os +import tempfile from pathlib import ( Path, ) @@ -33,6 +37,43 @@ # Metal-units reduced Planck constant used by pair_dpa4spin{,_kokkos}.cpp. HBAR_METAL = 6.5821191e-04 +_MODEL_ENV = "DPA4SPIN_LMP_MODEL" +_DEFAULT_MODEL = ( + Path(__file__).resolve().parents[2] + / "tests" + / "infer" + / "deeppot_dpa4c_spin_lmp.pt2" +) + +# 6-atom NiO geometry shared with source/tests/infer/gen_dpa4c_spin.py +# (box edge 6 A = DPA4C ghost cutoff). Ni is spin-active; O is not. +_BOX = np.array([0.0, 6.0, 0.0, 6.0, 0.0, 6.0, 0.0, 0.0, 0.0]) +_COORD = np.array( + [ + [1.0, 1.0, 1.0], + [3.2, 1.4, 1.1], + [1.3, 1.8, 1.0], + [0.4, 1.2, 1.6], + [3.6, 2.0, 1.3], + [3.4, 0.7, 1.7], + ], + dtype=np.float64, +) +# Atom 0: spin-active Ni at zero moment (cartesian fm must stay nonzero). +# Atoms 1-2: nonzero Ni moments. Atoms 3-5: non-magnetic O. +_SPIN = np.array( + [ + [0.0, 0.0, 0.0], + [-0.07, 0.09, 0.03], + [0.02, -0.06, 0.08], + [0.0, 0.0, 0.0], + [0.0, 0.0, 0.0], + [0.0, 0.0, 0.0], + ], + dtype=np.float64, +) +_TYPE = np.array([1, 1, 1, 2, 2, 2], dtype=np.int32) + def store_magnetic_force( force_mag: np.ndarray, @@ -168,3 +209,285 @@ def test_non_unit_energy_conversion_applies_only_to_cartesian() -> None: ) np.testing.assert_allclose(cartesian, [[46.121098, 0.0, 0.0]]) np.testing.assert_allclose(precession, force_mag * (2.0 / HBAR_METAL)) + + +# --------------------------------------------------------------------------- +# Behavioral coverage through the real ``dpa4spin`` adapter +# --------------------------------------------------------------------------- + + +def _resolve_model() -> Path: + env = os.environ.get(_MODEL_ENV) + return Path(env) if env else _DEFAULT_MODEL + + +def _require_dpa4spin_host(): + """Import PyLammps and require a host ``dpa4spin`` registration.""" + try: + from lammps import ( + lammps, + ) + except ImportError: + pytest.skip("PyLammps is not available") + model = _resolve_model() + if not model.is_file(): + pytest.skip( + f"native-spin dpa4spin fixture missing ({model}); " + f"set {_MODEL_ENV} or run source/tests/infer/gen_dpa4c_spin_lmp.py" + ) + probe = lammps(cmdargs=["-log", "none", "-screen", "none"]) + try: + if not probe.has_style("pair", "dpa4spin"): + pytest.skip("pair style dpa4spin is not registered in this LAMMPS build") + finally: + probe.close() + return lammps, model + + +def _extract_fm_cartesian_force(lmp) -> int: + """Read the scalar int from ``Pair::extract("fm_cartesian_force", dim)``.""" + extract = lmp.lib.lammps_extract_pair + # The Python binding defaults restype to c_int, which truncates pointers. + previous = extract.restype + extract.restype = ctypes.c_void_p + try: + c_ptr = extract(lmp.lmp, b"fm_cartesian_force") + finally: + extract.restype = previous + if not c_ptr: + raise AssertionError("extract(fm_cartesian_force) returned null") + return int(ctypes.cast(c_ptr, ctypes.POINTER(ctypes.c_int)).contents.value) + + +def _open_dpa4spin(lammps_cls, data_file: Path, model: Path, pair_args: str): + """Build a single-rank metal/spin system with ``pair_style dpa4spin``.""" + lmp = lammps_cls(cmdargs=["-log", "none", "-screen", "none"]) + try: + for cmd in ( + "units metal", + "dimension 3", + "boundary p p p", + "atom_style spin", + "atom_modify map yes", + f"read_data {data_file.resolve()}", + "mass 1 58.6934", + "mass 2 15.9994", + f"pair_style dpa4spin {model.resolve()} {pair_args}".rstrip(), + "pair_coeff * * Ni O", + "neighbor 2.0 bin", + "neigh_modify every 1 delay 0 check yes", + "compute fmprop all property/atom fmx fmy fmz", + "compute spprop all property/atom spx spy spz sp", + "run 0", + ): + lmp.command(cmd) + except Exception: + lmp.close() + raise + return lmp + + +def _gather_fm_sp(lmp) -> tuple[np.ndarray, np.ndarray]: + natoms = int(lmp.get_natoms()) + fm = np.array(lmp.gather("c_fmprop", 1, 3), dtype=np.float64).reshape(natoms, 3) + sp = np.array(lmp.gather("c_spprop", 1, 4), dtype=np.float64).reshape(natoms, 4) + return fm, sp + + +@pytest.fixture(scope="module") +def dpa4spin_data(tmp_path_factory: pytest.TempPathFactory) -> Path: + """Write the mixed zero-/nonzero-moment NiO data file once per module.""" + from write_lmp_data import ( + write_lmp_data_spin, + ) + + path = tmp_path_factory.mktemp("dpa4spin_mf") / "nio_zero_and_finite.lmp" + write_lmp_data_spin(_BOX, _COORD, _SPIN, _TYPE, path) + return path + + +def test_dpa4spin_extract_reports_magnetic_force_selection( + dpa4spin_data: Path, +) -> None: + """``magnetic_force`` selection must be visible via Pair::extract.""" + lammps_cls, model = _require_dpa4spin_host() + for args, expect in ( + ("", 0), + ("magnetic_force precession", 0), + ("magnetic_force cartesian", 1), + ): + lmp = _open_dpa4spin(lammps_cls, dpa4spin_data, model, args) + try: + assert _extract_fm_cartesian_force(lmp) == expect + finally: + lmp.close() + + +def test_dpa4spin_cartesian_fm_matches_energy_fd_and_keeps_zero_moment( + dpa4spin_data: Path, +) -> None: + """Cartesian ``fm`` is ``-dE/dS`` and stays nonzero at zero moment on Ni.""" + lammps_cls, model = _require_dpa4spin_host() + lmp = _open_dpa4spin(lammps_cls, dpa4spin_data, model, "magnetic_force cartesian") + try: + fm, sp = _gather_fm_sp(lmp) + assert sp[0, 3] == pytest.approx(0.0) + assert np.linalg.norm(fm[0]) > 1e-8, ( + "spin-active Ni at zero moment must keep a nonzero Cartesian fm" + ) + assert sp[1, 3] > 0.0 + assert np.linalg.norm(fm[1]) > 1e-8 + + # Energy finite difference w.r.t. Sx of Ni atom 2 (1-based id 2). + # ``set atom ... spin`` takes magnitude then direction (mag > 0). + s0 = sp[1, :3] * sp[1, 3] + eps = 1e-5 + + def energy_at(moment: np.ndarray) -> float: + mag = float(np.linalg.norm(moment)) + assert mag > 0.0 + u = moment / mag + lmp.command( + f"set atom 2 spin {mag:.16e} {u[0]:.16e} {u[1]:.16e} {u[2]:.16e}" + ) + lmp.command("run 0") + return float(lmp.get_thermo("pe")) + + ep = energy_at(s0 + np.array([eps, 0.0, 0.0])) + em = energy_at(s0 - np.array([eps, 0.0, 0.0])) + energy_at(s0) + fm_restored, _ = _gather_fm_sp(lmp) + fd = -(ep - em) / (2.0 * eps) + analytic = fm_restored[1, 0] + rel = abs(fd - analytic) / max(1e-12, abs(analytic) + abs(fd)) + assert rel < 5e-3 or abs(analytic - fd) < 1e-8, ( + f"cartesian fm_x vs energy FD: analytic={analytic}, fd={fd}, rel={rel}" + ) + finally: + lmp.close() + + +def test_dpa4spin_default_precession_scales_cartesian_by_moment_over_hbar( + dpa4spin_data: Path, +) -> None: + """Default precession ``fm`` equals cartesian ``fm`` times ``|m| / hbar``.""" + lammps_cls, model = _require_dpa4spin_host() + lmp_c = _open_dpa4spin(lammps_cls, dpa4spin_data, model, "magnetic_force cartesian") + try: + fm_c, sp = _gather_fm_sp(lmp_c) + finally: + lmp_c.close() + + lmp_p = _open_dpa4spin(lammps_cls, dpa4spin_data, model, "") + try: + assert _extract_fm_cartesian_force(lmp_p) == 0 + fm_p, _ = _gather_fm_sp(lmp_p) + finally: + lmp_p.close() + + # Zero-moment spin-active Ni: precession must vanish; cartesian must not. + np.testing.assert_allclose(fm_p[0], 0.0, atol=1e-14) + assert np.linalg.norm(fm_c[0]) > 1e-8 + + # Nonzero Ni moments: ratio must be |m|/hbar. + for i in (1, 2): + m = sp[i, 3] + expect = m / HBAR_METAL + for d in range(3): + if abs(fm_c[i, d]) < 1e-18: + assert abs(fm_p[i, d]) < 1e-14 + continue + ratio = fm_p[i, d] / fm_c[i, d] + np.testing.assert_allclose(ratio, expect, rtol=1e-8, atol=0.0) + + # Non-magnetic O: both conventions stay at zero. + np.testing.assert_allclose(fm_c[3:], 0.0, atol=1e-14) + np.testing.assert_allclose(fm_p[3:], 0.0, atol=1e-14) + + +def test_dpa4spin_kokkos_cartesian_skipped_without_gpu_backend() -> None: + """Exercise ``dpa4spin/kk`` when a GPU Kokkos build is present; else skip.""" + try: + from lammps import ( + lammps, + ) + except ImportError: + pytest.skip("PyLammps is not available") + + model = _resolve_model() + if not model.is_file(): + pytest.skip(f"fixture missing: {model}") + + probe = lammps(cmdargs=["-log", "none", "-screen", "none"]) + try: + has_kk = probe.has_style("pair", "dpa4spin/kk") + has_kokkos = probe.has_package("KOKKOS") + finally: + probe.close() + if not has_kk or not has_kokkos: + pytest.skip("dpa4spin/kk or KOKKOS package not available") + + # dpa4spin/kk rejects host Kokkos backends and graph-lower (uncompressed) + # archives; a GPU + compressed canonical model is required. + try: + lmp = lammps( + cmdargs=[ + "-k", + "on", + "g", + "1", + "-sf", + "kk", + "-log", + "none", + "-screen", + "none", + ] + ) + except Exception as exc: + pytest.skip(f"GPU Kokkos backend unavailable: {exc}") + + try: + if not lmp.has_style("pair", "dpa4spin/kk"): + pytest.skip("dpa4spin/kk not registered under Kokkos") + # Probe whether this archive is accepted; graph-lower must skip. + from write_lmp_data import ( + write_lmp_data_spin, + ) + + with tempfile.TemporaryDirectory() as tmp: + data = Path(tmp) / "nio.lmp" + write_lmp_data_spin(_BOX, _COORD, _SPIN, _TYPE, data) + try: + for cmd in ( + "units metal", + "dimension 3", + "boundary p p p", + "atom_style spin", + "atom_modify map yes", + f"read_data {data.resolve()}", + "mass 1 58.6934", + "mass 2 15.9994", + f"pair_style dpa4spin/kk {model.resolve()} " + "magnetic_force cartesian", + "pair_coeff * * Ni O", + "neighbor 2.0 bin", + "neigh_modify every 1 delay 0 check yes", + "compute fmprop all property/atom fmx fmy fmz", + "run 0", + ): + lmp.command(cmd) + except Exception as exc: + msg = str(exc) + if "GPU backend only" in msg or "compact canonical" in msg: + pytest.skip(f"dpa4spin/kk cannot run here: {msg}") + raise + assert _extract_fm_cartesian_force(lmp) == 1 + natoms = int(lmp.get_natoms()) + fm = np.array(lmp.gather("c_fmprop", 1, 3), dtype=np.float64).reshape( + natoms, 3 + ) + assert np.linalg.norm(fm[0]) > 1e-8 + assert np.linalg.norm(fm[1]) > 1e-8 + finally: + lmp.close() diff --git a/source/tests/infer/gen_dpa4c_spin_lmp.py b/source/tests/infer/gen_dpa4c_spin_lmp.py new file mode 100644 index 0000000000..2e9798f2fa --- /dev/null +++ b/source/tests/infer/gen_dpa4c_spin_lmp.py @@ -0,0 +1,91 @@ +#!/usr/bin/env python3 +# SPDX-License-Identifier: LGPL-3.0-or-later +"""Generate a minimal native-spin DPA4C graph ``.pt2`` for LAMMPS ``dpa4spin``. + +Unlike ``gen_dpa4c_spin.py``, this archive carries no fparam, aparam, or +charge/spin condition so ``pair_style dpa4spin`` can load it without extra +keywords. Used by ``source/lmp/tests/test_dpa4spin_magnetic_force_convention.py``. +""" + +from __future__ import annotations + +import copy +import json +import os +import sys +import zipfile + +import numpy as np + +sys.path.insert(0, os.path.join(os.path.dirname(__file__), "..")) + +from dpa4_fixtures import ( + jitter_zero_arrays, +) +from gen_common import ( + ensure_inductor_compiler, + load_custom_ops, +) + +CFG = { + "type_map": ["Ni", "O"], + "descriptor": { + "type": "dpa4c", + "rcut": 4.0, + "channels": 8, + "lmax": 2, + "n_radial": 4, + "precision": "float64", + "seed": 7, + }, + "fitting_net": {"neuron": [8, 8], "precision": "float64", "seed": 7}, + "spin": {"use_spin": [True, False], "scheme": "native"}, +} + +_JITTER_SEED = 20260817 +OUT = os.path.join(os.path.dirname(__file__), "deeppot_dpa4c_spin_lmp.pt2") + + +def main() -> None: + from deepmd.dpmodel.model.model import ( + get_model, + ) + from deepmd.pt_expt.utils.serialization import ( + deserialize_to_file, + ) + + ensure_inductor_compiler() + load_custom_ops() + model = get_model(copy.deepcopy(CFG)) + md = jitter_zero_arrays(model.serialize(), np.random.default_rng(_JITTER_SEED)) + data = { + "model": md, + "model_def_script": CFG, + "backend": "dpmodel", + "software": "deepmd-kit", + "version": "3.0.0", + } + print(f"Exporting {OUT}") # noqa: T201 + deserialize_to_file(OUT, data, do_atomic_virial=True, lower_kind="graph") + with zipfile.ZipFile(OUT) as zf: + meta = json.loads(zf.read("model/extra/metadata.json").decode()) + print( # noqa: T201 + json.dumps( + { + k: meta.get(k) + for k in ( + "spin_scheme", + "has_comm_artifact", + "dim_fparam", + "dim_aparam", + "dim_chg_spin", + "lower_input_kind", + ) + }, + indent=2, + ) + ) + + +if __name__ == "__main__": + main() From f6d1a1bcb72939ffc849a63eb9e7043fbb7823d2 Mon Sep 17 00:00:00 2001 From: OutisLi-Bot <44916828+OutisLi-Bot@users.noreply.github.com> Date: Thu, 15 Oct 2026 18:54:06 +0800 Subject: [PATCH 7/8] test(lammps): FD largest |fm| component without abs escape Probe the live cartesian fm component and require |analytic| > 1e-8 with a relative-only check so a scatter sign flip cannot pass. --- ...test_dpa4spin_magnetic_force_convention.py | 25 ++++++++++++++----- 1 file changed, 19 insertions(+), 6 deletions(-) diff --git a/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py index 93a9d2488b..af1bbe38b1 100644 --- a/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py +++ b/source/lmp/tests/test_dpa4spin_magnetic_force_convention.py @@ -338,10 +338,17 @@ def test_dpa4spin_cartesian_fm_matches_energy_fd_and_keeps_zero_moment( assert sp[1, 3] > 0.0 assert np.linalg.norm(fm[1]) > 1e-8 - # Energy finite difference w.r.t. Sx of Ni atom 2 (1-based id 2). + # Energy FD on the largest |fm| component of Ni atom 2 (1-based id 2). # ``set atom ... spin`` takes magnitude then direction (mag > 0). s0 = sp[1, :3] * sp[1, 3] eps = 1e-5 + comp = int(np.argmax(np.abs(fm[1]))) + analytic0 = float(fm[1, comp]) + assert abs(analytic0) > 1e-8, ( + f"largest |fm[1]| component must be live for FD: fm[1]={fm[1]}" + ) + ehat = np.zeros(3, dtype=np.float64) + ehat[comp] = 1.0 def energy_at(moment: np.ndarray) -> float: mag = float(np.linalg.norm(moment)) @@ -353,15 +360,21 @@ def energy_at(moment: np.ndarray) -> float: lmp.command("run 0") return float(lmp.get_thermo("pe")) - ep = energy_at(s0 + np.array([eps, 0.0, 0.0])) - em = energy_at(s0 - np.array([eps, 0.0, 0.0])) + ep = energy_at(s0 + eps * ehat) + em = energy_at(s0 - eps * ehat) energy_at(s0) fm_restored, _ = _gather_fm_sp(lmp) fd = -(ep - em) / (2.0 * eps) - analytic = fm_restored[1, 0] + analytic = float(fm_restored[1, comp]) + assert abs(analytic) > 1e-8, ( + f"analytic fm[{comp}] vanished after restore: analytic={analytic}" + ) + # Relative only: an absolute escape of ~1e-8 would let a scatter + # sign flip (+= -> -=) pass when |fm| is tiny; with |analytic| > 1e-8, + # |analytic - (-analytic)| = 2|analytic| fails both rel and a tight abs. rel = abs(fd - analytic) / max(1e-12, abs(analytic) + abs(fd)) - assert rel < 5e-3 or abs(analytic - fd) < 1e-8, ( - f"cartesian fm_x vs energy FD: analytic={analytic}, fd={fd}, rel={rel}" + assert rel < 5e-3, ( + f"cartesian fm[{comp}] vs energy FD: analytic={analytic}, fd={fd}, rel={rel}" ) finally: lmp.close() From 22c8a06034b5b6d0f6fe96844a56c5e11b17dfe0 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 15 Oct 2026 10:55:58 +0000 Subject: [PATCH 8/8] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- source/tests/infer/gen_dpa4c_spin_lmp.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/source/tests/infer/gen_dpa4c_spin_lmp.py b/source/tests/infer/gen_dpa4c_spin_lmp.py index 2e9798f2fa..e0b2fecf0d 100644 --- a/source/tests/infer/gen_dpa4c_spin_lmp.py +++ b/source/tests/infer/gen_dpa4c_spin_lmp.py @@ -7,7 +7,9 @@ keywords. Used by ``source/lmp/tests/test_dpa4spin_magnetic_force_convention.py``. """ -from __future__ import annotations +from __future__ import ( + annotations, +) import copy import json