Coordinates#
- class autode.opt.coordinates.base.OptCoordinates(input_array: Sequence | ndarray, units: str | Unit)#
Coordinates used to perform optimisations
- abstract property active_indexes: List[int]#
A list of indexes which are active in this coordinate set
- property active_mol_indexes: List[int]#
Active indexes that are actually atomic coordinates in the molecule
- abstract property cart_proj_g: ndarray | None#
The Cartesian gradient with any constraints projected out
- clear_tensors() None #
Helper function for clearing the energy, gradient and Hessian for these coordinates. Called if the coordinates have been perturbed, making these quantities not accurate any more for the new coordinates
- copy(order='C')#
Return a copy of the array.
- Parameters:
order ({'C', 'F', 'A', 'K'}, optional) – Controls the memory layout of the copy. ‘C’ means C-order, ‘F’ means F-order, ‘A’ means ‘F’ if a is Fortran contiguous, ‘C’ otherwise. ‘K’ means match the layout of a as closely as possible. (Note that this function and
numpy.copy()
are very similar but have different default values for their order= arguments, and this function always passes sub-classes through.)
See also
numpy.copy
Similar function with different default behavior
numpy.copyto
Notes
This function is the preferred method for creating an array copy. The function
numpy.copy()
is similar, but it defaults to using order ‘K’, and will not pass sub-classes through by default.Examples
>>> x = np.array([[1,2,3],[4,5,6]], order='F')
>>> y = x.copy()
>>> x.fill(0)
>>> x array([[0, 0, 0], [0, 0, 0]])
>>> y array([[1, 2, 3], [4, 5, 6]])
>>> y.flags['C_CONTIGUOUS'] True
For arrays containing Python objects (e.g. dtype=object), the copy is a shallow one. The new array will contain the same object which may lead to surprises if that object can be modified (is mutable):
>>> a = np.array([1, 'm', [2, 3, 4]], dtype=object) >>> b = a.copy() >>> b[2][0] = 10 >>> a array([1, 'm', list([10, 3, 4])], dtype=object)
To ensure all elements within an
object
array are copied, use copy.deepcopy:>>> import copy >>> a = np.array([1, 'm', [2, 3, 4]], dtype=object) >>> c = copy.deepcopy(a) >>> c[2][0] = 10 >>> c array([1, 'm', list([10, 3, 4])], dtype=object) >>> a array([1, 'm', list([2, 3, 4])], dtype=object)
- property e: PotentialEnergy | None#
Energy
- property g: ndarray | None#
Gradient of the energy
\[G = \nabla E \equiv \left\{\frac{\partial E}{\partial\boldsymbol{R}_{i}}\right\}\]where \(\boldsymbol{R}\) are a general vector of coordinates.
- property h: ndarray | None#
Hessian (second derivative) matrix of the energy
\[\begin{split}H = \begin{pmatrix} \frac{\partial^2 E} {\partial\boldsymbol{R}_{0}\partial\boldsymbol{R}_{0}} & \cdots \\ \vdots & \ddots \end{pmatrix}\end{split}\]where \(\boldsymbol{R}\) are a general vector of coordinates.
- property h_inv: ndarray | None#
Inverse of the Hessian matrix
\[H^{-1}\]
- h_or_h_inv_has_correct_shape(arr: ndarray | None)#
Does a Hessian or its inverse have the correct shape?
- abstract iadd(value: ndarray) OptCoordinates #
Inplace addition of some coordinates
- abstract property inactive_indexes: List[int]#
A list of indexes which are non-active in this coordinate set
- property indexes: List[int]#
Indexes of the coordinates in this set
- make_hessian_positive_definite() None #
Make the Hessian matrix positive definite by shifting eigenvalues
- property min_eigval: float#
Obtain the minimum eigenvalue of the molecular Hessian in the active space
- Returns:
The minimum eigenvalue
- Return type:
(float)
- abstract property n_constraints: int#
Number of constraints in these coordinates
- abstract property n_satisfied_constraints: int#
Number of constraints that are satisfied in these coordinates
- pred_quad_delta_e(new_coords: ndarray) float #
Calculate the estimated change in energy at the new coordinates based on the quadratic model (i.e. second order Taylor expansion)
- Parameters:
new_coords (np.ndarray) – The new coordinates
- Returns:
The predicted change in energy
- Return type:
(float)
- property raw: ndarray#
Raw numpy array of these coordinates
- property rfo_shift: float#
Get the RFO diagonal shift factor λ for the molecular Hessian that can be applied (H - λI) to obtain the RFO downhill step. The shift is only calculated in active subspace
- Returns:
The shift parameter
- Return type:
(float)
- abstract to(*args, **kwargs) OptCoordinates #
Transformation between these coordinates and another type
- update_g_from_cart_g(arr: Gradient | None) None #
Update the gradient from a Cartesian gradient, zeroing those atoms that are constrained
- update_h_from_cart_h(arr: Hessian | None) None #
Update the Hessian from a cartesian Hessian with shape 3N x 3N for N atoms, zeroing the second derivatives if required
- update_h_from_old_h(old_coords: OptCoordinates, hessian_update_types: List[Type[HessianUpdater]]) None #
Update the Hessian \(H\) from an old Hessian using an update scheme. Requires the gradient to be set, and the old set of coordinates with gradient to be available
- Parameters:
old_coords (OptCoordinates) – Old set of coordinates with gradient and hessian defined
hessian_update_types (list[type[HessianUpdater]]) – A list of hessian updater classes - the first updater that meets the mathematical conditions will be used
- class autode.opt.coordinates.cartesian.CartesianCoordinates(input_array, units='Å')#
Flat Cartesian coordinates shape = (3 × n_atoms, )
- property active_indexes: List[int]#
A list of indexes which are active in this coordinate set
- property cart_proj_g: ndarray | None#
The Cartesian gradient with any constraints projected out
- property expected_number_of_dof: int#
Expected number of degrees of freedom for the system
- iadd(value: ndarray) OptCoordinates #
Inplace addition of some coordinates
- property inactive_indexes: List[int]#
A list of indexes which are non-active in this coordinate set
- property n_constraints: int#
Number of constraints in these coordinates
- property n_satisfied_constraints: int#
Number of constraints that are satisfied in these coordinates
- to(value: str) OptCoordinates #
Transform between cartesian and internal coordinates e.g. delocalised internal coordinates or other units
- Returns:
Transformed coordinates
- Return type:
(autode.opt.coordinates.OptCoordinates)
- Raises:
(ValueError) – If the conversion cannot be performed
Delocalised internal coordinate implementation from: 1. https://aip.scitation.org/doi/pdf/10.1063/1.478397 and references cited therein. Also used is 2. https://aip.scitation.org/doi/pdf/10.1063/1.1515483
The notation follows the paper and is briefly summarised below:
- class autode.opt.coordinates.dic.DIC(input_array)#
Delocalised internal coordinates (DIC)
- property active_indexes: List[int]#
A list of indexes for the active modes in this coordinate set
- property cart_proj_g: ndarray | None#
The Cartesian gradient with any constraints projected out
- classmethod from_cartesian(x: CartesianCoordinates, primitives: PIC | None = None) DIC #
Convert cartesian coordinates to primitives then to delocalised internal coordinates (DICs), of which there should be 3N-6 for a polyatomic system with N atoms
all pairwise inverse distances
- Returns:
Delocalised internal coordinates
- Return type:
(autode.opt.coordinates.DIC)
- iadd(value: ndarray) OptCoordinates #
Set some new internal coordinates and update the Cartesian coordinates
\[x^(k+1) = x(k) + ({B^T})^{-1}(k)[s_{new} - s(k)]\]for an iteration k.
- Raises:
(RuntimeError) – If the transformation diverges
- property inactive_indexes: List[int]#
A list of indexes for the non-active modes in this coordinate set
- to(value: str) OptCoordinates #
Convert these DICs to another type of coordinate
- Returns:
Coordinates
- Return type:
(autode.opt.coordinates.OptCoordinates)
- class autode.opt.coordinates.dic.DICWithConstraints(input_array)#
Delocalised internal coordinates (DIC) with constraints. Uses Lagrangian multipliers to enforce the constraints with:
..math:
L(X, λ) = E(s) + \sum_{i=1}^m \lambda_i C_i(X)
where s are internal coordinates, and C the constraint functions. The optimisation space then is the n non-constrained internal coordinates and the m Lagrangian multipliers (lambda_i).
- property active_indexes: List[int]#
Generate a list of indexes for the active modes in this coordinate set
- property cart_proj_g: ndarray | None#
Obtain Cartesian gradient with constraints projected out
- classmethod from_cartesian(x: CartesianCoordinates, primitives: PIC | None = None) DICWithConstraints #
Generate delocalised internal coordinates with constraints with the Lagrangian multipliers initialised as zeroes
- Parameters:
x – Cartesian coordinates
primitives – Primitive internal coordinates. If undefined then use all pairwise inverse distances
- Returns:
DIC with constraints
- Return type:
- property g#
Gradient of the energy, contains the Lagrangian dL/d_λi terms where λi is the i-th lagrangian multiplier.
- property h#
The Hessian matrix, containing Lagrangian constraint terms
- Return type:
(np.ndarray)
- iadd(value: ndarray) OptCoordinates #
Add a step in internal coordinates (along with Lagrange multipliers) to this set of coordinates, and update the Cartesian coordinates
- Parameters:
value – Difference between current and new DICs, and the multipliers
- property inactive_indexes: List[int]#
Generate a list of mode indexes that are inactive in the optimisation space. This requires the m constrained modes being at the end of the coordinate set. It also includes the lagrange multipliers
- property raw: ndarray#
Raw numpy array of these coordinates including the multipliers
- update_lagrange_multipliers(arr: ndarray) None #
Update the lagrange multipliers by adding a set of values
- class autode.opt.coordinates.primitives.CompositeBonds(bonds: List[Tuple[int, int]], coeffs: List[float])#
Linear Combination of several bond distances
- __init__(bonds: List[Tuple[int, int]], coeffs: List[float])#
Linear combination of a list of bonds and the corresponding coefficients given as a list of real numbers
- Parameters:
bonds – A list of tuples (i, j) representing bonds
coeffs – A list of floating point coefficients in order
- class autode.opt.coordinates.primitives.ConstrainedCompositeBonds(bonds: List[Tuple[int, int]], coeffs: List[float], value: float)#
Constrained linear combindation of bonds
- __init__(bonds: List[Tuple[int, int]], coeffs: List[float], value: float)#
Linear combination of a list of bonds and the corresponding coefficients given as a list of real numbers
- Parameters:
bonds – A list of tuples (i, j) representing bonds
coeffs – A list of floating point coefficients in order
value – The target value for this coordinate
- class autode.opt.coordinates.primitives.ConstrainedPrimitive(*atom_indexes: int)#
A primitive internal coordinate constrained to a value
- delta(x: CartesianCoordinates) float #
Difference between the observed and required value
- is_constrained = True#
- is_satisfied(x: CartesianCoordinates, tol: float = 0.0001) bool #
Is this constraint satisfied to within an absolute tolerance
- class autode.opt.coordinates.primitives.ConstrainedPrimitiveBondAngle(m: int, o: int, n: int, value: float)#
- __init__(m: int, o: int, n: int, value: float)#
Angle (m-o-n) constrained to a value (in radians)
- Parameters:
n – Atom index
value – Required value of the constrained angle
- class autode.opt.coordinates.primitives.ConstrainedPrimitiveDistance(i: int, j: int, value: float)#
- __init__(i: int, j: int, value: float)#
Distance constrained to a value
- Parameters:
value – Required value of the constrained distance
- class autode.opt.coordinates.primitives.LinearAngleBase(m: int, o: int, n: int, r: int, axis: LinearBendType)#
- __init__(m: int, o: int, n: int, r: int, axis: LinearBendType)#
Linear Bend: m-o-n
- class autode.opt.coordinates.primitives.LinearBendType(value)#
For linear angles, there are two orthogonal directions
- BEND = 0#
- COMPLEMENT = 1#
- class autode.opt.coordinates.primitives.Primitive(*atom_indexes: int)#
Primitive internal coordinate
- __init__(*atom_indexes: int)#
A primitive internal coordinate that involves a number of atoms
- derivative(x: CartesianCoordinates) ndarray #
Calculate the derivatives with respect to cartesian coordinates
\[\frac{dq} {d\boldsymbol{X}_{i, k}} {\Bigg\rvert}_{X=X0}\]where \(q\) is the primitive coordinate and \(\boldsymbol{X}\) are the cartesian coordinates.
- Returns:
Derivative array of shape (N, )
- Return type:
(np.ndarray)
- is_constrained = False#
- second_derivative(x: CartesianCoordinates) ndarray #
Calculate the second derivatives with respect to cartesian coordinates
\[\frac{d^2 q} {d\boldsymbol{X}_{i, k}^2} {\Bigg\rvert}_{X=X0}\]where \(q\) is the primitive coordinate and \(\boldsymbol{X}\) are the cartesian coordinates.
- Returns:
Second derivative matrix of shape (N, N)
- Return type:
(np.ndarray)
- class autode.opt.coordinates.primitives.PrimitiveBondAngle(m: int, o: int, n: int)#
Bond angle between three atoms, calculated with the arccosine of the normalised dot product
- __init__(m: int, o: int, n: int)#
Bond angle m-o-n
- class autode.opt.coordinates.primitives.PrimitiveDihedralAngle(m: int, o: int, p: int, n: int)#
- __init__(m: int, o: int, p: int, n: int)#
Dihedral angle: m-o-p-n
- class autode.opt.coordinates.primitives.PrimitiveDistance(i: int, j: int)#
Distance between two atoms:
\[q = |\boldsymbol{X}_i - \boldsymbol{X}_j|\]
- class autode.opt.coordinates.primitives.PrimitiveDummyLinearAngle(m: int, o: int, n: int, axis: LinearBendType)#
Linear bend with a dummy atom
- __init__(m: int, o: int, n: int, axis: LinearBendType)#
Linear Bend: m-o-n
- class autode.opt.coordinates.primitives.PrimitiveImproperDihedral(m: int, o: int, p: int, n: int)#
Out-of-Plan (improper) dihedral angles
- class autode.opt.coordinates.primitives.PrimitiveInverseDistance(i: int, j: int)#
Inverse distance between two atoms:
\[q = \frac{1} {|\boldsymbol{X}_i - \boldsymbol{X}_j|}\]
- class autode.opt.coordinates.primitives.PrimitiveLinearAngle(m: int, o: int, n: int, r: int, axis: LinearBendType)#
Linear Angle w.r.t. a reference atom
Coordinates for dimer optimisations. Notation follows: [1] https://aip.scitation.org/doi/pdf/10.1063/1.2815812
- class autode.opt.coordinates.dimer.DimerCoordinates(input_array: Sequence | ndarray, units: str | Unit = 'Å amu^1/2')#
Mass weighted Cartesian coordinates for two points in space forming a dimer, such that the midpoint is close to a first order saddle point
- property delta: float#
Distance between the dimer point, Δ
- property did_rotation#
Rotated this iteration?
- property did_translation#
Translated this iteration?
- property dist: MWDistance#
Distance that the dimer was translated by from its last position
- property f_r: ndarray#
Rotational force F_R. eqn. 3 in ref. [1]
- property f_t: ndarray#
Translational force F_T, eqn. 2 in ref. [1]
- classmethod from_species(species1: Species, species2: Species) DimerCoordinates #
Initialise a set of DimerCoordinates from two species, i.e. those either side of the saddle point.
- property g0: ndarray#
Gradient at the midpoint of the dimer
- property g1: ndarray#
Gradient on the ‘left’ side of the dimer
- property g2: ndarray#
Gradient on the ‘right’ side of the dimer
- g_at(point: DimerPoint) ndarray #
- iadd(value: ndarray) OptCoordinates #
Inplace addition of some coordinates
- set_g_at(point: DimerPoint, arr: ndarray, mass_weighted: bool = True)#
Set the gradient vector at a particular point
- property tau: ndarray#
Direction between the two ends of the dimer (τ)
- property tau_hat: ndarray#
Normalised direction between the two ends of the dimer
- to(*args, **kwargs) OptCoordinates #
Transformation between these coordinates and another type
- property x0: ndarray#
Midpoint of the dimer
- property x1: ndarray#
Coordinates on the ‘left’ side of the dimer
- property x2: ndarray#
Coordinates on the ‘right’ side of the dimer
- x_at(point: DimerPoint, mass_weighted: bool = True) ndarray #
Coordinates at a point in the dimer