Skip to content

Commit

Permalink
Various TS check improvements (#715)
Browse files Browse the repository at this point in the history
This PR contains various commit from the infamous `ts2` branch.
Also, added an option for users to specify `skip_nmd: True` in the input
file to skip the TS normal mode displacement check if it is problematic.
  • Loading branch information
calvinp0 authored Nov 26, 2023
2 parents 976fde1 + e7aea63 commit 36c361d
Show file tree
Hide file tree
Showing 16 changed files with 107 additions and 61 deletions.
38 changes: 25 additions & 13 deletions arc/checks/ts.py
Original file line number Diff line number Diff line change
Expand Up @@ -46,6 +46,7 @@ def check_ts(reaction: 'ARCReaction',
output: Optional[dict] = None,
sp_level: Optional['Level'] = None,
freq_scale_factor: float = 1.0,
skip_nmd: bool = False,
verbose: bool = True,
):
"""
Expand All @@ -68,6 +69,7 @@ def check_ts(reaction: 'ARCReaction',
output (dict, optional): The Scheduler output dictionary.
sp_level (Level, optional): The single-point energy level of theory.
freq_scale_factor (float, optional): The frequency scaling factor.
skip_nmd (bool, optional): Whether to skip the normal mode displacement check.
verbose (bool, optional): Whether to print logging messages.
"""
checks = checks or list()
Expand All @@ -89,19 +91,26 @@ def check_ts(reaction: 'ARCReaction',
if reaction.ts_species.ts_checks['E0'] is None and not reaction.ts_species.ts_checks['e_elect']:
check_rxn_e_elect(reaction=reaction, verbose=verbose)

if 'freq' in checks or (not reaction.ts_species.ts_checks['normal_mode_displacement'] and job is not None):
check_normal_mode_displacement(reaction, job=job)

if 'rotors' in checks or (ts_passed_all_checks(species=reaction.ts_species, exemptions=['E0', 'warnings', 'IRC'])
if 'freq' in checks or (not reaction.ts_species.ts_checks['NMD'] and job is not None):
try:
check_normal_mode_displacement(reaction, job=job)
except (ValueError, KeyError) as e:
logger.error(f'Could not check normal mode displacement, got: \n{e}')
reaction.ts_species.ts_checks['NMD'] = True
if skip_nmd and not reaction.ts_species.ts_checks['NMD']:
logger.warning(f'Skipping normal mode displacement check for TS {reaction.ts_species.label}')
reaction.ts_species.ts_checks['NMD'] = True

if 'rotors' in checks or (ts_passed_checks(species=reaction.ts_species, exemptions=['E0', 'warnings', 'IRC'])
and job is not None):
invalidate_rotors_with_both_pivots_in_a_reactive_zone(reaction, job,
rxn_zone_atom_indices=rxn_zone_atom_indices)


def ts_passed_all_checks(species: 'ARCSpecies',
exemptions: Optional[List[str]] = None,
verbose: bool = False,
) -> bool:
def ts_passed_checks(species: 'ARCSpecies',
exemptions: Optional[List[str]] = None,
verbose: bool = False,
) -> bool:
"""
Check whether the TS species passes all checks other than ones specified in ``exemptions``.
Expand Down Expand Up @@ -239,6 +248,9 @@ def check_rxn_e0(reaction: 'ARCReaction',
chosen_ts=reaction.ts_species.chosen_ts, energy_type='E0')
if r_e0 >= ts_e0 or p_e0 >= ts_e0:
reaction.ts_species.ts_checks['E0'] = False
if r_e0 + 1 >= ts_e0 or p_e0 + 1 >= ts_e0:
logger.warning('TS energy gas relative to one fo the wells is lower than 1 kJ/mol, skipping this TS')
reaction.ts_species.ts_checks['E0'] = False
else:
reaction.ts_species.ts_checks['E0'] = True

Expand Down Expand Up @@ -294,7 +306,7 @@ def check_normal_mode_displacement(reaction: 'ARCReaction',
rmgdb.determine_family(reaction)
amplitudes = amplitudes or [0.1, 0.2, 0.4, 0.6, 0.8, 1]
amplitudes = [amplitudes] if isinstance(amplitudes, float) else amplitudes
reaction.ts_species.ts_checks['normal_mode_displacement'] = False
reaction.ts_species.ts_checks['NMD'] = False
rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=reaction) or list()
freqs, normal_modes_disp = parser.parse_normal_mode_displacement(path=job.local_path_to_output_file, raise_error=False)
if not len(normal_modes_disp):
Expand Down Expand Up @@ -335,13 +347,13 @@ def check_normal_mode_displacement(reaction: 'ARCReaction',
and not any(entry is None for entry in breaking) and not any(entry is None for entry in forming) \
and all(entry == breaking[0] for entry in breaking) and all(entry == forming[0] for entry in forming) \
and breaking[0] != forming[0]:
reaction.ts_species.ts_checks['normal_mode_displacement'] = True
reaction.ts_species.ts_checks['NMD'] = True
done = True
break
if not got_expected_changing_bonds and not reaction.ts_species.ts_checks['normal_mode_displacement']:
if not got_expected_changing_bonds and not reaction.ts_species.ts_checks['NMD']:
reaction.ts_species.ts_checks['warnings'] += 'Could not compare normal mode displacement to expected ' \
'breaking/forming bonds due to a missing RMG template; '
reaction.ts_species.ts_checks['normal_mode_displacement'] = True
reaction.ts_species.ts_checks['NMD'] = True
break
if not len(rmg_reactions):
# Just check that some bonds break/form, and that this is not a torsional saddle point.
Expand All @@ -351,7 +363,7 @@ def check_normal_mode_displacement(reaction: 'ARCReaction',
reaction.ts_species.ts_checks['warnings'] += warning + '; '
if any(bond not in dmat_bonds_2 for bond in dmat_bonds_1) \
or any(bond not in dmat_bonds_1 for bond in dmat_bonds_2):
reaction.ts_species.ts_checks['normal_mode_displacement'] = True
reaction.ts_species.ts_checks['NMD'] = True
break
if done:
break
Expand Down
60 changes: 30 additions & 30 deletions arc/checks/ts_test.py
Original file line number Diff line number Diff line change
Expand Up @@ -237,39 +237,39 @@ def test_check_ts(self):
"""Test the check_ts() function."""
self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'freq', 'TS_C3_intraH_8.out')
self.rxn_2a.ts_species.populate_ts_checks()
self.assertFalse(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2a.ts_species.ts_checks['NMD'])
ts.check_ts(reaction=self.rxn_2a, job=self.job1)
self.assertTrue(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(self.rxn_2a.ts_species.ts_checks['NMD'])

self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite', 'keto_enol_ts.out')
self.rxn_7.ts_species.populate_ts_checks()
self.assertFalse(self.rxn_7.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_7.ts_species.ts_checks['NMD'])
ts.check_ts(reaction=self.rxn_7, job=self.job1)
self.assertTrue(self.rxn_7.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(self.rxn_7.ts_species.ts_checks['NMD'])

self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'freq', 'TS_nC3H7-iC3H7.out')
self.rxn_8.ts_species.populate_ts_checks()
self.assertFalse(self.rxn_8.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_8.ts_species.ts_checks['NMD'])
ts.check_ts(reaction=self.rxn_8, job=self.job1)
self.assertTrue(self.rxn_8.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(self.rxn_8.ts_species.ts_checks['NMD'])

def test_did_ts_pass_all_checks(self):
"""Test the did_ts_pass_all_checks() function."""
spc = ARCSpecies(label='TS', is_ts=True)
spc.populate_ts_checks()
self.assertFalse(ts.ts_passed_all_checks(spc))
self.assertFalse(ts.ts_passed_checks(spc))

self.ts_checks = {'E0': False,
'e_elect': False,
'IRC': False,
'freq': False,
'normal_mode_displacement': False,
'NMD': False,
'warnings': '',
}
for key in ['E0', 'e_elect', 'IRC', 'freq']:
spc.ts_checks[key] = True
self.assertFalse(ts.ts_passed_all_checks(spc))
self.assertTrue(ts.ts_passed_all_checks(spc, exemptions=['normal_mode_displacement', 'warnings']))
self.assertFalse(ts.ts_passed_checks(spc))
self.assertTrue(ts.ts_passed_checks(spc, exemptions=['NMD', 'warnings']))
spc.ts_checks['e_elect'] = False

def test_check_rxn_e_elect(self):
Expand Down Expand Up @@ -377,89 +377,89 @@ def test_check_rxn_e0(self):
def test_check_normal_mode_displacement(self):
"""Test the check_normal_mode_displacement() function."""
self.rxn_2a.ts_species.populate_ts_checks()
self.assertFalse(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2a.ts_species.ts_checks['NMD'])
self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS_intra_H_migration_CBS-QB3.out')
self.rxn_2a.determine_family(rmg_database=self.rmgdb)
ts.check_normal_mode_displacement(reaction=self.rxn_2a, job=self.job1)
self.assertTrue(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(self.rxn_2a.ts_species.ts_checks['NMD'])
self.rxn_2a.ts_species.populate_ts_checks()

self.rxn_2b.ts_species.populate_ts_checks()
self.assertFalse(self.rxn_2b.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2b.ts_species.ts_checks['NMD'])
ts.check_normal_mode_displacement(reaction=self.rxn_2b, job=self.job1)
self.assertFalse(self.rxn_2b.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2b.ts_species.ts_checks['NMD'])

# Wrong TS for intra H migration [CH2]CC <=> C[CH]C
self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS_C3_intraH_1.out') # A wrong TS.
self.rxn_2a.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_2a, job=self.job1)
self.assertFalse(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2a.ts_species.ts_checks['NMD'])

self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS_C3_intraH_2.out') # A wrong TS.
self.rxn_2a.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_2a, job=self.job1)
self.assertFalse(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2a.ts_species.ts_checks['NMD'])

self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS_C3_intraH_3.out') # ** The correct TS. **
self.rxn_2a.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_2a, job=self.job1)
self.assertTrue(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(self.rxn_2a.ts_species.ts_checks['NMD'])

self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS_C3_intraH_4.out') # A wrong TS.
self.rxn_2a.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_2a, job=self.job1)
self.assertFalse(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2a.ts_species.ts_checks['NMD'])

self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS_C3_intraH_5.out') # A wrong TS.
self.rxn_2a.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_2a, job=self.job1)
self.assertFalse(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2a.ts_species.ts_checks['NMD'])

self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS_C3_intraH_6.out') # A wrong TS.
self.rxn_2a.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_2a, job=self.job1)
self.assertFalse(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2a.ts_species.ts_checks['NMD'])

self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS_C3_intraH_7.out') # A wrong TS.
self.rxn_2a.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_2a, job=self.job1)
self.assertFalse(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2a.ts_species.ts_checks['NMD'])

self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'freq',
'TS_C3_intraH_8.out') # Correct TS (freq run, not composite).
self.rxn_2a.ts_species.populate_ts_checks()
self.assertFalse(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_2a.ts_species.ts_checks['NMD'])
ts.check_normal_mode_displacement(reaction=self.rxn_2a, job=self.job1)
self.assertTrue(self.rxn_2a.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(self.rxn_2a.ts_species.ts_checks['NMD'])

# CCO[O] + CC <=> CCOO + [CH2]C, incorrect TS:
self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS0_composite_2043.out')
self.rxn_4.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_4, job=self.job1)
self.assertFalse(self.rxn_4.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_4.ts_species.ts_checks['NMD'])

# CCO[O] + CC <=> CCOO + [CH2]C, correct TS:
self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS0_composite_2102.out')
self.rxn_4.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_4, job=self.job1)
self.assertTrue(self.rxn_4.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(self.rxn_4.ts_species.ts_checks['NMD'])

# NCC + H <=> CH3CHNH2 + H2, correct TS:
self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
'TS0_composite_2044.out')
self.rxn_5.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=self.rxn_5, job=self.job1)
self.assertTrue(self.rxn_5.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(self.rxn_5.ts_species.ts_checks['NMD'])

# NH2 + N2H3 <=> NH + N2H4:
self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'freq', 'TS_NH2+N2H3.out')
Expand Down Expand Up @@ -490,7 +490,7 @@ def test_check_normal_mode_displacement(self):
rxn_6.ts_species.mol_from_xyz()
rxn_6.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=rxn_6, job=self.job1)
self.assertTrue(rxn_6.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(rxn_6.ts_species.ts_checks['NMD'])

# [CH2]CC=C <=> CCC=[CH] butylene intra H migration:
self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'composite',
Expand Down Expand Up @@ -526,17 +526,17 @@ def test_check_normal_mode_displacement(self):
rxn_7.ts_species.mol_from_xyz()
rxn_7.ts_species.populate_ts_checks()
ts.check_normal_mode_displacement(reaction=rxn_7, job=self.job1)
self.assertTrue(rxn_7.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(rxn_7.ts_species.ts_checks['NMD'])

@work_in_progress
def test_check_normal_mode_displacement_wip(self):
"""Test the check_normal_mode_displacement() function."""
self.job1.local_path_to_output_file = os.path.join(ts.ARC_PATH, 'arc', 'testing', 'freq',
'TS_NH3+H=NH2+H2.out') # NH3 + H <=> NH2 + H2
self.rxn_3.ts_species.populate_ts_checks()
self.assertFalse(self.rxn_3.ts_species.ts_checks['normal_mode_displacement'])
self.assertFalse(self.rxn_3.ts_species.ts_checks['NMD'])
ts.check_normal_mode_displacement(reaction=self.rxn_3, job=self.job1)
self.assertTrue(self.rxn_3.ts_species.ts_checks['normal_mode_displacement'])
self.assertTrue(self.rxn_3.ts_species.ts_checks['NMD'])

def test_invalidate_rotors_with_both_pivots_in_a_reactive_zone(self):
"""Test the invalidate_rotors_with_both_pivots_in_a_reactive_zone() function."""
Expand Down
2 changes: 1 addition & 1 deletion arc/job/adapters/ts/heuristics.py
Original file line number Diff line number Diff line change
Expand Up @@ -302,7 +302,7 @@ def execute_incore(self):
rxn.ts_species.ts_guesses.append(ts_guess)
save_geo(xyz=xyz,
path=self.local_path,
filename=f'Heuristics {method_index}',
filename=f'Heuristics_{method_index}',
format_='xyz',
comment=f'Heuristics {method_index}, family: {family_label}',
)
Expand Down
Loading

0 comments on commit 36c361d

Please sign in to comment.