Skip to content

Geometry Optimization

Geometry optimizations and transition-state searches are configured using the OptimizationInput building block. The initial geometry is specified through the initial_molecule argument, and the underlying energy method through method.

Running the input through sierra.run() produces an OptimizationResult. The OptimizationResult most notably contains the fields converged, final_molecule and energies. These results can easily be utilized for more complex workflows.

Simple geometry optimization

In this example a water molecule geometry is optimized using the GFN1-xTB method:

import sierra
from sierra.inputs import *

water = Molecule(pubchem="water")
method = XTBMethod(model="gfn1")
input = OptimizationInput(initial_molecule=water, method=method)
result = sierra.run(input)

print(f"Converged: {result.converged}")
#> Converged: True

print("Optimized molecule:")
#> Optimized molecule:
print(result.final_molecule.write())
"""
3
0 1
O     0.013698339491   0.010141575412  -0.007165048852
H     0.268966472642   0.893269716708   0.262654097613
H     0.601535188188  -0.248811292744  -0.717989047365
"""

print("Energies at the each optimization step:")
#> Energies at the each optimization step:
for i, e in enumerate(result.energies):
    print(f"  {i+1} {e:15.9f}")
    #>   1    -5.768440813
    #>   2    -5.768678192
    #>   3    -5.768772267
    #>   4    -5.768781558
    #>   5    -5.768783910
    #>   6    -5.768783913

Constrained geometry optimization

Constraints on bond-lengths, angles, and dihedrals can be applied during the geometry optimization. The CoordinateConstraint building block can either:

  • Maintain the original value of the coordinate with freeze=True
  • Drive the coordinate towards a value specified by the value argument
import sierra
from sierra.inputs import *
from codex.models.components import CoordinateConstraint

water = Molecule(pubchem="water")
method = XTBMethod(model="gfn1")
constraints = [
    CoordinateConstraint(indices=[0, 1], value="2.0 bohr"),
    CoordinateConstraint(indices=[0, 1, 2], freeze=True),
]
input = OptimizationInput(
    initial_molecule=water, method=method, constraints=constraints
)
result = sierra.run(input)

final_molecule = result.final_molecule

print("Constrained coordinates:")
#> Constrained coordinates:
for c in constraints:
    print(f"Coordinate {c.indices}: {final_molecule.measure(c.indices)}")
    #> Coordinate [0, 1]: 1.9999999971231277
    #> Coordinate [0, 1, 2]: 38.01000739722804

print(f"Unconstrained OH bond length: {final_molecule.measure([0, 2])} Bohr")
#> Unconstrained OH bond length: 1.8308089792371645 Bohr

OptimizationInput

Configuration for performing an optimization.

Fields

initial_molecule

Specify the initial molecule at the start of the optimization.

method

The method for the geometry optimization.

constraints

Specify the coordinate constraints to place on the optimization.

transition_search

If True run a transition state search rather than a energy minimization. The transition state optimisation follows the most negative eigenvector of the Hessian uphill, and the remaining eigenvectors downhill. Starting from a geometry close to the target transition state is crucial for finding the correct transition state.

  • Type: bool
  • Default: False
details

Specify additional details to supply to the optimizer. A list of supported options is given at the bottom of this page.

  • Type: Dict[str, Any]

OptimizationResult

The results of an optimization calculation.

Fields

All of the fields in OptimizationInput and the following:

converged

If the optimization converged or not

  • Type: bool
energies

A list of energies at each step of the geometry optimization

  • Type: Array[float]
extras

Additional key/value pairs generated during the optimization

  • Type: Mapping[str, Any]
  • Default:{}
final_molecule

The optimized molecule or the result of the transition state search

trajectory

A list of results for each step in the geometry optimization

Available options for input details

max_iter

Maximum number of optimization steps.

  • The type is int
  • The default is 300
  • The value must be positive
energy_change

Convergence threshold on the energy (in a.u.) for geometry optimization.

  • The type is float
  • The default is 1e-6
  • The value must be positive
gradient_max

Convergence threshold on the max element of gradient (in a.u.) for geometry optimization.

  • The type is float
  • The default is 3e-4 (9e-3 when OrbNet Denali is used)
  • The value must be positive
gradient_rms

Convergence threshold on the root-mean-square of gradient (in a.u.) for geometry optimization.

  • The type is float
  • The default is 4.5e-4 (1.35e-2 when OrbNet Denali is used)
  • The value must be positive
step_max

Convergence threshold on the max element of the step (in a.u.) for geometry optimization.

  • The type is float
  • The default is 1.2e-3
  • The value must be positive
step_rms

Convergence threshold on the root-mean-square of the step (in a.u.) for geometry optimization.

  • The type is float
  • The default is 1.8e-3
  • The value must be positive
trust_radius_min

Minimum trust radius.

  • The type is LengthQuantity
  • The default is 1e-6 bohr (1e-3 bohr when OrbNet Denali is used)
  • The value must be positive