Optimization with CIOpt

CIOpt is TeraChem's dedicated optimizer for locating conical intersections and minimal-energy crossing points between electronic states. It is selected with

run ciopt

and the optimization algorithm is chosen via the cioptalgorithm keyword. CIOpt operates on excited-state methods that provide non-adiabatic couplings: CIS, TDA-TDDFT, CASSCF, CASCI (including FOMO-CASCI), CIS-NO, STEX, REKS, and (h)h-TDA.

CIOpt can locate:

Two-state conical intersections — both minimum-energy crossings (MECI, the lowest-energy point on the seam) and minimum-distance crossings (MDCI, the seam point closest to a user-supplied reference geometry).
Three-state conical intersections — ME3CI and MD3CI variants, where the optimizer locates a point at which three electronic states are simultaneously degenerate.

A separate, older conical-intersection optimizer is available through DL-Find with run conical and min_method_ci; that path is documented under DL-Find. For new work the run ciopt path documented here is preferred.

Choosing an algorithm

cioptalgorithm selects both the kind of crossing being located (MECI vs MDCI, two-state vs three-state) and the formulation used to enforce the degeneracy constraint:

Algorithm value	Crossing type	Formulation
`meci_brs`	2-state MECI	Branching-space updating (default)
`me3ci_brs`	3-state ME3CI	Branching-space updating
`mdci_brs`	2-state MDCI	Branching-space updating
`md3ci_brs`	3-state MD3CI	Branching-space updating
`meci_penalty`	2-state MECI	Penalty function
`me3ci_penalty`	3-state ME3CI	Penalty function
`mdci_penalty`	2-state MDCI	Penalty function
`md3ci_penalty`	3-state MD3CI	Penalty function
`meci_aug`	2-state MECI	Augmented Lagrangian
`me3ci_aug`	3-state ME3CI	Augmented Lagrangian
`mdci_aug`	2-state MDCI	Augmented Lagrangian
`minimize`	Single-state minimum	Standard minimization through CIOpt

The three formulations differ in how the energy-gap constraint is imposed:

Branching-space (BRS) methods project the gradient into the branching plane using the energy-difference and derivative-coupling vectors, following the weighting scheme of Keal, Koslowski, and Thiel.¹ The weights are controlled by cioptc3 and cioptc4.
Penalty-function methods minimize an objective in which the energy gap appears through a smoothed penalty term. The penalty strength is controlled by cioptsigma1/cioptsigma2 and the smoothing parameter by cioptalpha.
Augmented-Lagrangian methods iteratively tighten a Lagrange multiplier on the energy gap; ciopttighten and ciopttightener control the tightening schedule.

Specifying the states

The two (or three) electronic states involved in the crossing are specified with the nacstate1, nacstate2, and nacstate3 keywords (state indices are 0-based, so nacstate1 0 is the ground state). For methods that support multiple spin multiplicities, the multiplicity of each state is set with nacstatemult1, nacstatemult2, and nacstatemult3. The non-adiabatic coupling between the specified states must be available from the underlying electronic structure method, which typically requires turning on the appropriate gradient and CPHF/CPSA options (e.g. cphfiter, cphftol).

For an MDCI search, the reference geometry from which "distance" is measured defaults to the input coordinates; it can be overridden with cioptrefgeom.

Example: MECI with the penalty function and CASSCF

This input locates the S0/S1 MECI of ethylene at the CASSCF(2,2)/6-31G* level using the penalty-function formulation:

meci_penalty_casscf.in

basis             6-31gs
coordinates       ethylene_brs.xyz
method            rhf
run               ciopt

cioptalgorithm    meci_penalty
cioptthresh       0.0005

charge            0
spinmult          1

casscf            yes
closed            7
active            2
cassinglets       2

nacstate1         0
nacstate2         1
nacstatemult1     1

precision         double
convthre          1.0e-7
threall           1.0e-14

end

Example: MECI with branching-space updating and CASSCF

The same MECI search using the BRS formulation:

meci_brs_casscf.in

basis             6-31gs
coordinates       ethylene_brs.xyz
method            rhf
run               ciopt

cioptalgorithm    meci_brs
cioptthresh       0.001

casscf            yes
closed            7
active            2
cassinglets       2

nacstate1         0
nacstate2         1
nacstatemult1     1

tvarnac           yes

end

Example: MDCI with FOMO-CASCI

This input locates the seam point closest to the input geometry, using FOMO-CASCI as the underlying method and BRS for the constraint:

mdci_brs_fomo.in

basis             6-31gs
coordinates       geom0.xyz
method            rhf
run               ciopt

cioptalgorithm    mdci_brs
cioptthresh       0.001
cioptc4           0.9

fon               yes
fon_temperature   0.15
purify            no

casci             yes
closed            7
active            2
cassinglets       2

nacstate1         0
nacstate2         1
nacstatemult1     1
nacstatemult2     1

end

Summary of CIOpt keywords

Algorithm and convergence

Keyword	Type	Default	Description
`cioptalgorithm`	string	`meci_brs`	Algorithm and crossing type, as listed in the table above
`cioptthresh`	float	0.001	Overall convergence threshold
`cioptgradtolinf`	float	4.5e-4	Infinity-norm gradient tolerance
`cioptgradtol2`	float	4.5e-4	2-norm gradient tolerance
`cioptgradtolperp`	float	5.0e-3	Tolerance for the gradient component perpendicular to the branching plane
`cioptsteptolinf`	float	1.8e-3	Infinity-norm step tolerance
`cioptsteptol2`	float	1.2e-3	2-norm step tolerance
`cioptfxtol`	float	1.0e-6	Tolerance on objective function decrease

State selection

Keyword	Type	Default	Description
`nacstate1`	int	0	Index (0-based) of the first state
`nacstate2`	int	1	Index (0-based) of the second state
`nacstate3`	int	2	Index of the third state (3-state algorithms only)
`nacstatemult1`	int	1	Spin multiplicity of state 1
`nacstatemult2`	int	1	Spin multiplicity of state 2
`nacstatemult3`	int	1	Spin multiplicity of state 3
`cioptrefgeom`	filename	input coords	Reference geometry for MDCI searches

Penalty function (`*_penalty` algorithms)

Keyword	Type	Default	Description
`cioptalpha`	float	0.02	Smoothing parameter \(\alpha\) in the penalty objective. A value of 0 gives a non-smooth objective.
`cioptsigma1`	float	4.0	Penalty strength for the first energy-gap constraint
`cioptsigma2`	float	4.0	Penalty strength for the second constraint (three-state algorithms)
`cioptsigmamult`	float	4.0	Multiplier applied to all \(\sigma\) values

Branching-space updating (`*_brs` algorithms)

Keyword	Type	Default	Description
`cioptc3`	float	1.0	Weighting on the energy-difference gradient component
`cioptc4`	float	0.9	Weighting on the derivative-coupling gradient component
`cioptbpdim`	int	0	Dimension of the branching space (used with `cioptuseweighting`)
`cioptuseweighting`	int	0	Enable weighting for projected-gradient methods
`cioptweightingstart`	float	0.002	Energy gap at which weighting begins
`cioptweightingthresh`	float	0.4	Threshold controlling weighting behavior

Augmented Lagrangian (`*_aug` algorithms)

Keyword	Type	Default	Description
`ciopttighten`	int	0	If nonzero, tighten the multiplier each macroiteration
`ciopttightener`	float	0.66	Tightening factor applied each iteration

Optimizer and line search

Keyword	Type	Default	Description
`cioptoptimizer`	string	`bfgs`	Underlying optimizer: `bfgs` or `nelder_mead`
`cioptlinesearch`	string	`linesearch`	Step control: `linesearch` or `trustregion`
`cioptlbfgsm`	int	20	L-BFGS history length. A value of 0 reduces to gradient descent.
`cioptftol`	float	1.0e-4	Wolfe condition tolerance on function decrease
`cioptgtol`	float	0.9	Wolfe condition tolerance on gradient

Mass weighting and units

Keyword	Type	Default	Description
`cioptunits`	string	`angstrom`	Units for input/output coordinates: `angstrom` or `bohr`
`cioptnomw`	int	0	If nonzero, disable mass weighting
`cioptalphamw`	float	0.1	Mass-weighting parameter
`cioptdmwpre`	float	2.0	Mass-weighting prefactor

Output and diagnostics

Keyword	Type	Default	Description
`cioptgh`	int	0	If nonzero, return the gradient-difference (g) and derivative-coupling (h) vectors
`cioptflipc`	int	0	If nonzero, flip the sign of the derivative-coupling vector
`ciopttheta`	float	0.0	Rotation angle in the g/h plane (diagnostic)
`cioptprc`	int	0	Plane-rotation check flag
`cioptprt`	float	0.4	Plane-rotation threshold
`cioptrebuild`	int	0	Force rebuild of the optimizer state

T. W. Keal, A. Koslowski, and W. Thiel, Theor. Chem. Acc. 118, 837 (2007). ↩