5. Coordinate Systems

Choice of coordinate system is often crucial, for the successful outcome of geometry optimizations. Pysisyphus supports different coordinate systems. By default, molecular optimizations are carried out in redundant internal coordinates (RIC).

5.1. Supported Coordinate Systems

5.1.1. Cartesian Coordinates

  • Strong coupling

  • Unambiguously defined, if translation and rotation (TR) are removed

  • Redundant set, if TR are not removed

  • type: cart

5.1.2. Redundant Internal Coordinates

  • Comprising bond stretches, (linear) bends, (improper) torsions and other primitives

  • Less coupling, compared to Cartesians

  • Easy estimation of diagonal model Hessians

  • Support constraints

  • Require sophisticated setup algorithm

  • Iterative internal-Cartesian backtransformation, which may fail

  • Usually highly redundant set

  • type: redund

5.1.3. Delocalized Internal Coordinates (DLC)

  • Complicated linear combinations of primitive internals

  • No coupling, coordinates are orthogonal to each other, at least at the geometrey they were defined at

  • Non redundant set

  • More efficient compared to RIC for bigger systems (if initial DLC generation is feasible)

  • Same comments apply, as for RICs

  • type: dlc

5.1.4. Translation & Rotation Internal Coordinates (TRIC)

  • Especially suited to optimize solvated/non-covalently bound systems

  • Translation and rotation coordinates are assigned to every fragment

  • Avoids error-prone assignment of interfragment coordinates

  • See 10.1063/1.4952956 for a full discussion

  • By default the B-Matrix is recalculated in every step of the internal-Cartesian backtransformation when TRIC is enabled

  • type: tric

5.2. Supported File Formats

All formats can be read, at least to a certain extent. Coordinate files are expected to be in Å, if not otherwise dictated by the respective format.

Suffix

Write

Comment

.xyz

Plain XYZ, single geometry.

.trj

Plain XYZ, multiple geometries.

.pdb

Protein Data Bank.

.molden

Restricted to [Geometries] block.

.zmat

Z-Matrix, see below for an example.

.cjson

As saved by Avogadro.

.crd

CHARMM card format

.sdf

Structure-data file

5.2.1. Z-Matrix example

C
N 1 1.35
H 1 1.0 2 105.0
H 2 1.4 1 105.0 3 150.0
H 2 1.4 1 110.0 3 -160.0

Indexing in Z-matrices is 1-based, so the first atom has index 1. Variable substitution is not supported in the current Z-matrix parser.

5.3. YAML Input

See below for an explanation of possible inputs in the geom section. Input related to internal coordinates is given mostly in the coord_kwargs subgroup, whereas atom-related input (isotops, freeze_atoms) is given one level above. See the example below:

geom:
 type: cart                     # Coordinate system (cart/redund/dlc/tric)
 fn: [input]                    # File name or inline input
 union: False                   # Define same set of primitives at multiple geometries
 isotopes: null                 # Specify different isotopes
 freeze_atoms: null             # Freeze Cartesians of certain atoms
 coord_kwargs:                  # Keywords that are passed to the internal coordinate class
  define_prims: null            # Additionally define these primitives
  constrain_prims: null         # Primitive internals to be constrained
  freeze_atoms_exclude: False   # Whether to set up internal coordinates for frozen atoms
preopt:
 geom:                          # geom block in preopt takes same keywords as above
  ...                           # no 'fn' key here!
endopt:
 geom:                          # geom block in endopt takes same keywords as above
  ...                           # no 'fn' key here!

Employed coordinates and coordinate systems in preopt and endopt are similary controlled by a geom block. Same keywords are supported, as for the geom block, at the top level, except the fn key.

5.3.1. Inline input

Inline coordinates are supported for XYZ format, and expected in Å. Take care of proper indentation. The example below would yield RIC for the hydrogen molecule (1 bond).

geom:
 type: redund
 fn: |
  2

  H 0.0 0.0 0.0
  H 0.0 0.0 0.7

5.3.2. Types of Primitive Coordinates

Pysisyphus implements many different (primitive) internal coordinates. Every coordinate is defined by its type and a set of atom indices, e.g., 2 indices for a bond, 3 indices for a bend and 4 indices for a dihedral.

Specification of a type is necessary, as there are many different kinds of bonds, bends and dihedrals/out-of-plane. One can't just assume, that a coordinate comprised of 3 atom indices is always a regular bend, as it may also be a linear bend or a translational coordinate (TRANSLATION_X, 14), describin the mean Cartesian X coordinate of 3 atoms.

Atom indices start at 0!

# Primitive types
BOND = 0
AUX_BOND = 1
HYDROGEN_BOND = 2
INTERFRAG_BOND = 3
AUX_INTERFRAG_BOND = 4
BEND = 5
LINEAR_BEND = 6
LINEAR_BEND_COMPLEMENT = 7
PROPER_DIHEDRAL = 8
IMPROPER_DIHEDRAL = 9
OUT_OF_PLANE = 10
LINEAR_DISPLACEMENT = 11
LINEAR_DISPLACEMENT_COMPLEMENT = 12
# TRANSLATION = 13  # Dummy coordinate
TRANSLATION_X = 14
TRANSLATION_Y = 15
TRANSLATION_Z = 16
# ROTATION = 17  # Dummy coordinate
ROTATION_A = 18
ROTATION_B = 19
ROTATION_C = 20
# CARTESIAN = 21  # Dummy coordinate
CARTESIAN_X = 22
CARTESIAN_Y = 23
CARTESIAN_Z = 24
BONDED_FRAGMENT = 25
DUMMY_TORSION = 26
DISTANCE_FUNCTION = 27
# atan2 based coordinates
BEND2 = 28
TORSION2 = 29

As some of these types are quite unwieldy, several shortcuts are supported, that can be used in place of the types above.

# Additional shortcuts
# Using Cartesians in the framework of internal coordinates is mainly
# useful if one wants to constrain certain atoms.
"X": [PT.CARTESIAN_X],
"Y": [PT.CARTESIAN_Y],
"Z": [PT.CARTESIAN_Z],
"XY": [PT.CARTESIAN_X, PT.CARTESIAN_Y],
"XZ": [PT.CARTESIAN_X, PT.CARTESIAN_Z],
"YZ": [PT.CARTESIAN_Y, PT.CARTESIAN_Z],
"XYZ": [PT.CARTESIAN_X, PT.CARTESIAN_Y, PT.CARTESIAN_Z],
"ATOM": [PT.CARTESIAN_X, PT.CARTESIAN_Y, PT.CARTESIAN_Z],
# Primitive aliases
"B": [PT.BOND],
"A": [PT.BEND],
"A2": [PT.BEND2],
"D": [PT.PROPER_DIHEDRAL],
"DIHEDRAL": [PT.PROPER_DIHEDRAL],
"TORSION": [PT.PROPER_DIHEDRAL],
"D2": [PT.PROPER_DIHEDRAL2],
"DIHEDRAL2": [PT.PROPER_DIHEDRAL2],
"TORSION2": [PT.PROPER_DIHEDRAL2],
# Translation & Rotation coordinates
"TRANSLATION": [PT.TRANSLATION_X, PT.TRANSLATION_Y, PT.TRANSLATION_Z],
"ROTATION": [PT.ROTATION_A, PT.ROTATION_B, PT.ROTATION_C],
# Miscellaneous
"DIST_FUNC": [PT.DISTANCE_FUNCTION]

5.3.3. Define Additional Primitives

Pysisyphus tries its best, to automatically come up with a reasonable set of internal coordinates, but sometimes the algorithm misses an important one. Especially at transition state guesses, where increased atom distances are common, bonds may be missed.

In such cases, additional coordinates can be requested explicitly. If additional coordinates are requested, a nested list is expected [[coord0], [coord1], ...].

# General structure (list of coordinate lists)
define_prims: [[PrimType or Shortcut], *[atom indices], ...]

# Examples

# Additional bond between atoms 4 and 7 (0-based indexing).
# All three lines below result in the same bond; the latter two use shortcuts.
define_prims: [[0, 4, 7]]
define_prims: [[B, 4, 7]]
define_prims: [[BOND, 4, 7]]

# Wrong specification (forgot outer list/brackets):
define_prims: [0, 4, 7]

# Also define an additional dihedral, beside the bond
define_prims: [[0, 4, 7], ["D", 0, 1, 2, 3]]

5.3.4. Freeze Atoms

All three Cartesian coordinates (X, Y, Z) of certain atoms can be frozen, so they always remain at their initial value. By setting freeze_atoms_exclude in coord_kwargs, frozen atoms can be excluded from the internal coordinate setup. By default frozen atoms are included.

geom:
 type: [type]
 fn: [fn]
 freeze_atoms: [*atom indices]
 coord_kwargs:
  freeze_atoms_exclude: False

# Example; fully freeze Cartesians of first and second atom.
geom:
 type: cart
 fn: input.xyz
 freeze_atoms: [0, 1]
 coord_kwargs:
  freeze_atoms: True

5.3.5. Constraints

Constraints beyond frozen atoms are currently only supported in conjunction with RIC (`type: redund`). It is not (yet) possible to modify the value of the specified coordinate via YAML input; the internal coordinate is constrained at its initial value. The same syntax as for define_prims is used. If the coordinate of the requested constraint is not already defined, it will be defined subsequently. There is no need to also add the constrained coordinate to define_prims.

# General structure (nested list of coordinates)
constrain_prims: [[[PrimType or Shortcut], *[atom indices]], ...]

# Examples

# Constrain Cartesian coordinate of atom 0.
# Both lines result in the same constraint.
constrain_prims: [[XYZ, 0]]
constrain_prims: [[ATOM, 0]]

# Constrain only Cartesian X and Y component of atom 0.
constrain_prims: [[XY, 0]]

# Constraint bond between atoms 4 and 7 (0-based indexing).
# All three lines below result in the same constraint; the latter two use shortcuts.
constrain_prims: [[0, 4, 7]]
constrain_prims: [[B, 4, 7]]
constrain_prims: [[BOND, 4, 7]]

Constraining the Cartesian coordinates (X, Y and Z) of one atom does not affect the final energy of an optimization. But constraining more than one atome does.

Harmonic restraints to selected primitive internals can be specified in the calc: section (see the Restraint documentation).

5.3.6. Isotopes

Different isotope masses can be requested. The system works similar to Gaussians system. A list of pairs is expected, where the first number specifies the atom and the second number is either an integer or a float. If it is an integer, the isotope mass closest to this integer is looked up in an internal database. Floats are used as is.

# General structure (nested list of coordinates)
isotopes: [[[atom index], [new mass, integer/float], ...]

# Modify the mass of atom with index 2 (hydrogen in this case)
# Both lines give identical results (deuterium).
# In the second line, the mass is given directly.
isotopes: [[2, 2]]
isotopes: [[2, 2.014101778]]

Different isotope masses affect calculated frequencies and IRCs. Atoms can be fixed in IRC calculations by specifying a very high mass. YAML does not recognize 1e9 as float, take care to add a dot (1.e9).

# Fix atom 0 in IRC calculation.
isotopes: [[0, 1.e9]]

5.4. Geometry & RedundantCoords

The central class in pysisyphus, handling coordinates and delegating calculations to (external QC) codes, is the Geometry class, similar to ASE's Atoms.

class pysisyphus.Geometry.Geometry(atoms, coords, fragments=None, coord_type='cart', coord_kwargs=None, isotopes=None, freeze_atoms=None, comment='', name='')[source]
__init__(atoms, coords, fragments=None, coord_type='cart', coord_kwargs=None, isotopes=None, freeze_atoms=None, comment='', name='')[source]

Object representing atoms in a coordinate system.

The Geometry represents atoms and their positions in coordinate system. By default cartesian coordinates are used, but internal coordinates are also possible.

Parameters:
  • atoms (iterable) -- Iterable of length N, containing element symbols.

  • coords (1d iterable) -- 1d iterable of length 3N, containing the cartesian coordinates of N atoms.

  • fragments (dict, optional) -- Dict with different keys denoting different fragments. The values contain lists of atom indices.

  • coord_type ({"cart", "redund"}, optional) -- Type of coordinate system to use. Right now cartesian (cart) and redundand (redund) are supported.

  • coord_kwargs (dict, optional) -- Dictionary containing additional arguments that get passed to the constructor of the internal coordinate class.

  • isotopes (iterable of pairs, optional) -- Iterable of pairs consisting of 0-based atom index and either an integer or a float. If an integer is given the closest isotope mass will be selected. Given a float, this float will be directly used as mass.

  • freeze_atoms (iterable of integers) -- Specifies which atoms should remain fixed at their initial positions.

  • comment (str, optional) -- Comment string.

  • name (str, optional) -- Verbose name of the geometry, e.g. methanal or water. Used for printing

align_principal_axes()[source]

Align the principal axes to the cartesian axes.

https://math.stackexchange.com/questions/145023

property all_energies

Return energies of all states that were calculated.

This will also set self.energy, which may NOT be the ground state, but the state correspondig to the 'root' attribute of the calculator.

approximate_radius()[source]

Approximate molecule radius from the biggest atom distance along an axis.

as_ase_atoms()[source]
as_g98_list()[source]

Returns data for fake Gaussian98 standard orientation output.

Returns:

g98_list -- List with one row per atom. Every row contains [center number, atomic number, atomic type (always 0 for now), X Y Z coordinates in Angstrom.

Return type:

list

as_xyz(comment='', atoms=None, cart_coords=None)[source]

Current geometry as a string in XYZ-format.

Parameters:
  • comment (str, optional) -- Will be written in the second line (comment line) of the XYZ-string.

  • cart_coords (np.array, 1d, shape (3 * atoms.size, )) -- Cartesians for dumping instead of self._coords.

Returns:

xyz_str -- Current geometry as string in XYZ-format.

Return type:

str

assert_cart_coords(coords)[source]
assert_compatibility(other)[source]

Assert that two Geometries can be substracted from each other.

Parameters:

other (Geometry) -- Geometry for comparison.

atom_indices()[source]

Dict with atom types as key and corresponding indices as values.

Returns:

inds_dict -- Unique atom types as keys, corresponding indices as values.

Return type:

dict

property atom_types
atom_xyz_iter()[source]
property atomic_numbers
property bond_sets
calc_double_ao_overlap(geom2)[source]
calc_energy_and_forces()[source]

Force a calculation of the current energy and forces.

property cart_coords
property cart_forces
property cart_gradient
property cart_hessian
center()[source]
property center_of_mass

Returns the center of mass.

Returns:

R -- Center of mass.

Return type:

np.array, shape(3, )

center_of_mass_at(coords3d)[source]

Returns the center of mass at given coords3d.

Parameters:

coords3d (np.array, shape(N, 3)) -- Cartesian coordiantes.

Returns:

R -- Center of mass.

Return type:

np.array, shape(3, )

property centroid

Geometric center of the Geometry.

Returns:

R -- Geometric center of the Geometry.

Return type:

np.array, shape(3, )

clear()[source]

Reset the object state.

property comment
coord_types = {'cart': None, 'cartesian': <class 'pysisyphus.intcoords.CartesianCoords.CartesianCoords'>, 'dlc': <class 'pysisyphus.intcoords.DLC.DLC'>, 'hdlc': <class 'pysisyphus.intcoords.DLC.HDLC'>, 'hredund': <class 'pysisyphus.intcoords.RedundantCoords.HybridRedundantCoords'>, 'mwcartesian': <class 'pysisyphus.intcoords.CartesianCoords.MWCartesianCoords'>, 'redund': <class 'pysisyphus.intcoords.RedundantCoords.RedundantCoords'>, 'tmtric': <class 'pysisyphus.intcoords.RedundantCoords.TMTRIC'>, 'tric': <class 'pysisyphus.intcoords.RedundantCoords.TRIC'>}
property coords

1d vector of atomic coordinates.

Returns:

coords -- 1d array holding the current coordinates.

Return type:

np.array

property coords3d

Coordinates in 3d.

Returns:

coords3d -- Coordinates of the Geometry as 2D array.

Return type:

np.array

property coords_by_type

Coordinates in 3d by atom type and their corresponding indices.

Returns:

  • cbt (dict) -- Dictionary with the unique atom types of the Geometry as keys. It's values are the 3d coordinates of the corresponding atom type.

  • inds (dict) -- Dictionary with the unique atom types of the Geometry as keys. It's values are the original indices of the 3d coordinates in the whole coords3d array.

copy(coord_type=None, coord_kwargs=None)[source]

Returns a new Geometry object with same atoms and coordinates.

Parameters:
  • coord_type (str) -- Desired coord_type, defaults to current coord_type.

  • coord_kwargs (dict, optional) -- Any desired coord_kwargs that will be passed to the RedundantCoords object.

Returns:

geom -- New Geometry object with the same atoms and coordinates.

Return type:

Geometry

copy_all(coord_type=None, coord_kwargs=None)[source]
property covalent_radii
del_atoms(inds, **kwargs)[source]
describe()[source]
dump_xyz(fn, cart_coords=None, **kwargs)[source]
eckart_projection(mw_hessian, return_P=False, full=False)[source]
property energy

Energy of the current atomic configuration.

Returns:

energy -- Energy of the current atomic configuration.

Return type:

float

fd_coords3d_gen(step_size=0.001)[source]

Iterator returning 3d Cartesians for finite-differences.

property forces

Energy of the current atomic configuration.

Returns:

force -- 1d array containing the forces acting on the atoms. Negative of the gradient.

Return type:

np.array

get_energy_and_cart_forces_at(cart_coords)[source]
get_energy_and_cart_hessian_at(cart_coords)[source]
get_energy_and_forces_at(coords)[source]

Calculate forces and energies at the given coordinates.

The results are not saved in the Geometry object.

get_energy_at(coords)[source]
get_energy_at_cart_coords(cart_coords)[source]
get_fragments(regex)[source]
get_imag_frequencies(hessian=None, thresh=1e-06)[source]
get_normal_modes(cart_hessian=None, full=False)[source]

Normal mode wavenumbers, eigenvalues and Cartesian displacements Hessian.

get_restart_info()[source]
get_sphere_radius(offset=4)[source]
get_subgeom(indices, coord_type='cart', sort=False)[source]

Return a Geometry containing a subset of the current Geometry.

Parameters:
  • indices (iterable of ints) -- Atomic indices that the define the subset of the current Geometry.

  • coord_type (str, ("cart", "redund"), optional) -- Coordinate system of the new Geometry.

Returns:

sub_geom -- Subset of the current Geometry.

Return type:

Geometry

get_subgeom_without(indices, **kwargs)[source]
get_temporary_coords(coords)[source]
get_thermoanalysis(energy=None, cart_hessian=None, T=298.15, p=101325, point_group='c1')[source]
get_trans_rot_projector(full=False)[source]
property gradient

Negative of the force.

Returns:

gradient -- 1d array containing the negative of the current forces.

Return type:

np.array

property hessian

Matrix of second derivatives of the energy in respect to atomic displacements.

Returns:

hessian -- 2d array containing the second derivatives of the energy with respect to atomic/coordinate displacements depending on the type of coordiante system.

Return type:

np.array

property inertia_tensor
property is_analytical_2d
jmol(atoms=None, cart_coords=None)[source]

Show geometry in jmol.

property layers
mass_weigh_hessian(hessian)[source]
property masses
property masses_rep
property mm_inv

Inverted mass matrix.

Returns a diagonal matrix containing the inverted atomic masses.

property mm_sqrt_inv

Inverted square root of the mass matrix.

modes3d()[source]
property moving_atoms
moving_atoms_jmol()[source]
property mw_coords

Mass-weighted coordinates.

Returns:

mw_coords -- 1d array containing the mass-weighted cartesian coordiantes.

Return type:

np.array

property mw_gradient

Mass-weighted gradient.

Returns:

mw_gradient -- Returns the mass-weighted gradient.

Return type:

np.array

property mw_hessian

Mass-weighted hessian.

Returns:

mw_hessian -- 2d array containing the mass-weighted hessian M^(-1/2) H M^(-1/2).

Return type:

np.array

principal_axes_are_aligned()[source]

Check if the principal axes are aligned with the cartesian axes.

Returns:

aligned -- Wether the principal axes are aligned or not.

Return type:

bool

reparametrize()[source]
reset_coords(new_typed_prims=None)[source]
rmsd(geom)[source]
rotate(copy=False, rng=None)[source]
set_calculator(calculator, clear=True)[source]

Reset the object and set a calculator.

set_coord(ind, coord)[source]

Set a coordinate by index.

Parameters:
  • ind (int) -- Index in of the coordinate to set in the self.coords array.

  • coord (float) -- Coordinate value.

set_coords(coords, cartesian=False, update_constraints=False)[source]
set_h5_hessian(fn)[source]
set_restart_info(restart_info)[source]
set_results(results)[source]

Save the results from a dictionary.

Parameters:

results (dict) -- The keys in this dict will be set as attributes in the current object, with the corresponding item as value.

standard_orientation()[source]
property sum_formula
tmp_xyz_handle(atoms=None, cart_coords=None)[source]
property total_mass
unweight_mw_hessian(mw_hessian)[source]

Unweight a mass-weighted hessian.

Parameters:

mw_hessian (np.array) -- Mass-weighted hessian to be unweighted.

Returns:

hessian -- 2d array containing the hessian.

Return type:

np.array

property vdw_radii
vdw_volume(**kwargs)[source]
without_hydrogens()[source]
zero_frozen_forces(cart_forces)[source]
pysisyphus.Geometry.get_trans_rot_projector(cart_coords, masses, full=False)[source]
pysisyphus.Geometry.get_trans_rot_vectors(cart_coords, masses, rot_thresh=1e-06)[source]

Vectors describing translation and rotation.

These vectors are used for the Eckart projection by constructing a projector from them.

See Martin J. Field - A Pratcial Introduction to the simulation of Molecular Systems, 2007, Cambridge University Press, Eq. (8.23), (8.24) and (8.26) for the actual projection.

See also https://chemistry.stackexchange.com/a/74923.

Parameters:
  • cart_coords (np.array, 1d, shape (3 * atoms.size, )) -- Atomic masses in amu.

  • masses (iterable, 1d, shape (atoms.size, )) -- Atomic masses in amu.

Returns:

ortho_vecs -- 2d array containing row vectors describing translations and rotations.

Return type:

np.array(6, 3*atoms.size)

pysisyphus.Geometry.inertia_tensor(coords3d, masses)[source]

Inertita tensor.

x² xy xz |
(x y z)^T . (x y z) = | xy y² yz |
xz yz z² |

Basic infrastructure to work with internal coordinates is provided by RedundantCoords.

class pysisyphus.intcoords.RedundantCoords.HybridRedundantCoords(*args, **kwargs)[source]
class pysisyphus.intcoords.RedundantCoords.RedundantCoords(atoms, coords3d, masses=None, bond_factor=1.3, typed_prims=None, define_prims=None, constrain_prims=None, freeze_atoms=None, freeze_atoms_exclude=False, internals_with_frozen=False, define_for=None, bonds_only=False, check_bends=True, rebuild=True, bend_min_deg=15, dihed_max_deg=175, lb_min_deg=175, weighted=False, min_weight=0.3, svd_inv_thresh=0.000316, recalc_B=False, tric=False, hybrid=False, hbond_angles=False, rm_for_frag=None)[source]
property B

Wilson B-Matrix

property B_inv

Generalized inverse of the Wilson B-Matrix.

property B_inv_prim

Generalized inverse of the primitive Wilson B-Matrix.

property B_prim

Wilson B-Matrix

property Bt_inv

Transposed generalized inverse of the Wilson B-Matrix.

property Bt_inv_prim

Transposed generalized inverse of the primitive Wilson B-Matrix.

property C

Diagonal matrix. Entries for constraints are set to one.

property P

Projection matrix onto B. See [1] Eq. (4).

backtransform_hessian(redund_hessian, int_gradient=None)[source]

Transform Hessian in internal coordinates to Cartesians.

property bend_atom_indices
property bend_indices
property bond_atom_indices
property bond_indices
property bond_typed_prims
property cartesian_indices
clear()[source]
property constrained_indices
property coords
property coords3d
property dihedral_atom_indices
property dihedral_indices
eval(coords3d, attr=None)[source]
get_K_matrix(int_gradient=None)[source]
get_index_of_typed_prim(typed_prim)[source]

Index in self.typed_prims for the supplied typed_prim.

get_prim_internals_by_indices(indices)[source]
inv_B(B)[source]
inv_Bt(B)[source]
property linear_bend_indices
log(message)[source]
log_int_grad_msg(int_gradient)[source]
property outofplane_indices
property prim_coords
property prim_indices_set
property prim_internals
property primitives
print_typed_prims()[source]
project_hessian(H, shift=1000)[source]

Expects a hessian in internal coordinates. See Eq. (11) in [1].

project_vector(vector)[source]

Project supplied vector onto range of B.

return_inds(slice_)[source]
property rotation_indices
set_inds_from_typed_prims(typed_prims)[source]
set_primitive_indices(atoms, coords3d)[source]
transform_forces(cart_forces)[source]

Combination of Eq. (9) and (11) in [1].

transform_hessian(cart_hessian, int_gradient=None)[source]

Transform Cartesian Hessian to internal coordinates.

transform_int_step(int_step, update_constraints=False, pure=False)[source]
property translation_indices
property typed_prims
class pysisyphus.intcoords.RedundantCoords.TMTRIC(atoms, *args, **kwargs)[source]
class pysisyphus.intcoords.RedundantCoords.TRIC(*args, **kwargs)[source]