Species#
- class autode.species.species.Species(name: str, atoms: List[Atom] | Atoms | None, charge: float | int, mult: float | int, solvent_name: str | None = None)#
- __init__(name: str, atoms: List[Atom] | Atoms | None, charge: float | int, mult: float | int, solvent_name: str | None = None)#
A molecular species. A collection of atoms with a charge and spin multiplicity in a solvent (None is gas phase)
- Parameters:
atoms (list(autode.atoms.Atom) | None) – List of atoms in the species, or None
charge (int) – Charge on the species
mult (int) – Spin multiplicity of the species. 2S+1, where S is the number of unpaired electrons
- Keyword Arguments:
solvent_name (str | None) – Name of the solvent, or None for a species in the gas phase
- property atomic_masses: List[float]#
Atom masses of all the atoms in this species
- property atomic_symbols: List[str]#
Atomic symbols of all atoms in this species
- property bond_matrix: ndarray#
Numpy boolean array containing which atoms are bonded, also known as an adjacency matrix.
- calc_g_cont(*args, **kwargs) None #
Calculate the Gibbs free (G) contribution for this species using Species.calc_thermo()
- calc_h_cont(*args, **kwargs) None #
Calculate the enthalpic (H) contribution for this species using Species.calc_thermo()
- calc_hessian(method: Method, keywords: Keywords | None = None, numerical: bool = False, use_central_differences: bool = False, coordinate_shift: float | ~autode.values.Distance = Distance(0.002 Å), n_cores: int | None = None) None #
Calculate the Hessian
numerical = True
- Parameters:
numerical – Whether to do a numerical frequency calculation using analytic gradients
use_central_differences –
- Use central differences to calculate the
numerical Hessian. If True then use
df/dx = [f(x+h) - f(x-h)] / 2h otherwise use single sided differences (faster but less accurate) df/dx = [f(x+h) - f(x)] / h
coordinate_shift – Shift applied to each Cartesian coordinate (h) in the calculation of the numerical Hessian
n_cores – Number of cores to use for the calculation. If None then default to Config.n_cores
- calc_thermo(method: Method | None = None, calc: ~autode.calculations.calculation.Calculation | None = None, temp: ~autode.values.Temperature | float = Temperature(298.15 kelvin), keywords: ~typing.Sequence[str] | str | None = None, **kwargs) None #
Calculate the free energy and enthalpy contributions using the ideal gas approximation
- Parameters:
calc (autode.calculation.Calculation)
keywords (autode.wrappers.keywords.Keywords)
temp (float | autode.values.Temperature) – Temperature in K
- Keyword Arguments:
lfm_method (LFMethod | str) – Method to treat low frequency modes. {‘igm’, ‘truhlar’, ‘grimme’}. Defaults to Config.lfm_method
ss (str) – Standard state to use. Defaults to Config.standard_state
- Raises:
(autode.exceptions.CalculationException | ValueError) –
See also
autode.thermochemistry.igm.calculate_thermo_cont()
for additional kwargs
- centre() None #
Translate this molecule so the centroid (~COM) is at the origin
- property charge: int#
Total charge on this species
- property conformers: Conformers#
Conformers of this species
- property coordinates: Coordinates | None#
Numpy array of coordinates
- copy() TypeSpecies #
Copy this whole species
- energies#
All energies calculated at a geometry (autode.values.Energies)
- property energy: PotentialEnergy | None#
Last computed potential energy. Setting with a float assumes electornic Hartree units. Example:
>>> import autode as ade >>> species = ade.Species(name='H', atoms=[ade.Atom('H')], charge=0, mult=1) >>> species.energy is None True >>> species.energy = -0.5 >>> species.energy Energy(-0.5 Ha) >>> species.single_point(method=ade.methods.ORCA()) >>> species.energy Energy(-0.50104 Ha)
Energies are instances of autode.values.Energy so can be converted to different units simply:
>>> species.energy.to('kcal mol-1') Energy(-314.40567 kcal mol-1) >>> species.energy.to('eV') Energy(-13.63394 eV)
All previsouly calculated energies of a species are availble with the energies attribute:
>>> species.energies [Energy(-0.5 Ha), Energy(-0.50104 Ha)]
- property enthalpy: Enthalpy | None#
Enthalpy (H) of this species, calculated using the last energy and enthalpy contribution. Example:
>>> import autode as ade >>> h2 = ade.Molecule(smiles='[H][H]') >>> orca = ade.methods.ORCA() >>> >>> h2.optimise(method=orca) >>> h2.calc_h_cont(method=orca) >>> h2.enthalpy Enthalpy(-1.15069 Ha)
The enthalpy contribution is seperated, so performing a single point provides a new enthalpy using the electronic energy at the single-point level of theory:
>>> h2.single_point(method=orca) >>> h2.enthalpy Enthalpy(-1.15497 Ha)
- explicitly_solvate(num: int = 10, solvent: str | Species | None = None) None #
Explicitly solvate this Molecule
- Keyword Arguments:
solvent (str | autode.species.Species | None)
- Raises:
(ValueError) – If the solvent is not defined as a string or a Species and the solvent of this species is not defined
- find_lowest_energy_conformer(lmethod: Method | None = None, hmethod: Method | None = None, allow_connectivity_changes: bool = False) None #
Find the lowest energy conformer of this species. Populates species.conformers and sets species.atoms and species.energy. By default will only optimise at a low-level method
- Keyword Arguments:
hmethod (autode.wrappers.ElectronicStructureMethod) – High-level method to use.
allow_connectivity_changes (bool) – Allow changes in connectivity, although not (by definition) a conformer it is useful to allow
- Raises:
(RuntimeError) – If no conformers (with energies) can be generated
- property formula: str#
Molecular formula of this species. Example:
>>> import autode as ade >>> blank_mol = ade.Molecule() >>> blank_mol.formula '' >>> h2 = ade.Molecule(smiles='[H][H]') >>> h2.formula 'H2'
- property free_energy: FreeEnergy | None#
Free energy (G or A) of this species, calculated using the last energy and free energy contribution
- property frequencies: List[Frequency] | None#
Frequencies from Hessian diagonalisation, in cm-1 by default and are projected from rotation and translation
- property g_cont: FreeEnergyCont | None#
Return the Gibbs (free) contribution to the energy
- property graph: MolecularGraph | None#
Molecular graph with atoms(V) and bonds(E)
Note: Graphs are lazily evaluated, i.e. if one has not been generated for this species before and it does have atoms then a graph will be generated. Subsequent accesses of this property will use the cached internal/private self._graph attribute
- property h_cont: EnthalpyCont | None#
Return the enthalpic contribution to the energy
- has_identical_composition_as(species: Species) bool #
Does this species have the same chemical identity as another?
- property has_reasonable_coordinates: bool#
Does this species have a ‘reasonable’ set of coordinates? I.e. No atom-atom distances that are particularly short or long. Also checks that all the atoms don’t lie in a single plane, which is possible for a failed 3D embedding of a structure.
- has_same_connectivity_as(other: Species) bool #
Determine if this species have the same connectivity as another
- Returns:
Does another species have the same connectivity?
- Return type:
(bool)
- property has_valid_spin_state: bool#
Does this species have a valid spin state given the atomic composition and charge state?
>>> import autode as ade >>> h = ade.Molecule(atoms=[ade.Atom('H')], charge=0, mult=1) >>> h.has_valid_spin_state False >>> hydride = ade.Molecule(atoms=[ade.Atom('H')], charge=-1, mult=1) >>> hydride.has_valid_spin_state True
- property hessian: Hessian | None#
Hessian (d^2E/dx^2) at this geometry (autode.values.Hessian | None) shape = (3*n_atoms, 3*n_atoms)
- property imaginary_frequencies: List[Frequency] | None#
Imaginary frequencies of a molecule
- Imaginary frequencies, or
None if there are none
- Return type:
(list(autode.values.Frequency) | None)
- property is_explicitly_solvated: bool#
- property is_implicitly_solvated: bool#
- is_linear(tol: float | None = None, angle_tol: ~autode.values.Angle = Angle(1.0 °)) bool #
Determine if a species is linear i.e all atoms are colinear
0 to n (n > 1). Present for compatibility and overrides angle_tol if not None
- Keyword Arguments:
angle_tol (autode.values.Angle) – Tolerance on the angle considered to be linear
- is_planar(tol: float | ~autode.values.Distance = Distance(0.0001 Å)) bool #
Determine if a species is planar i.e all atoms are coplanar
- property mult: int#
Total spin multiplicity on this species (2S + 1)
- property n_conformers: int#
Number of conformers of this species
- new_species(name='species', with_constraints: bool = False) Species #
A new version of this species, identical properties without any energies, gradients, hessian, conformers or constraints.
- Parameters:
with_constraints (bool) – Should the constraints from this species be copied into the new one
- Return type:
(autode.species.Species)
- normal_mode(mode_number: int) Coordinates | None #
Vibrational normal mode indexed from 0, the first 6 are translation and rotation and have zero displacements. The first vibrational mode has mode_number = 6.
- Return type:
- optimise(method: Method | None = None, reset_graph: bool = False, calc: Calculation | None = None, keywords: Sequence[str] | str | None = None, n_cores: int | None = None) None #
Optimise the geometry using a method
- Parameters:
reset_graph (bool) – Reset the molecular graph
calc (autode.calculation.Calculation) – Different e.g. constrained optimisation calculation
keywords (list(str) | None) – Calculation keywords to use, if None then use the default for the method. Does not include solvent-specific ones
n_cores (int | None) – Number of cores to use for the calculation, if None then will default to autode.Config.n_cores
- Raises:
- property partial_charges: List[float]#
Partial charges on all the atoms present in this species
- populate_conformers(*args, **kwargs)#
Populate self.conformers
- print_xyz_file(title_line: str | None = None, filename: str | None = None, additional_title_line: str | None = None, with_solvent: bool = True, append: bool = False) None #
Print a standard xyz file from this molecule’s atoms
name of this molecule
- Keyword Arguments:
additional_title_line – Additional elements to add to the title line
with_solvent – If the solvent is explicit then include the solvent atoms in the .xyz file
append – Should the structure be appended to the existing file
- property radius: Distance#
Calculate an approximate radius of this species. Does not consider any VdW radii of the outer most atoms i.e. purely determined on nuclear positions
- reorder_atoms(mapping: dict) None #
Reorder the atoms in this species (in place) using a mapping. For example, to reorder the atoms in a HF molecule:
>>> import autode as ade >>> hf = ade.Species(name='HF', charge=0, mult=1, ... atoms=[ade.Atom('H'), ade.Atom('F', x=1)]) >>> hf.atoms Atoms([Atom(H, 0.0 0.0 0.0), Atom(F, 1.0, 0.0, 0.0)]) >>> hf.reorder_atoms(mapping={0: 1, 1: 0})
- Raises:
(ValueError) – If the mapping is invalid
- reset_graph() None #
Reset the molecular graph of this species by its connectivity
- rotate(axis: ndarray | Sequence, theta: Angle | float, origin: ndarray | Sequence | None = None) None #
Rotate the molecule by around an axis
- Parameters:
theta (Angle | float) – Angle to rotate anticlockwise by if float then assume radian units
origin (np.ndarray | list(float) | None) – Origin of the rotation
- single_point(method: Method, keywords: Keywords | Sequence[str] | str | None = None, n_cores: int | None = None) None #
Calculate the single point energy of the species using a method
- Parameters:
keywords (list(str) | None) – Calculation keywords to use, if None then use the default for the method
n_cores (int | None) – Number of cores to use for the calculation, if None then use autode.Config.n_cores
- Raises:
- property sn: int#
Calculate the symmetry number (σ_R) of the atoms. Only implemented for ‘small’ molecules <50 atoms
References: [1] Theor Chem Account (2007) 118:813 [2] . Phys. Chem. B (2010) 114:16304
- property solvent: Solvent | None#
Solvent which this species is immersed in
- Solvent or None if the species is
in the gas phase
- Return type:
(autode.solvent.Solvent | None)
- property solvent_name: str | None#
Name of the solvent or None
- Return type:
(str | None)
- property sorted_atomic_symbols: List[str]#
Atomic symbols of all atoms sorted alphabetically
- property symmetry_number: int#
Calculate the symmetry number (σ_R) of the atoms. Only implemented for ‘small’ molecules <50 atoms
References: [1] Theor Chem Account (2007) 118:813 [2] . Phys. Chem. B (2010) 114:16304
- translate(vec: Sequence[float]) None #
Translate the molecule by vector