fix: skip zero-force fabrication when ABACUS MD dump lacks FORCE columns

dpdata/formats/abacus/md.py

@@ -72,7 +72,7 @@ def get_coords_from_dump(dumplines, natoms):
     )
     cells = np.zeros([nframes_dump, 3, 3])
     stresses = np.zeros([nframes_dump, 3, 3])
-    forces = np.zeros([nframes_dump, total_natoms, 3])
+    forces = np.zeros([nframes_dump, total_natoms, 3]) if calc_force else None
     coords = np.zeros([nframes_dump, total_natoms, 3])
     iframe = 0
     for iline in range(nlines):
@@ -115,7 +115,7 @@ def get_coords_from_dump(dumplines, natoms):
                 if not newversion:
                     coords[iframe, iat] *= celldm
 
-                if calc_force:
+                if calc_force and forces is not None:
                     forces[iframe, iat] = np.array(
                         [
                             float(i)
@@ -187,7 +187,8 @@ def get_frame(fname):
         if np.isnan(iene):
             coords = np.delete(coords, i - ndump, axis=0)
             cells = np.delete(cells, i - ndump, axis=0)
-            force = np.delete(force, i - ndump, axis=0)
+            if force is not None:
+                force = np.delete(force, i - ndump, axis=0)
             stress = np.delete(stress, i - ndump, axis=0)
             energy = np.delete(energy, i - ndump, axis=0)
             unconv_stru += "%d " % i  # noqa: UP031
@@ -207,7 +208,8 @@ def get_frame(fname):
     #    data['cells'][:, :, :] = cell
     data["coords"] = coords
     data["energies"] = energy
-    data["forces"] = force
+    if force is not None:
+        data["forces"] = force
     data["virials"] = stress
     if not isinstance(data["virials"], np.ndarray):
         del data["virials"]

tests/abacus.md.noforce/INPUT

@@ -0,0 +1,2 @@
+suffix autotest
+md_dumpfreq 1

tests/abacus.md.noforce/OUT.autotest/MD_dump

@@ -0,0 +1,22 @@
+MDSTEP:  0
+LATTICE_CONSTANT: 10.200000000000
+LATTICE_VECTORS
+  0.500000000000  0.500000000000  0.000000000000
+  0.500000000000  0.000000000000  0.500000000000
+  0.000000000000  0.500000000000  0.500000000000
+INDEX    LABEL    POSITIONS
+  0  Si  0.000000000000  0.000000000000  0.000000000000
+  1  Si  0.241000000000  0.255000000000  0.251000000000
+
+
+MDSTEP:  1
+LATTICE_CONSTANT: 10.200000000000
+LATTICE_VECTORS
+  0.500000000000  0.500000000000  0.000000000000
+  0.500000000000  0.000000000000  0.500000000000
+  0.000000000000  0.500000000000  0.500000000000
+INDEX    LABEL    POSITIONS
+  0  Si  0.010000000000  0.010000000000  0.010000000000
+  1  Si  0.251000000000  0.265000000000  0.261000000000
+
+

tests/abacus.md.noforce/OUT.autotest/running_md.log

@@ -0,0 +1,8 @@
+ STEP OF MOLECULAR DYNAMICS : 0
+ -------------------------------------------
+ final etot is -211.771846025 eV
+
+ -------------------------------------------
+ STEP OF MOLECULAR DYNAMICS : 1
+ -------------------------------------------
+ final etot is -211.781119663 eV

tests/abacus.md.noforce/STRU

@@ -0,0 +1,19 @@
+ATOMIC_SPECIES
+Si 1 ../tools/PP_ORB/Si_ONCV_PBE-1.0.upf
+
+LATTICE_CONSTANT
+10.2
+
+LATTICE_VECTORS
+0.5 0.5 0 #latvec1
+0.5 0 0.5 #latvec2
+0 0.5 0.5 #latvec3
+
+ATOMIC_POSITIONS
+Cartesian
+
+Si #label
+0 #magnetism
+2 #number of atoms
+0  0  0  m  1  1  1  v  -0  0  0
+0.241  0.255  0.251  m  1  1  1  v  0  -0  -0

tests/test_abacus_md.py

@@ -249,6 +249,17 @@ def test_to_system(self):
                     iline += 1
                 self.assertEqual(iline, 30)
 
+    def test_no_force_columns(self):
+        """Test that MD dumps without FORCE columns do not fabricate zero forces."""
+        system_noforce = dpdata.LabeledSystem("abacus.md.noforce", fmt="abacus/md")
+        # When FORCE is absent from the dump, forces should not be in data
+        self.assertNotIn("forces", system_noforce.data)
+        # Other data should still be present
+        self.assertIn("coords", system_noforce.data)
+        self.assertIn("energies", system_noforce.data)
+        self.assertIn("cells", system_noforce.data)
+        self.assertEqual(len(system_noforce.data["coords"]), 2)
+
 
 if __name__ == "__main__":
     unittest.main()