diff --git a/source/lmp/fix_dplr.cpp b/source/lmp/fix_dplr.cpp index e66cb643f9..cdd91fb1c6 100644 --- a/source/lmp/fix_dplr.cpp +++ b/source/lmp/fix_dplr.cpp @@ -125,6 +125,8 @@ FixDPLR::FixDPLR(LAMMPS *lmp, int narg, char **arg) comm_reverse = 3; } +/* ---------------------------------------------------------------------- */ + int FixDPLR::setmask() { int mask = 0; #if LAMMPS_VERSION_NUMBER < 20210210 @@ -134,9 +136,14 @@ int FixDPLR::setmask() { mask |= POST_INTEGRATE; mask |= PRE_FORCE; mask |= POST_FORCE; + mask |= MIN_PRE_EXCHANGE; + mask |= MIN_PRE_FORCE; + mask |= MIN_POST_FORCE; return mask; } +/* ---------------------------------------------------------------------- */ + void FixDPLR::init() { // double **xx = atom->x; // double **vv = atom->v; @@ -152,8 +159,12 @@ void FixDPLR::init() { // } } +/* ---------------------------------------------------------------------- */ + void FixDPLR::setup_pre_force(int vflag) { pre_force(vflag); } +/* ---------------------------------------------------------------------- */ + void FixDPLR::setup(int vflag) { // if (strstr(update->integrate_style,"verlet")) post_force(vflag); @@ -167,6 +178,12 @@ void FixDPLR::setup(int vflag) { } } +/* ---------------------------------------------------------------------- */ + +void FixDPLR::min_setup(int vflag) { setup(vflag); } + +/* ---------------------------------------------------------------------- */ + void FixDPLR::get_valid_pairs(vector > &pairs) { pairs.clear(); @@ -215,6 +232,8 @@ void FixDPLR::get_valid_pairs(vector > &pairs) { } } +/* ---------------------------------------------------------------------- */ + void FixDPLR::post_integrate() { double **x = atom->x; double **v = atom->v; @@ -237,6 +256,8 @@ void FixDPLR::post_integrate() { } } +/* ---------------------------------------------------------------------- */ + void FixDPLR::pre_force(int vflag) { double **x = atom->x; int *type = atom->type; @@ -356,6 +377,8 @@ void FixDPLR::pre_force(int vflag) { // } } +/* ---------------------------------------------------------------------- */ + void FixDPLR::post_force(int vflag) { if (vflag) { v_setup(vflag); @@ -515,6 +538,20 @@ void FixDPLR::post_force(int vflag) { } } +/* ---------------------------------------------------------------------- */ + +void FixDPLR::min_pre_exchange() { post_integrate(); } + +/* ---------------------------------------------------------------------- */ + +void FixDPLR::min_pre_force(int vflag) { pre_force(vflag); } + +/* ---------------------------------------------------------------------- */ + +void FixDPLR::min_post_force(int vflag) { post_force(vflag); } + +/* ---------------------------------------------------------------------- */ + int FixDPLR::pack_reverse_comm(int n, int first, double *buf) { int m = 0; int last = first + n; diff --git a/source/lmp/fix_dplr.h b/source/lmp/fix_dplr.h index 953e9ad183..f7b0361eec 100644 --- a/source/lmp/fix_dplr.h +++ b/source/lmp/fix_dplr.h @@ -40,11 +40,15 @@ class FixDPLR : public Fix { ~FixDPLR() override{}; int setmask() override; void init() override; - void setup_pre_force(int) override; void setup(int) override; + void setup_pre_force(int) override; + void min_setup(int) override; void post_integrate() override; void pre_force(int) override; void post_force(int) override; + void min_pre_exchange() override; + void min_pre_force(int) override; + void min_post_force(int) override; int pack_reverse_comm(int, int, double *) override; void unpack_reverse_comm(int, int *, double *) override; double compute_scalar(void) override; diff --git a/source/lmp/tests/test_dplr.py b/source/lmp/tests/test_dplr.py index b92e0ddbfc..0e74e8fc45 100644 --- a/source/lmp/tests/test_dplr.py +++ b/source/lmp/tests/test_dplr.py @@ -126,6 +126,36 @@ ] ) +expected_x_min_step1 = np.array( + [ + [1.26321372, 1.28623039, 0.97992808], + [0.94931355, 3.25164736, 1.32847644], + [0.64743204, 1.32472127, 1.76763885], + [1.30026330, 0.31505758, 0.71705476], + [1.91176075, 3.50957943, 1.43047464], + [0.50553664, 4.04638297, 0.90796723], + [1.22952996, 1.15353755, 1.01171484], + [0.97744342, 3.38790804, 1.31885315], + ] +) + +expected_e_min_step1 = -40.46708779 + +expected_e_kspace_min_step1 = 0.16209504 + +expected_f_min_step1 = np.array( + [ + [-0.04230911, -0.09057170, -0.00644688], + [0.03377777, 0.01767955, -0.18004346], + [-0.10280108, 0.22936656, 0.17841170], + [0.13239183, 0.02472290, -0.12924794], + [0.11956397, -0.13352759, 0.08242015], + [-0.14062339, -0.04766971, 0.05490644], + [0.00000000, 0.00000000, 0.00000000], + [0.00000000, 0.00000000, 0.00000000], + ] +) + box = np.array([0, 20, 0, 20, 0, 20, 0, 0, 0]) coord = np.array( [ @@ -230,10 +260,30 @@ def test_pair_deepmd_lr(lammps): lammps.fix(f"0 all dplr model {pb_file.resolve()} type_associate 1 3 bond_type 1") lammps.fix_modify("0 virial yes") lammps.run(0) - assert lammps.eval("pe") == pytest.approx(expected_e_lr) + for ii in range(2): + assert lammps.atoms[6 + ii].position == pytest.approx(expected_WC[ii]) assert lammps.eval("elong") == pytest.approx(expected_e_kspace) + assert lammps.eval("pe") == pytest.approx(expected_e_lr) for ii in range(6): assert lammps.atoms[ii].force == pytest.approx(expected_f_lr[ii]) - for ii in range(2): - assert lammps.atoms[6 + ii].position == pytest.approx(expected_WC[ii]) lammps.run(1) + + +def test_min_dplr(lammps): + lammps.pair_style(f"deepmd {pb_file.resolve()}") + lammps.pair_coeff("* *") + lammps.bond_style("zero") + lammps.bond_coeff("*") + lammps.special_bonds("lj/coul 1 1 1 angle no") + lammps.kspace_style("pppm/dplr 1e-5") + lammps.kspace_modify(f"gewald {beta:.2f} diff ik mesh {mesh:d} {mesh:d} {mesh:d}") + lammps.fix(f"0 all dplr model {pb_file.resolve()} type_associate 1 3 bond_type 1") + lammps.fix_modify("0 virial yes") + lammps.min_style("cg") + lammps.minimize("0 1.0e-6 2 2") + for ii in range(8): + assert lammps.atoms[ii].position == pytest.approx(expected_x_min_step1[ii]) + assert lammps.eval("pe") == pytest.approx(expected_e_min_step1) + assert lammps.eval("elong") == pytest.approx(expected_e_kspace_min_step1) + for ii in range(8): + assert lammps.atoms[ii].force == pytest.approx(expected_f_min_step1[ii])