arc.species.species

A module for representing stationary points (chemical species and transition states). If the species is a transition state (TS), its ts_guesses attribute will have one or more TSGuess objects.

class arc.species.species.ARCSpecies(adjlist: str = '', bdes: Optional[list] = None, bond_corrections: Optional[dict] = None, charge: Optional[int] = None, checkfile: Optional[str] = None, compute_thermo: Optional[bool] = None, include_in_thermo_lib: Optional[bool] = True, consider_all_diastereomers: bool = True, directed_rotors: Optional[dict] = None, e0_only: bool = False, external_symmetry: Optional[int] = None, fragments: Optional[List[List[int]]] = None, force_field: str = 'MMFF94s', inchi: str = '', is_ts: bool = False, irc_label: Optional[str] = None, label: Optional[str] = None, mol: Optional[Molecule] = None, multiplicity: Optional[int] = None, multi_species: Optional[str] = None, number_of_radicals: Optional[int] = None, occ: Optional[int] = None, optical_isomers: Optional[int] = None, preserve_param_in_scan: Optional[list] = None, rmg_species: Optional[Species] = None, run_time: Optional[timedelta] = None, rxn_label: Optional[str] = None, rxn_index: Optional[int] = None, smiles: str = '', species_dict: Optional[dict] = None, ts_number: Optional[int] = None, xyz: Optional[Union[list, dict, str]] = None, yml_path: Optional[str] = None, keep_mol: bool = False, project_directory: Optional[str] = None)[source]

A class for representing stationary points.

Structures (rotors_dict is initialized in conformers.find_internal_rotors; pivots/scan/top values are 1-indexed):

rotors_dict: {0: {'pivots': ``List[int]``,  # 1-indexed
                  'top': ``List[int]``,  # 1-indexed
                  'scan': ``List[int]``,  # 1-indexed
                  'torsion': ``List[int]``,  # 0-indexed
                  'number_of_running_jobs': ``int``,
                  'success': Optional[``bool``],  # ``None`` by default
                  'invalidation_reason': ``str``,
                  'times_dihedral_set': ``int``,
                  'scan_path': <path to scan output file>,
                  'max_e': ``float``,  # relative to the minimum energy, in kJ/mol,
                  'trsh_counter': ``int``,
                  'trsh_methods': ``List[str]``,
                  'symmetry': ``int``,
                  'dimensions': ``int``,
                  'original_dihedrals': ``list``,
                  'cont_indices': ``list``,
                  'directed_scan_type': ``str``,
                  'directed_scan': ``dict``,  # keys: tuples of dihedrals as strings,
                                              # values: dicts of energy, xyz, is_isomorphic, trsh
                 }
              1: {}, ...
             }
Parameters:
  • label (str, optional) – The species label.

  • is_ts (bool, optional) – Whether the species represents a transition state.

  • rmg_species (Species, optional) – An RMG Species object to be converted to an ARCSpecies object.

  • mol (Molecule, optional) – An RMG Molecule. Atom order corresponds to the order in .initial_xyz

  • xyz (list, str, dict, optional) – Entries are either string-format coordinates, file paths, or ARC’s dict format. (If there’s only one entry, it could be given directly, not in a list). The file paths could direct to either a .xyz file, ARC conformers (w/ or w/o energies), or an ESS log/input files.

  • multiplicity (int, optional) – The species’ electron spin multiplicity. Can be determined from the adjlist/smiles/xyz (If unspecified, assumed to be either a singlet or a doublet).

  • charge (int, optional) – The species’ net charge. Assumed to be 0 if unspecified.

  • smiles (str, optional) – A SMILES representation for the species 2D graph.

  • adjlist (str, optional) – An RMG adjacency list representation for the species 2D graph.

  • inchi (str, optional) – An InChI representation for the species 2D graph.

  • bond_corrections (dict, optional) – The bond additivity corrections (BAC) to be used. Determined from the structure if not directly given.

  • compute_thermo (bool, optional) – Whether to calculate thermodynamic properties for this species.

  • include_in_thermo_lib (bool, optional) – Whether to include in the output RMG library.

  • e0_only (bool, optional) – Whether to only run statmech (w/o thermo) to compute E0.

  • species_dict (dict, optional) – A dictionary to create this object from (used when restarting ARC).

  • yml_path (str, optional) – Path to an Arkane YAML file representing a species (for loading the object).

  • ts_number (int, optional) – An auto-generated number associating the TS ARCSpecies object with the corresponding ARCReaction object.

  • rxn_label (str, optional) – The reaction string (relevant for TSs).

  • rxn_index (int, optional) – The reaction index which is the respective key to the Scheduler rxn_dict.

  • external_symmetry (int, optional) – The external symmetry of the species (not including rotor symmetries).

  • optical_isomers (int, optional) – Whether (=2) or not (=1) the species has chiral center/s.

  • run_time (timedelta, optional) – Overall species execution time.

  • checkfile (str, optional) – The local path to the latest checkfile by Gaussian for the species.

  • number_of_radicals (int, optional) – The number of radicals (inputted by the user, ARC won’t attempt to determine it). Defaults to None. Important, e.g., if a Species is a bi-rad singlet, in which case the job should be unrestricted, but the multiplicity does not have the required information to make that decision (r vs. u).

  • force_field (str, optional) – The force field to be used for conformer screening. The default is MMFF94s. Other optional force fields are MMFF94, UFF, or GAFF (not recommended, slow). If ‘fit’ is specified for this parameter, some initial MMFF94s conformers will be generated, then force field parameters will be fitted for this molecule and conformers will be re-run with the fitted force field (recommended for drug-like species and species with many heteroatoms). Another option is specifying ‘cheap’, and the “old” RDKit embedding method will be used.

  • bdes (list, optional) – Specifying for which bonds should bond dissociation energies be calculated. Entries are bonded atom indices tuples (1-indexed). An ‘all_h’ string entry is also allowed, triggering BDE calculations for all hydrogen atoms in the molecule.

  • directed_rotors (dict) – Execute a directed internal rotation scan (i.e., a series of sp or constrained opt jobs) for the pivots of interest. The optional primary keys are: - ‘brute_force_sp’ - ‘brute_force_opt’ - ‘cont_opt’ - ‘ess’ The brute force methods will generate all the geometries in advance and submit all relevant jobs simultaneously. The continuous method will wait for the previous job to terminate, and use its geometry as the initial guess for the next job. Another set of three keys is allowed, adding _diagonal to each of the above keys. the secondary keys are therefore: - ‘brute_force_sp_diagonal’ - ‘brute_force_opt_diagonal’ - ‘cont_opt_diagonal’ Specifying ‘_diagonal’ will increment all the respective dihedrals together, resulting in a 1D scan instead of an ND scan. Values are nested lists. Each value is a list where the entries are either pivot lists (e.g., [1, 5]) or lists of pivot lists (e.g., [[1, 5], [6, 8]]), or a mix (e.g., [[4, 8], [[6, 9], [3, 4]]]). The requested directed scan type will be executed separately for each list entry in the value. A list entry that contains only two pivots will result in a 1D scan, while a list entry with N pivots will consider all of them, and will result in an ND scan if ‘_diagonal’ is not specified. ARC will generate geometries using the rotor_scan_resolution argument in settings.py Note: An ‘all’ string entry is also allowed in the value list, triggering a directed internal rotation scan for all torsions in the molecule. If ‘all’ is specified within a second level list, then all the dihedrals will be considered together. Currently, ARC does not automatically identify torsions to be treated as ND, and this attribute must be specified by the user. An additional supported key is ‘ess’, in which case ARC will allow the ESS to take care of spawning the ND continuous constrained optimizations (not yet implemented).

  • consider_all_diastereomers (bool, optional) – Whether to consider all different chiralities (tetrahedral carbon centers, nitrogen inversions, and cis/trans double bonds) when generating conformers. True to consider all. If no 3D coordinates are given for the species, all diastereomers will be considered, otherwise the chirality specified by the given coordinates will be preserved.

  • preserve_param_in_scan (list, optional) – Entries are length two iterables of atom indices (1-indexed) between which distances and dihedrals of these pivots must be preserved. Used for identification of rotors which break a TS.

  • fragments (Optional[List[List[int]]]) – Fragments represented by this species, i.e., as in a VdW well or a TS. Entries are atom index lists of all atoms in a fragment, each list represents a different fragment.

  • occ (int, optional) – The number of occupied orbitals (core + val) from a molpro CCSD sp calc.

  • irc_label (str, optional) – The label of an original ARCSpecies object (a TS) for which an IRC job was spawned. The present species object instance represents a geometry optimization job of the IRC result in one direction.

  • project_directory (str, optional) – The path to the project directory.

  • multi_species – (str, optional): The multi-species set this species belongs to. Used for running a set of species simultaneously in a single ESS input file. A species marked as multi_species will only have one conformer considered (n_confs set to 1).

label

The species’ label.

Type:

str

original_label

The species’ label prior to modifications (removing forbidden characters).

Type:

str

multiplicity

The species’ electron spin multiplicity. Can be determined from the adjlist/smiles/xyz (If unspecified, assumed to be either a singlet or a doublet).

Type:

int

charge

The species’ net charge. Assumed to be 0 if unspecified.

Type:

int

number_of_radicals

The number of radicals (inputted by the user, ARC won’t attempt to determine it). Defaults to None. Important, e.g., if a Species is a bi-rad singlet, in which case the job should be unrestricted, but the multiplicity does not have the required information to make that decision (r vs. u).

Type:

int

e_elect

The total electronic energy (without ZPE) at the chosen sp level, in kJ/mol.

Type:

float

e0

The 0 Kelvin energy (total electronic energy plus ZPE) at the chosen sp level, in kJ/mol.

Type:

float

is_ts

Whether the species represents a transition state. True if it does.

Type:

bool

number_of_rotors

The number of potential rotors to scan.

Type:

int

rotors_dict

A dictionary of rotors. structure given below.

Type:

dict

conformers

A list of selected conformers XYZs (dict format).

Type:

list

conformer_energies

A list of conformers E0 (in kJ/mol).

Type:

list

cheap_conformer

A string format xyz of a cheap conformer (not necessarily the best/lowest one).

Type:

str

most_stable_conformer

The index of the best/lowest conformer in self.conformers.

Type:

int

recent_md_conformer

A length three list containing the coordinates of the recent conformer generated by MD, the energy in kJ/mol, and the number of MD runs. Used to detect when the MD algorithm converges on a single structure.

Type:

list

initial_xyz

The initial geometry guess.

Type:

dict

final_xyz

The optimized species geometry.

Type:

dict

_radius

The species radius in Angstrom.

Type:

float

opt_level

Level of theory for geometry optimization. Saved for archiving.

Type:

str

_number_of_atoms

The number of atoms in the species/TS.

Type:

int

mol

An RMG Molecule object used for BAC determination. Atom order corresponds to the order in .initial_xyz

Type:

Molecule

mol_list

A list of localized structures generated from ‘mol’, if possible.

Type:

list

rmg_species

An RMG Species object to be converted to an ARCSpecies object.

Type:

Species

bond_corrections

The bond additivity corrections (BAC) to be used. Determined from the structure if not directly given.

Type:

dict

run_time

Overall species execution time.

Type:

timedelta

t1

The T1 diagnostic parameter from Molpro.

Type:

float

neg_freqs_trshed

A list of negative frequencies this species was troubleshooted for.

Type:

list

compute_thermo

Whether to calculate thermodynamic properties for this species.

Type:

bool

include_in_thermo_lib

Whether to include in the output RMG library.

Type:

bool

e0_only

Whether to only run statmech (w/o thermo) to compute E0.

Type:

bool

thermo

The thermodata calculated by ARC.

Type:

HeatCapacityModel

rmg_thermo

The thermodata generated by RMG for comparison.

Type:

HeatCapacityModel

long_thermo_description

A description for the species entry in the thermo library outputted.

Type:

str

ts_guesses

A list of TSGuess objects for each of the specified methods.

Type:

list

successful_methods

Methods used to generate a TS guess that successfully generated an XYZ guess.

Type:

list

unsuccessful_methods

Methods used to generate a TS guess that were unsuccessfully.

Type:

list

chosen_ts

The TSGuess index corresponding to the chosen TS conformer used for optimization.

Type:

int

chosen_ts_list

The TSGuess index corresponding to the TS guesses that were tried out.

Type:

List[int]

chosen_ts_method

The TS method that was actually used for optimization.

Type:

str

ts_checks

Checks that a TS species went through.

Type:

Dict[str, bool]

rxn_zone_atom_indices

0-indexed atom indices of the active reaction zone.

Type:

List[int]

ts_conf_spawned

Whether conformers were already spawned for the Species (representing a TS) based on its TSGuess objects.

Type:

bool

tsg_spawned

If this species is a TS, this attribute describes whether TS guess jobs were already spawned.

Type:

bool

ts_guesses_exhausted

Whether all TS guesses were tried out with no luck (True if no convergence achieved).

Type:

bool

ts_number

An auto-generated number associating the TS ARCSpecies object with the corresponding ARCReaction object.

Type:

int

ts_report

A description of all methods used for guessing a TS and their ranking.

Type:

str

rxn_label

The reaction string (relevant for TSs).

Type:

str

rxn_index

The reaction index which is the respective key to the Scheduler rxn_dict.

Type:

int

arkane_file

Path to the Arkane Species file generated in processor.

Type:

str

yml_path

Path to an Arkane YAML file representing a species (for loading the object).

Type:

str

keep_mol

Label to prevent the generation of a new Molecule object.

Type:

bool

checkfile

The local path to the latest checkfile by Gaussian for the species.

Type:

str

external_symmetry

The external symmetry of the species (not including rotor symmetries).

Type:

int

optical_isomers

Whether (=2) or not (=1) the species has chiral center/s.

Type:

int

transport_data

A placeholder for updating transport properties after Lennard-Jones calculation (using OneDMin).

Type:

TransportData

force_field

The force field to be used for conformer screening. The default is MMFF94s. Other optional force fields are MMFF94, UFF, or GAFF (not recommended, slow). If ‘fit’ is specified for this parameter, some initial MMFF94s conformers will be generated, then force field parameters will be fitted for this molecule and conformers will be re-run with the fitted force field (recommended for drug-like species and species with many heteroatoms). Another option is specifying ‘cheap’, and the “old” RDKit embedding method will be used.

Type:

str

conf_is_isomorphic

Whether the lowest conformer is isomorphic with the 2D graph representation of the species. True if it is. Defaults to None. If True, an isomorphism check will be strictly enforced for the final optimized coordinates.

Type:

bool

conformers_before_opt

Conformers XYZs of a species before optimization.

Type:

tuple

bdes

Specifying for which bonds should bond dissociation energies be calculated. Entries are bonded atom indices tuples (1-indexed). An ‘all_h’ string entry is also allowed, triggering BDE calculations for all hydrogen atoms in the molecule.

Type:

list

directed_rotors

Execute a directed internal rotation scan (i.e., a series of constrained optimizations). Data is in 3 levels of nested lists, converted from pivots to four-atom torsion indices.

Type:

dict

consider_all_diastereomers

Whether to consider all different chiralities (tetrahydral carbon centers, nitrogen inversions, and cis/trans double bonds) when generating conformers. True to consider all.

Type:

bool, optional

zmat

The species internal coordinates (Z Matrix).

Type:

dict

preserve_param_in_scan

Entries are length two iterables of atom indices (1-indexed) between which distances and dihedrals of these pivots must be preserved.

Type:

list

fragments

Fragments represented by this species, i.e., as in a VdW well or a TS. Entries are atom index lists of all atoms in a fragment, each list represents a different fragment.

Type:

Optional[List[List[int]]]

occ

The number of occupied orbitals (core + val) from a molpro CCSD sp calc.

Type:

int

irc_label

The label of an original ARCSpecies object (a TS) for which an IRC job was spawned. The present species object instance represents a geometry optimization job of the IRC result in one direction. If a species is a transition state, then this attribute contains the labels of the two corresponding “IRC species”, separated by a blank space.

Type:

str

project_directory

The path to the project directory.

Type:

str

multi_species

(str): The multi-species set this species belongs to. Used for running a set of species simultaneously in a single ESS input file.

as_dict(reset_atom_ids: bool = False) dict[source]

A helper function for dumping this object as a dictionary in a YAML file for restarting ARC.

Parameters:

reset_atom_ids (bool, optional) – Whether to reset the atom IDs in the .mol Molecule attribute. Useful when copying the object to avoid duplicate atom IDs between different object instances.

Returns:

The dictionary representation of the object instance.

Return type:

dict

check_xyz_isomorphism(mol: Optional[Molecule] = None, xyz: Optional[dict] = None, allow_nonisomorphic_2d: Optional[bool] = False, verbose: Optional[bool] = True) bool[source]

Check whether the perception of self.final_xyz or xyz is isomorphic with self.mol. If it is not isomorphic, compliant coordinates will be checked (equivalent to checking isomorphism without bond order information, only does not necessitate a molecule object, directly checks bond lengths).

Parameters:
  • mol (Molecule, optional) – A molecule to check instead of self.mol.

  • xyz (dict, optional) – The coordinates to check instead of self.final_xyz.

  • allow_nonisomorphic_2d (bool, optional) – Whether to allow non-isomorphic representations to pass this test.

  • verbose (bool, optional) – Whether to log isomorphism findings and errors.

Returns: bool

Whether the perception of self.final_xyz is isomorphic with self.mol, True if it is.

cluster_tsgs()[source]

Cluster TSGuesses.

copy()[source]

Get a copy of this object instance.

Returns:

A copy of this object instance.

Return type:

ARCSpecies

determine_multiplicity(smiles: str, adjlist: str, mol: Optional[Molecule])[source]

Determine the spin multiplicity of the species.

Parameters:
  • smiles (str) – The SMILES descriptor .

  • adjlist (str) – The adjacency list descriptor.

  • mol (Molecule) – The respective RMG Molecule object.

determine_multiplicity_from_descriptors(smiles: str, adjlist: str, mol: Optional[Molecule])[source]

Determine the spin multiplicity of the species from the chemical descriptors.

Parameters:
  • smiles (str) – The SMILES descriptor .

  • adjlist (str) – The adjacency list descriptor.

  • mol (Molecule) – The respective RMG Molecule object.

determine_multiplicity_from_xyz()[source]

Determine the spin multiplicity of the species from the xyz.

determine_rotors(verbose: bool = False) None[source]

Determine possible unique rotors in the species to be treated as hindered rotors, taking into account all localized structures. The resulting rotors are saved in a {'pivots': [1, 3], 'top': [3, 7], 'scan': [2, 1, 3, 7]} format in self.rotors_dict. Also updates self.number_of_rotors.

determine_symmetry() None[source]

Determine the external symmetry and chirality (optical isomers) of the species.

from_dict(species_dict)[source]

A helper function for loading this object from a dictionary in a YAML file for restarting ARC

from_yml_file(label: Optional[str] = None) bool[source]

Load important species attributes such as label and final_xyz from an Arkane YAML file. Actual QM data parsing is done later when processing thermo and kinetics.

Parameters:

label (str, optional) – The specie label.

Raises:

ValueError – If the adjlist cannot be read.

Returns:

Whether self.mol should be regenerated

Return type:

bool

generate_conformers(n_confs: int = 10, e_confs: float = 5, plot_path: Optional[str] = None) None[source]

Generate conformers.

Parameters:
  • n_confs (int, optional) – The max number of conformers to store in the .conformers attribute that will later be DFT’ed at the conformers_level.

  • e_confs (float, optional) – The energy threshold in kJ/mol above the lowest energy conformer below which all (unique) generated conformers will be stored in the .conformers attribute.

  • plot_path (str, optional) – A folder path in which the plot will be saved. If None, the plot will not be shown (nor saved).

get_cheap_conformer()[source]

Cheaply (limiting the number of possible conformers) get a reasonable conformer, this could very well not be the best (lowest energy) one.

get_xyz(generate: bool = True, return_format: str = 'dict') Optional[Union[dict, str]][source]

Get the highest quality xyz the species has. If it doesn’t have any 3D information, and if generate is True, cheaply generate it. Returns None if no xyz can be retrieved nor generated.

Parameters:
  • generate (bool, optional) – Whether to cheaply generate an FF conformer if no xyz is found. True to generate. If generate is False and the species has no xyz data, the method will return None.

  • return_format (str, optional) – Whether to output a ‘dict’ or a ‘str’ representation of the respective xyz.

Returns:

The xyz coordinates in the requested representation.

Return type:

Optional[Union[dict, str]]

initialize_directed_rotors()[source]

Initialize self.directed_rotors, correcting the input to nested lists and to torsions instead of pivots. Also modifies self.rotors_dict to remove respective 1D rotors dicts if appropriate and adds ND rotors dicts. Principally, from this point we don’t really need self.directed_rotors, self.rotors_dict should have all relevant rotors info for the species, including directed and ND rotors.

Raises:

SpeciesError – If the pivots don’t represent a dihedral in the species.

is_diatomic() Optional[bool][source]

Determine whether the species is diatomic.

Returns:

Whether the species is diatomic.

Return type:

Optional[bool]

is_isomorphic(other: Union[ARCSpecies, Species, Molecule]) Optional[bool][source]

Determine whether the species is isomorphic with other.

Parameters:

other (Union[ARCSpecies, Species, Molecule]) – An ARCSpecies, RMG Species, or RMG Molecule object instance to compare isomorphism with.

Returns:

Whether the species is isomorphic with other.

Return type:

Optional[bool]

is_monoatomic() Optional[bool][source]

Determine whether the species is monoatomic.

Returns:

Whether the species is monoatomic.

Return type:

Optional[bool]

label_atoms()[source]

Labels atoms in order. The label is stored in the atom.label property.

make_ts_report()[source]

A helper function to write content into the .ts_report attribute

mol_from_xyz(xyz: Optional[dict] = None, get_cheap: bool = False) None[source]

Make sure atom order in self.mol corresponds to xyz. Important for TS searches and for identifying rotor indices. This works by generating a molecule from xyz and using the 2D structure to confirm that the perceived molecule is correct. If xyz is not given, the species xyz attribute will be used.

Parameters:
  • xyz (dict, optional) – Alternative coordinates to use.

  • get_cheap (bool, optional) – Whether to generate conformers if the species has no xyz data.

property number_of_atoms

The number of atoms in the species

property number_of_heavy_atoms: int

The number of heavy (non hydrogen) atoms in the species

populate_ts_checks()[source]

Populate (or restart) the .ts_checks attribute with default (None) values.

process_completed_tsg_queue_jobs(yml_path: str)[source]

Process YAML files which are the output of running a TS guess job in the queue.

Parameters:

yml_path (str) – The path to the output YAML file.

process_xyz(xyz_list: Union[list, str, dict])[source]

Process the user’s input and add either to the .conformers attribute or to .ts_guesses.

Parameters:

xyz_list (list, str, dict) – Entries are either string-format, dict-format coordinates or file paths. (If there’s only one entry, it could be given directly, not in a list) The file paths could direct to either a .xyz file, ARC conformers (w/ or w/o energies), or an ESS log/input files, making this method extremely flexible. Internal coordinates (either string or dict) are also allowed and will be converted into cartesian coordinates.

property radius: float

Determine the largest distance from the coordinate system origin attributed to one of the molecule’s atoms in 3D space.

Returns:

The radius in Angstrom.

Return type:

float

scissors(sort_atom_labels: bool = False) list[source]

Cut chemical bonds to create new species from the original one according to the .bdes attribute, preserving the 3D geometry other than the scissioned bond. If one of the scission-resulting species is a hydrogen atom, it will be returned last, labeled as ‘H’. Other species labels will be <original species label>_BDE_index1_index2_X, where “X” is either “A” or “B”, and the indices are 1-indexed.

Parameters:

sort_atom_labels (bool, optional) – Boolean flag, determines whether sorting is required.

Returns: list

The scission-resulting species.

set_dihedral(scan: list, index: int = 1, deg_increment: Optional[float] = None, deg_abs: Optional[float] = None, count: bool = True, xyz: Optional[dict] = None, chk_rotor_list: bool = True)[source]

Set the dihedral angle value of the torsion scan. Either increment by a given value or set to an absolute value. All bonded atoms are rotated as groups. The result is saved to self.initial_xyz.

Parameters:
  • scan (list) – The atom indices representing the dihedral.

  • index (int, optional) – Whether the atom indices are 1-indexed (pass 1) or 0-indexed (pass 0).

  • deg_increment (float, optional) – The dihedral angle increment in degrees.

  • deg_abs (float, optional) – The absolute desired dihedral angle in degrees.

  • count (bool, optional) – Whether to increment the rotor’s times_dihedral_set parameter. True to increment.

  • xyz (dict, optional) – An alternative xyz to use instead of self.final_xyz.

  • chk_rotor_list (bool, optional) – Whether to check if the changing dihedral is in the rotor list.

Raises:
  • InputError – If both deg_increment and deg_abs are None.

  • RotorError – If the rotor could not be identified based on the pivots.

  • TypeError – If deg_increment or deg_abs are of wrong type.

set_mol_list()[source]

Set the .mol_list attribute from self.mol by generating resonance structures, preserving atom order. The mol_list attribute is used for identifying rotors and generating conformers.

set_transport_data(lj_path: str, opt_path: str, bath_gas: str, opt_level: Level, freq_path: Optional[str] = '', freq_level: Optional[Level] = None)[source]

Set the species.transport_data attribute after a Lennard-Jones calculation (via OneDMin).

Parameters:
  • lj_path (str) – The path to a oneDMin job output file.

  • opt_path (str) – The path to an opt job output file.

  • bath_gas (str) – The oneDMin job bath gas.

  • opt_level (Level) – The optimization level of theory.

  • freq_path (str, optional) – The path to a frequencies job output file.

  • freq_level (Level, optional) – The frequencies level of theory.

class arc.species.species.TSGuess(index: Optional[int] = None, method: Optional[str] = None, method_index: Optional[int] = None, method_direction: Optional[str] = None, constraints: Optional[Dict[List[int], int]] = None, t0: Optional[datetime] = None, execution_time: Optional[Union[str, timedelta]] = None, success: Optional[bool] = None, family: Optional[str] = None, xyz: Optional[Union[dict, str]] = None, rmg_reaction: Optional[Reaction] = None, arc_reaction: Optional = None, ts_dict: Optional[dict] = None, energy: Optional[float] = None, cluster: Optional[List[int]] = None, project_directory: Optional[str] = None)[source]

A class for representing TS a guess.

Parameters:
  • index (int, optional) – A running index of all TSGuess objects belonging to an ARCSpecies object.

  • method (str, optional) – The method/source used for the xyz guess.

  • method_index (int, optional) – A sub-index, used for cases where a single method generates several guesses. Counts separately for each direction, ‘F’ and ‘R’.

  • method_direction (str, optional) – The reaction direction used for generating the guess (‘F’ or ‘R’).

  • constraints (dict, optional) – Any constraints to be used when first optimizing this guess (i.e., keeping bond lengths of the reactive site constant).

  • family (str, optional) – The RMG family that corresponds to the reaction, if applicable.

  • success (bool, optional) – Whether the TS guess method succeeded in generating an XYZ guess or not.

  • xyz (dict, str, optional) – The 3D coordinates guess.

  • rmg_reaction (Reaction, optional) – An RMG Reaction object.

  • arc_reaction (ARCReaction, optional) – An ARC Reaction object.

  • ts_dict (dict, optional) – A dictionary to create this object from (used when restarting ARC).

  • energy (float, optional) – Relative energy of all TS conformers in kJ/mol.

  • t0 (datetime.datetime, optional) – Initial time of spawning the guess job.

  • execution_time (datetime.timedelta, optional) – Overall execution time for the TS guess method.

  • project_directory (str, optional) – The path to the project directory.

initial_xyz

The 3D coordinates guess.

Type:

dict

opt_xyz

The 3D coordinates after optimization at the ts_guesses level.

Type:

dict

method

The method/source used for the xyz guess.

Type:

str

method_index

A sub-index, used for cases where a single method generates several guesses. Counts separately for each direction, ‘F’ and ‘R’.

Type:

int

method_direction

The reaction direction used for generating the guess (‘F’ or ‘R’).

Type:

str

family

The RMG family that corresponds to the reaction, if applicable.

Type:

str

rmg_reaction

An RMG Reaction object.

Type:

Reaction

arc_reaction

An ARC Reaction object.

Type:

ARCReaction

t0

Initial time of spawning the guess job.

Type:

datetime.datetime, optional

execution_time

Overall execution time for the TS guess method.

Type:

str

success

Whether the TS guess method succeeded in generating an XYZ guess or not.

Type:

bool

energy

Relative energy of all TS conformers in kJ/mol.

Type:

float

index

A running index of all TSGuess objects belonging to an ARCSpecies object.

Type:

int

imaginary_freqs

The imaginary frequencies of the TS guess after optimization.

Type:

List[float]

conformer_index

An index corresponding to the conformer jobs spawned for each TSGuess object. Assigned only if self.success is True.

Type:

int

successful_irc

Whether the IRS run(s) identified this to be the correct TS by isomorphism of the wells.

Type:

bool

successful_normal_mode

Whether a normal mode check was successful.

Type:

bool

errors

Problems experienced with this TSGuess. Used for logging.

Type:

str

cluster

Indices of TSGuess object instances clustered together.

Type:

List[int]

almost_equal_tsgs(other: TSGuess) bool[source]

Determine whether two TSGuess object instances represent the same geometry.

Parameters:

other (TSGuess) – The other TSGuess object instance to compare to.

Returns:

Whether the two TSGuess object instances represent the same geometry.

Return type:

bool

as_dict(for_report: bool = False) dict[source]

A helper function for dumping this object as a dictionary.

Parameters:

for_report (bool, optional) – Whether to generate a concise dictionary representation for the final_ts_guess_report.

Returns:

The dictionary representation.

Return type:

dict

from_dict(ts_dict: dict)[source]

A helper function for loading this object from a dictionary in a YAML file for restarting ARC

get_xyz(return_format: str = 'dict') Optional[Union[dict, str]][source]

Get the highest quality xyz the TSGuess has. Returns None if no xyz can be retrieved.

Parameters:

return_format (str, optional) – Whether to output a ‘dict’ or a ‘str’ representation of the respective xyz.

Returns:

The xyz coordinates in the requested representation.

Return type:

Optional[Union[dict, str]]

property initial_xyz

The initial coordinate guess

property opt_xyz

The optimized coordinates

process_xyz(xyz: Union[dict, str], project_directory: Optional[str] = None)[source]

Process the user’s input. If xyz represents a file path, parse it.

Parameters:
  • xyz (dict, str) – The coordinates in a dict/string form or a path to a file containing the coordinates.

  • project_directory (str, optional) – The path to the project directory.

Raises:

InputError – If xyz is of the wrong type.

tic()[source]

Initialize self.t0.

tok()[source]

Assign the time difference between now and self.t0 into self.execution_time.

arc.species.species.are_coords_compliant_with_graph(xyz: dict, mol: Molecule) bool[source]

Check whether the Cartesian coordinates represent the same 2D connectivity as the graph. Bond orders are not considered here, this function checks whether the coordinates represent a bond length below 120% of the single bond length of the respective elements.

Parameters:
  • xyz (dict) – The Cartesian coordinates.

  • mol (Molecule) – The 2D graph connectivity information.

Returns:

Whether the coordinates are compliant with the 2D connectivity graph. True if they are.

Return type:

bool

arc.species.species.check_atom_balance(entry_1: Union[dict, str, Molecule], entry_2: Union[dict, str, Molecule], verbose: Optional[bool] = True) bool[source]

Check whether the two entries are in atom balance.

Parameters:
  • entry_1 (Union[dict, str, Molecule]) – Either an xyz (dict or str) or an RMG Molecule object.

  • entry_2 (Union[dict, str, Molecule]) – Either an xyz (dict or str) or an RMG Molecule object.

  • verbose (Optional[bool]) – Whether to log the differences if found.

Raises:

SpeciesError – If both entries are empty.

Returns:

Whether entry1 and entry2 are in atomic balance. True id they are.

Return type:

bool

arc.species.species.check_label(label: str, is_ts: bool = False, verbose: bool = False) Tuple[str, Optional[str]][source]

Check whether a species (or reaction) label is legal, modify it if needed.

Parameters:
  • label (str) – A label.

  • is_ts (bool, optional) – Whether the species label belongs to a transition state.

  • verbose (bool, optional) – Whether to log errors.

Raises:
  • TypeError – If the label is not a string type.

  • SpeciesError – If label is illegal and cannot be automatically fixed.

Returns: Tuple[str, Optional[str]]
  • A legal label.

  • The original label if the label was modified, else None.

arc.species.species.check_xyz(xyz: dict, multiplicity: int, charge: int) bool[source]

Checks xyz for electronic consistency with the spin multiplicity and charge.

Parameters:
  • xyz (dict) – The species coordinates.

  • multiplicity (int) – The species spin multiplicity.

  • charge (int) – The species net charge.

Returns:

Whether the input arguments are all in agreement. True if they are.

Return type:

bool

arc.species.species.colliding_atoms(xyz: dict, mol: Optional[Molecule] = None, threshold: float = 0.55) bool[source]

Check whether atoms are too close to each other. A default threshold of 55% the covalent radii of two atoms is used. For example: - C-O collide at 55% * 1.42 A = 0.781 A - N-N collide at 55% * 1.42 A = 0.781 A - C-N collide at 55% * 1.47 A = 0.808 A - C-H collide at 55% * 1.07 A = 0.588 A - H-H collide at 55% * 0.74 A = 0.588 A

Parameters:
  • xyz (dict) – The Cartesian coordinates.

  • mol (Molecule, optional) – The corresponding Molecule object instance with formal charge information.

  • threshold (float, optional) – The collision threshold to use.

Returns:

True if they are colliding, False otherwise.

Return type:

bool

arc.species.species.cyclic_index_i_minus_1(i: int) int[source]

A helper function for cyclic indexing rotor scans

arc.species.species.cyclic_index_i_plus_1(i: int, length: int) int[source]

A helper function for cyclic indexing rotor scans

arc.species.species.determine_occ(xyz, charge)[source]

Determines the number of occupied orbitals for an MRCI calculation.

Todo

arc.species.species.determine_rotor_symmetry(label: str, pivots: Union[List[int], str], rotor_path: str = '', energies: Optional[Union[list, ndarray]] = None, return_num_wells: bool = False, log: bool = True) Tuple[int, float, Optional[int]][source]

Determine the rotor symmetry number from a potential energy scan. The worst resolution for each peak and valley is determined. The first criterion for a symmetric rotor is that the highest peak and the lowest peak must be within the worst peak resolution (and the same is checked for valleys). A second criterion for a symmetric rotor is that the highest and lowest peaks must be within 10% of the highest peak value. This is only applied if the highest peak is above 2 kJ/mol.

Parameters:
  • label (str) – The species label (used for error messages).

  • pivots (list) – A list of two atom indices representing the pivots.

  • rotor_path (str) – The path to an ESS output rotor scan file.

  • energies (list, optional) – The list of energies in the scan in kJ/mol.

  • return_num_wells (bool, optional) – Whether to also return the number of wells, True to return, default is False.

  • log (bool, optional) – Whether to log info, error, and warning messages.

Raises:
  • InputError – If both or none of the rotor_path and energy arguments are given,

  • or if rotor_path does not point to an existing file.

Returns:

int: The symmetry number float: The highest torsional energy barrier in kJ/mol. int (optional): The number of peaks, only returned if return_len_peaks is True.

Return type:

Tuple[int, float, int]

arc.species.species.determine_rotor_type(rotor_path: str) str[source]

Determine whether this rotor should be treated as a HinderedRotor of a FreeRotor according to its maximum peak.

arc.species.species.enumerate_bonds(mol: Molecule) dict[source]

A helper function for calling Molecule.enumerate_bonds. First, get the Kekulized molecule (get the Kekule version with alternating single and double bonds if the molecule is aromatic), since we don’t have implementation for aromatic bond additivity corrections.

Parameters:

mol (Molecule) – The 2D graph representation of the molecule.

Returns:

Keys are bond types (elements and bond order symbol), values are number of occurrences per bond type.

Return type:

dict

arc.species.species.split_mol(mol: Molecule) Tuple[List[Molecule], List[List[int]]][source]

Split an RMG Molecule object by connectivity gaps while retaining the relative atom order.

Parameters:

mol (Molecule) – The Molecule to split.

Returns:

  • Entries are molecular fragments resulting from the split.

  • Entries are lists with indices that correspond to the original atoms that were assigned to each fragment.

Return type:

Tuple[List[Molecule], List[List[int]]]