diff --git a/dpgen2/op/run_caly_model_devi.py b/dpgen2/op/run_caly_model_devi.py index 9e191cf2..9258afa5 100644 --- a/dpgen2/op/run_caly_model_devi.py +++ b/dpgen2/op/run_caly_model_devi.py @@ -1,3 +1,4 @@ +import logging from collections import ( defaultdict, ) @@ -15,10 +16,14 @@ OPIO, Artifact, BigParameter, + FatalError, OPIOSign, Parameter, ) +from dpgen2.exploration.selector.distance_conf_filter import ( + safe_dist_dict as full_safe_dist_dict, +) from dpgen2.utils import ( set_directory, ) @@ -266,6 +271,15 @@ def parse_traj(traj_file): "H": 0.813, } + def safe_radius(symbol): + """Return a CALYPSO safe radius, preferring legacy tuned values.""" + if symbol in safe_dist_dict: + return safe_dist_dict[symbol] + try: + return full_safe_dist_dict[symbol] + except KeyError: + raise FatalError(f"no safe distance known for element {symbol!r}") from None + trajs: List[Atoms] = read(traj_file, index=":", format="traj") # type: ignore dthresh = 0.72 numb_traj = len(trajs) @@ -307,11 +321,7 @@ def parse_traj(traj_file): for a in range(len(atype)): for b in range(a + 1, len(atype)): dd = dist_dict[a][b] - dr = ( - (safe_dist_dict[atype[a]] + safe_dist_dict[atype[b]]) - * 0.529 - / 1.2 - ) + dr = (safe_radius(atype[a]) + safe_radius(atype[b])) * 0.529 / 1.2 if dd < dr: frame_is_reasonable = False @@ -321,6 +331,13 @@ def parse_traj(traj_file): else: selected_traj = None + if selected_traj == []: + logging.warning( + "All frames in CALYPSO trajectory %s were rejected by " + "safe-distance filtering.", + traj_file, + ) + return selected_traj diff --git a/tests/op/test_run_caly_model_devi.py b/tests/op/test_run_caly_model_devi.py index 647604a8..bcc686ad 100644 --- a/tests/op/test_run_caly_model_devi.py +++ b/tests/op/test_run_caly_model_devi.py @@ -271,6 +271,51 @@ def test_00_parse_traj(self): atoms_list_10 = parse_traj(self.traj_file_10) self.assertTrue(atoms_list_10 is None), self.atoms_abnormal + def test_parse_traj_supports_elements_outside_legacy_table(self): + """Use the shared safe-distance table for elements such as oxygen.""" + oxygen_traj = self.work_dir / "oxygen.traj" + oxygen = Atoms( + symbols=["O", "O"], + scaled_positions=[[0, 0, 0], [0.5, 0.5, 0.5]], + cell=np.eye(3) * 10, + ) + write(oxygen_traj, oxygen, format="traj") + + selected = parse_traj(oxygen_traj) + + self.assertIsNotNone(selected) + self.assertEqual(len(selected), 1) + + def test_parse_traj_applies_shared_nitrogen_safe_distance(self): + """Keep and reject nitrogen frames on opposite sides of the limit.""" + # The shared table gives an N-N threshold of 1.1638 Angstrom. These + # distances also distinguish it from the previous ASE fallback + # threshold of 1.1833 Angstrom. + accepted_traj = self.work_dir / "nitrogen-accepted.traj" + accepted = Atoms( + symbols=["N", "N"], + positions=[[0, 0, 0], [0, 0, 1.17]], + cell=np.eye(3) * 12, + pbc=True, + ) + write(accepted_traj, accepted, format="traj") + + selected = parse_traj(accepted_traj) + + self.assertIsNotNone(selected) + self.assertEqual(len(selected), 1) + + rejected_traj = self.work_dir / "nitrogen-rejected.traj" + rejected = accepted.copy() + rejected.positions[1, 2] = 1.15 + write(rejected_traj, rejected, format="traj") + + with self.assertLogs(level="WARNING") as captured: + selected = parse_traj(rejected_traj) + + self.assertEqual(selected, []) + self.assertIn("All frames in CALYPSO trajectory", captured.output[0]) + def test_01_atoms2lmpdump(self): dump_str = atoms2lmpdump(self.atoms_normal_2, 1, self.type_map) self.assertEqual(dump_str, self.ref_dump_str)