Skip to article frontmatterSkip to article content
Site not loading correctly?

This may be due to an incorrect BASE_URL configuration. See the MyST Documentation for reference.

Potential energy surfaces

VeloxChem enables the exploration of potential energy surfaces through efficient geometry optimizations and transition‑state searches, using the geomeTRIC module as a robust engine that provides stable structure updates for molecular systems.

Ground state optimization

VeloxChem evaluates analytic ground‑state gradients and employs quasi‑Newton optimization algorithms that approximate the local Hessian to achieve rapid convergence toward a minimum.

Python script

import veloxchem as vlx

molecule = vlx.Molecule.read_name("bithiophene")
molecule.set_dihedral_in_degrees([3, 4, 6, 10], 160)
basis = vlx.MolecularBasis.read(molecule, "def2-svp")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "b3lyp"
scf_drv.dispersion = True
results = scf_drv.compute(molecule, basis)

opt_drv = vlx.OptimizationDriver(scf_drv)
opt_results = opt_drv.compute(molecule, basis, results)
Loading...
<Figure size 650x400 with 1 Axes>

Text file

Please refer to the keyword list for a complete set of options. Dispersion can be activated in the @method settings section by using the keyword dispersion.

@jobs
task: optimize
@end

@method settings
xcfun: b3lyp
basis: def2-svp
dispersion: yes # use dft-d4 correction
@end

@molecule
charge: 0
multiplicity: 1
xyz:
...
@end

Constrained optimization

Internal coordinates (distances, angles, dihedrals) can be constrained during the molecular structure optimization with use of either the set, freeze, or scan directives.

Set or freeze internal coordinate

The set directive will aim at converging an internal coordinate to a desired value while freeze will keep it at its initial value.

Python script

import veloxchem as vlx

molecule = vlx.Molecule.read_xyz_string("""4
Hydrogen peroxide
O  -0.65564532 -0.06106286 -0.03621403
O   0.65564532  0.06106286 -0.03621403
H  -0.97628735  0.65082652  0.57474201
H   0.97628735 -0.65082652  0.57474201
""")

basis = vlx.MolecularBasis.read(molecule, "def2-svp")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "b3lyp"
scf_results = scf_drv.compute(molecule, basis)

opt_drv = vlx.OptimizationDriver(scf_drv)
opt_drv.constraints = ["set dihedral 3 1 2 4 90.0", "freeze distance 1 2"]
opt_results = opt_drv.compute(molecule, basis, results)
Loading...

Text file

@jobs
task: optimize
@end

@method settings
xcfun: b3lyp
basis: def2-svp
@end

@optimize
set dihedral 3 1 2 4 90.0
freeze distance 1 2 
@end

@molecule
charge: 0
multiplicity: 1
xyz:
...
@end

Scan coordinate

The scan directive will perform a relaxed scan of an internal coordinate from an initial to a final value in a given number of steps.

scan distance 6 1 1.4 1.5 9
scan angle    6 1 2 100 110 9
scan dihedral 6 1 2 3 0 360 19

Python script

import veloxchem as vlx

molecule = vlx.Molecule.read_xyz_string("""4
Hydrogen peroxide
O  -0.65564532 -0.06106286 -0.03621403
O   0.65564532  0.06106286 -0.03621403
H  -0.97628735  0.65082652  0.57474201
H   0.97628735 -0.65082652  0.57474201
""")
basis = vlx.MolecularBasis.read(molecule, "def2-svp")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "b3lyp"
scf_results = scf_drv.compute(molecule, basis)

opt_drv = vlx.OptimizationDriver(scf_drv)
opt_drv.constraints = ["scan dihedral 3 1 2 4 0 180 7"]
opt_results = opt_drv.compute(molecule, basis, scf_results)
<Figure size 650x400 with 1 Axes>

Text file

@jobs
task: optimize
@end

@method settings
xcfun: b3lyp
basis: def2-svp
@end

@optimize
scan dihedral 3 1 2 4 180 0 7
@end

@molecule
charge: 0
multiplicity: 1
xyz:
...
@end

Transition state optimization

Transition‑state optimization targets first‑order saddle points on the Born–Oppenheimer potential energy surface, characterized by vanishing gradients and a single imaginary Hessian eigenvalue. The resulting transition‑state geometry defines the critical point connecting reactant and product minima and provides the reference structure for calculating reaction pathways, activation energies, and vibrational signatures.

The optimization of a transition-state structure requires an initial guess of the structure that is close to the actual one.

Obtaining a start guess of the transition state structure

The transition-state guesser can be used to efficiently calculate an initial guess for a transition state structure. More information can be found in the eChem book

Python script

import veloxchem as vlx

reac1 = vlx.Molecule.read_smiles("C=CO")
reac2 = vlx.Molecule.read_smiles("C=O")

prod = vlx.Molecule.read_smiles("C(=O)CCO")

ts_guesser = vlx.TransitionStateGuesser()
ts_guesser.do_qm_scan = True
ts_guess_results = ts_guesser.find_transition_state([reac1, reac2], prod)
Loading...

Optimizing the transition state structure

Python script

import veloxchem as vlx

ts_guess=vlx.Molecule.read_xyz_string(ts_guess_results['max_qm_xyz'])
basis = vlx.MolecularBasis.read(ts_guess, "def2-svp")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "B3LYP"
scf_results = scf_drv.compute(ts_guess, basis)

opt_drv = vlx.OptimizationDriver(scf_drv)
opt_drv.transition = True
opt_results = opt_drv.compute(ts_guess, basis, scf_results)
Loading...

Text file

@jobs
task: optimize
@end

@method settings
xcfun: b3lyp
basis: def2-svp
dispersion: yes
@end

@optimize
transition: yes
@end

@molecule
charge: -1
multiplicity: 1
xyz:
...
@end

Verifying the transition state structure

Vibrational frequencies

Python script

ts_opt = vlx.Molecule.read_xyz_string(opt_results['final_geometry'])

vib_drv = vlx.VibrationalAnalysis(scf_drv)
vib_drv.do_ir = False

vib_results = vib_drv.compute(ts_opt, basis)
Loading...

Internal reaction coordinate

An intrinsic reaction coordinate (IRC) calculation can be performed to verify that a transition state structure connects the expected reactant and product structures.

Python script

import veloxchem as vlx

basis = vlx.MolecularBasis.read(ts_opt, "def2-svp")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "B3LYP"
scf_results = scf_drv.compute(ts_opt, basis)

opt_drv = vlx.OptimizationDriver(scf_drv)
opt_drv.irc = True
opt_results = opt_drv.compute(ts_opt, basis, scf_results)
opt_drv.plot_irc(opt_results)
<Figure size 650x400 with 1 Axes>

Text file

@jobs
task: optimize
@end

@method settings
xcfun: b3lyp
basis: def2-svp
dispersion: yes
@end

@optimize
irc: yes
@end

@molecule
charge: -1
multiplicity: 1
xyz:
...
@end

Excited state optimization

Python script

import veloxchem as vlx

molecule = vlx.Molecule.read_name("bithiophene")

basis = vlx.MolecularBasis.read(molecule, "def2-svp")

scf_drv = vlx.ScfRestrictedDriver()
scf_drv.xcfun = "b3lyp"
scf_results = scf_drv.compute(molecule, basis)

rsp_drv = vlx.LinearResponseEigenSolver()
rsp_drv.nstates = 2
rsp_results = rsp_drv.compute(molecule, basis, scf_results)

grad_drv = vlx.TddftGradientDriver(scf_drv)
grad_drv.state_deriv_index = 1

opt_drv = vlx.OptimizationDriver(grad_drv)
opt_results = opt_drv.compute(molecule, basis, scf_drv, rsp_drv, rsp_results)
Loading...

Text file

To optimize an excited state, you need to use the task optimize and to specify which state you want to optimize with the state_deriv_index: keyword in the @gradient section.

@jobs
task: optimize
@end

@response
property: absorption
nstates: 2
@end

@gradient
state_deriv_index: 1
@end

@method settings
xcfun: b3lyp
basis: def2-svp
@end

@molecule
charge: 0
multiplicity: 1
xyz:
...
@end
Figure: Optimization of the molecular structure in the excited state, S_1.

Figure: Optimization of the molecular structure in the excited state, S1_1.

Model chemistry
Environment
Properties
UV/vis absorption/emission