8. mcscf — Multi-configurational self-consistent field

CASCI and CASSCF

Simple usage:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol).run()
>>> mc = mcscf.CASCI(mf, 6, 6)
>>> mc.kernel()[0]
-108.980200816243354
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.kernel()[0]
-109.044401882238134
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> cas_list = [5,6,8,9] # pick orbitals for CAS space, 1-based indices
>>> mo = mcscf.sort_mo(mc, mf.mo_coeff, cas_list)
>>> mc.kernel(mo)[0]
-109.007378939813691

mcscf.CASSCF() or mcscf.CASCI() returns a proper instance of CASSCF/CASCI class. There are some parameters to control the CASSCF/CASCI method.

verbose
: int
Print level. Default value equals to Mole.verbose.
max_memory
: float or int
Allowed memory in MB. Default value equals to Mole.max_memory.
ncas
: int
Active space size.
nelecas
: tuple of int
Active (nelec_alpha, nelec_beta)
ncore
: int or tuple of int
Core electron number. In UHF-CASSCF, it’s a tuple to indicate the different core eletron numbers.
natorb
: bool
Whether to restore the natural orbital during CASSCF optimization. Default is not.
canonicalization
: bool
Whether to canonicalize orbitals. Default is True.
fcisolver
: an instance of FCISolver

The pyscf.fci module provides several FCISolver for different scenario. Generally, fci.direct_spin1.FCISolver can be used for all RHF-CASSCF. However, a proper FCISolver can provide better performance and better numerical stability. One can either use fci.solver() function to pick the FCISolver by the program or manually assigen the FCISolver to this attribute, e.g.

>>> from pyscf import fci
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> mc.fcisolver = fci.solver(mol, singlet=True)
>>> mc.fcisolver = fci.direct_spin1.FCISolver(mol)

You can control FCISolver by setting e.g.:

>>> mc.fcisolver.max_cycle = 30
>>> mc.fcisolver.conv_tol = 1e-7

For more details of the parameter for FCISolver, See fci.

By replacing this fcisolver, you can easily use the CASCI/CASSCF solver with other FCI replacements, such as DMRG, QMC. See dmrgscf and fciqmcscf.

The Following attributes are used for CASSCF

conv_tol
: float
Converge threshold. Default is 1e-7
conv_tol_grad
: float
Converge threshold for CI gradients and orbital rotation gradients. Default is 1e-4
max_stepsize
: float

The step size for orbital rotation. Small step size is prefered. Default is 0.03. (NOTE although the default step size is small enough for many systems, it happens that the orbital optimizor crosses the barriar of local minimum and converge to the neighbour solution, e.g. the CAS(4,4) for C2H4 in the test files. In these systems, adjusting max_stepsize, max_ci_stepsize and max_cycle_micro, max_cycle_micro_inner and ah_start_tol may be helpful)

>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.max_stepsize = .01
>>> mc.max_cycle_micro = 1
>>> mc.max_cycle_macro = 100
>>> mc.max_cycle_micro_inner = 1
>>> mc.ah_start_tol = 1e-6
max_ci_stepsize
: float
The max size for approximate CI updates. The approximate updates are used in 1-step algorithm, to estimate the change of CI wavefunction wrt the orbital rotation. Small step size is prefered. Default is 0.01.
max_cycle_macro
: int
Max number of macro iterations. Default is 50.
max_cycle_micro
: int
Max number of micro iterations in each macro iteration. Depending on systems, increasing this value might reduce the total macro iterations. Generally, 2 - 3 steps should be enough. Default is 2.
max_cycle_micro_inner
: int
Max number of steps for the orbital rotations allowed for the augmented hessian solver. It can affect the actual size of orbital rotation. Even with a small max_stepsize, a few max_cycle_micro_inner can accumulate the rotation and leads to a significant change of the CAS space. Depending on systems, increasing this value migh reduce the total number of macro iterations. The value between 2 - 8 is preferred. Default is 4.
frozen
: int or list
If integer is given, the inner-most orbitals are excluded from optimization. Given the orbital indices (0-based) in a list, any doubly occupied core orbitals, active orbitals and external orbitals can be frozen.
ah_level_shift
: float, for AH solver.
Level shift for the Davidson diagonalization in AH solver. Default is 0.
ah_conv_tol
: float, for AH solver.
converge threshold for Davidson diagonalization in AH solver. Default is 1e-8.
ah_max_cycle
: float, for AH solver.
Max number of iterations allowd in AH solver. Default is 20.
ah_lindep
: float, for AH solver.
Linear dependence threshold for AH solver. Default is 1e-16.
ah_start_tol
: flat, for AH solver.
In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 1e-4.
ah_start_cycle
: int, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 3.

ah_conv_tol, ah_max_cycle, ah_lindep, ah_start_tol and ah_start_cycle can affect the accuracy and performance of CASSCF solver. Lower ah_conv_tol and ah_lindep can improve the accuracy of CASSCF optimization, but slow down the performance.

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.conv_tol = 1e-10
>>> mc.ah_conv_tol = 1e-5
>>> mc.kernel()
-109.044401898486001
>>> mc.ah_conv_tol = 1e-10
>>> mc.kernel()
-109.044401887945668
chkfile
: str
Checkpoint file to save the intermediate orbitals during the CASSCF optimization. Default is the checkpoint file of mean field object.

Saved results

e_tot
: float
Total MCSCF energy (electronic energy plus nuclear repulsion)
ci
: ndarray
CAS space FCI coefficients
converged
: bool, for CASSCF only
It indicates CASSCF optimization converged or not.
mo_energy: ndarray,
Diagonal elements of general Fock matrix
mo_coeff
: ndarray, for CASSCF only
Optimized CASSCF orbitals coefficients Note the orbitals are NOT natural orbitals by default. There are two inbuilt methods to convert the mo_coeff to natural orbitals. 1. Set .natorb attribute. It can be used before calculation. 2. call .cas_natorb_ method after the calculation to in-place convert the orbitals

8.1. CASSCF active space solver

8.1.1. DMRG solver

8.1.2. FCIQMC solver

8.1.3. State-average FCI solver

8.1.4. State-average with mixed solver

8.2. Symmetry broken

8.3. Initial guess

8.4. Canonical orbitals

There are two relevant parameters for orbital canonicalization. They are mc.canonicalization and mc.natorb (assuming the MCSCF object is mc). In the CASCI/CASSCF calculations, the resultant orbitals are stored in the attribute mc.mo_coeff. These orbitals may be identical or partially identical to the initial orbitals, depending on the values of mc.canonicalization and mc.natorb.

mc.canonicalization controls whether the resultant CASCI/CASSCF orbitals are canonicalized with respect to the general Fock matrix. General Fock matrix is defined as

\[\begin{split}\mathbf{F} &= \mathbf{h}_{core} + \mathbf{J} - \mathbf{K} \\ J_{pq} &= \sum_{rs} (pq|rs) \gamma_{sr} \\ K_{pq} &= \sum_{qr} (pq|rs) \gamma_{qr} \\\end{split}\]

\(\gamma\) is the total density matrix which is the summation of doubly occupied core density matrix and correlated density matrix in active space. If mc.canonicalization is enabled, the CASCI/CASSCF program will call the mc.canonicalize() function to diagonalize the core and external space wrt the general Fock matrix. The eigenvalues of the core and external subspace are stored in attribute mc.mo_energy. By default, mc.canonicalization is enabled because the canonicalized MCSCF orbitals can simplify the algorithm of MRPT methods.

mc.natorb controls whether the CASCI/CASSCF active space orbitals are transformed to natural orbitals of the correlated density matrix. When this parameter is enabled, the natural orbitals will be stored in the active part of the attribute mc.mo_coeff and the CI coefficients mc.ci (if applicable) will be transformed accordingly. By default mc.natorb is disabled and it is important for the MCSCF solver. Generally, the value of mc.natorb does not affect (the default) FCI solver because an independent CASCI calculation following a previous MCSCF calculation should give the same solutions no matter mc.natorb is enabled or not. But this is not true for some external large active space solvers such as DMRG, selected CI methods. The CASCI calculation may produce different answers depending on the value of mc.natorb. Therefore, it is recommended to disable mc.natorb in your calculation.

Following presents what the mc.mo_coeff would be like for different combinations of mc.canonicalization and mc.natorb in a CASCI calculation:

  • mc.canonicalization = False and mc.natorb = False:

The resultant orbitals mc.mo_coeff are identical to the input orbitals. If the CASCI was initialized with a RHF calculation, mc.mo_coeff points to RHF orbitals.

  • mc.canonicalization = True and mc.natorb = False:

Core part and external part of mc.mo_coeff are canonicalized orbitals, which diagonalize the core and external blocks of general Fock matrix. The orbitals in active space are identical to the active orbitals in the input.

  • mc.canonicalization = False and mc.natorb = True

Core and external part of mc.mo_coeff are identical to the core and external part of the input orbitals. Active space orbitals are transformed to the natural orbitals of the correlated density matrix.

  • mc.canonicalization = True and mc.natorb = True

mc.mo_coeff are completely different to the input orbitals.

There is another parameter mc.sorting_mo_energy which may affect the ordering of MCSCF orbitals when mc.canonicalization or mc.natorb is enabled. Generally, the canonical orbitals in the core and external space are sorted by the orbital energies (from low to high) and the natural orbitals in the active space are sorted by natural occupations (from large to small). This ordering may not be held if point group symmetry is enabled in the calculation. When a system has high spatial symmetry and point group symmetry is enabled, each SCF orbital will be assigned to an irreducible representation label. In the MCSCF calculation and the canonicalization, the irreducible representation label of the orbitals will not be changed. They are always the same to the symmetry labels of the input orbitals. Although the orbitals are still ordered within each irreducible representation, the orbital energies (or occupancies) for all orbitals are not strictly sorted. Setting mc.sorting_mo_energy = Trye (though not recommended) can force the orbitals to be sorted regardless whether the point group symmetry is enabled. In certain scenario, you may want to enable mc.natorb and mc.sorting_mo_energy. examples/dmrg/31-cr2_scan/cr2-scan.py provides one example that you need to enable the two parameters. In that example, the dissociation curve of Cr dimer was scanned by heat-bath selected-CI method in which the active space of selected-CI-CASSCF was gradually enlarged in a series of CASSCF calculations. Since the selected-CI algorithm depends on the initial single determinant, the orbital ordering do have matters to the final CASSCF results. Thus mc.natorb and mc.sorting_mo_energy have to be enabled to make sure that the each selected-CI starts from the similar initial reference for each point on the dissociation curve. At some critical points, the difference in the orbital ordering in the active space can lead to discontinuous potential energy curve.

8.5. Program reference

8.5.1. CASCI

class pyscf.mcscf.casci.CASCI(mf_or_mol, ncas, nelecas, ncore=None)[source]
Args:
mf_or_mol
: SCF object or Mole object
SCF or Mole to define the problem size.
ncas
: int
Number of active orbitals.
nelecas
: int or a pair of int
Number of electrons in active space.
Kwargs:
ncore
: int
Number of doubly occupied core orbitals. If not presented, this parameter can be automatically determined.
Attributes:
verbose
: int
Print level. Default value equals to Mole.verbose.
max_memory
: float or int
Allowed memory in MB. Default value equals to Mole.max_memory.
ncas
: int
Active space size.
nelecas
: tuple of int
Active (nelec_alpha, nelec_beta)
ncore
: int or tuple of int
Core electron number. In UHF-CASSCF, it’s a tuple to indicate the different core eletron numbers.
natorb
: bool
Whether to restore the natural orbital in CAS space. Default is not. Be very careful to set this parameter when CASCI/CASSCF are combined with DMRG solver because this parameter changes the orbital ordering which DMRG relies on.
canonicalization
: bool
Whether to canonicalize orbitals. Default is True.
sorting_mo_energy
: bool
Whether to sort the orbitals based on the diagonal elements of the general Fock matrix. Default is False.
fcisolver
: an instance of FCISolver

The pyscf.fci module provides several FCISolver for different scenario. Generally, fci.direct_spin1.FCISolver can be used for all RHF-CASSCF. However, a proper FCISolver can provide better performance and better numerical stability. One can either use fci.solver() function to pick the FCISolver by the program or manually assigen the FCISolver to this attribute, e.g.

>>> from pyscf import fci
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> mc.fcisolver = fci.solver(mol, singlet=True)
>>> mc.fcisolver = fci.direct_spin1.FCISolver(mol)

You can control FCISolver by setting e.g.:

>>> mc.fcisolver.max_cycle = 30
>>> mc.fcisolver.conv_tol = 1e-7

For more details of the parameter for FCISolver, See fci.

Saved results

e_tot
: float
Total MCSCF energy (electronic energy plus nuclear repulsion)
e_cas
: float
CAS space FCI energy
ci
: ndarray
CAS space FCI coefficients
mo_coeff
: ndarray
When canonicalization is specified, the orbitals are canonical orbitals which make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is specified, the active segment of the mo_coeff is natural orbitls.
mo_energy
: ndarray
Diagonal elements of general Fock matrix (in mo_coeff representation).

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASCI(mf, 6, 6)
>>> mc.kernel()[0]
-108.980200816243354
as_scanner(mc)

Generating a scanner for CASCI PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total CASCI energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters of MCSCF object are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, ...) during calculation.

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mf = scf.RHF(gto.Mole().set(verbose=0))
>>> mc_scanner = mcscf.CASCI(mf, 4, 4).as_scanner()
>>> mc_scanner(gto.M(atom='N 0 0 0; N 0 0 1.1'))
>>> mc_scanner(gto.M(atom='N 0 0 0; N 0 0 1.5'))
canonicalize(mc, mo_coeff=None, ci=None, eris=None, sort=False, cas_natorb=False, casdm1=None, verbose=3, with_meta_lowdin=True)

Canonicalized CASCI/CASSCF orbitals of effecitve Fock matrix. Effective Fock matrix is built with one-particle density matrix (see also mcscf.casci.get_fock()). For state-average CASCI/CASSCF object, the canonicalized orbitals are based on the state-average density matrix. To obtain canonicalized orbitals for an individual state, you need to pass “casdm1” of the specific state to this function.

Args:
mc: a CASSCF/CASCI object or RHF object
Kwargs:
mo_coeff (ndarray): orbitals that span the core, active and external
space.
ci (ndarray): CI coefficients (or objects to represent the CI
wavefunctions in DMRG/QMC-MCSCF calculations).
eris: Integrals for the MCSCF object. Input this object to reduce the
overhead of computing integrals. It can be generated by mc.ao2mo() method.
sort (bool): Whether the canonicalized orbitals are sorted based on
orbital energy (diagonal part of the effective Fock matrix) within each subspace (core, active, external). If the point group symmetry is not available in the system, the orbitals are always sorted. When the point group symmetry is available, sort=False will keep the symmetry label of input orbitals and only sort the orbitals in each symmetry block while sort=True will reorder all orbitals in each subspace and the symmetry labels may be changed.
cas_natorb (bool): Whether to transform the active orbitals to natual
orbitals
casdm1 (ndarray): 1-particle density matrix in active space. This
density matrix is used to build effective fock matrix. Without input casdm1, the density matrix is computed with the input ci coefficients/object. If neither ci nor casdm1 were given, density matrix is computed by mc.fcisolver.make_rdm1() method. For state-average CASCI/CASCF calculation, this results in a set of canonicalized orbitals of state-average effective Fock matrix. To canonicalize the orbitals for one particular state, you can assign the density matrix of that state to the kwarg casdm1.
Returns:
A tuple, (natural orbitals, CI coefficients, orbital energies) The orbital energies are the diagonal terms of effective Fock matrix.
canonicalize_(mo_coeff=None, ci=None, eris=None, sort=False, cas_natorb=False, casdm1=None, verbose=None, with_meta_lowdin=True)[source]

Canonicalized CASCI/CASSCF orbitals of effecitve Fock matrix. Effective Fock matrix is built with one-particle density matrix (see also mcscf.casci.get_fock()). For state-average CASCI/CASSCF object, the canonicalized orbitals are based on the state-average density matrix. To obtain canonicalized orbitals for an individual state, you need to pass “casdm1” of the specific state to this function.

Args:
mc: a CASSCF/CASCI object or RHF object
Kwargs:
mo_coeff (ndarray): orbitals that span the core, active and external
space.
ci (ndarray): CI coefficients (or objects to represent the CI
wavefunctions in DMRG/QMC-MCSCF calculations).
eris: Integrals for the MCSCF object. Input this object to reduce the
overhead of computing integrals. It can be generated by mc.ao2mo() method.
sort (bool): Whether the canonicalized orbitals are sorted based on
orbital energy (diagonal part of the effective Fock matrix) within each subspace (core, active, external). If the point group symmetry is not available in the system, the orbitals are always sorted. When the point group symmetry is available, sort=False will keep the symmetry label of input orbitals and only sort the orbitals in each symmetry block while sort=True will reorder all orbitals in each subspace and the symmetry labels may be changed.
cas_natorb (bool): Whether to transform the active orbitals to natual
orbitals
casdm1 (ndarray): 1-particle density matrix in active space. This
density matrix is used to build effective fock matrix. Without input casdm1, the density matrix is computed with the input ci coefficients/object. If neither ci nor casdm1 were given, density matrix is computed by mc.fcisolver.make_rdm1() method. For state-average CASCI/CASCF calculation, this results in a set of canonicalized orbitals of state-average effective Fock matrix. To canonicalize the orbitals for one particular state, you can assign the density matrix of that state to the kwarg casdm1.
Returns:
A tuple, (natural orbitals, CI coefficients, orbital energies) The orbital energies are the diagonal terms of effective Fock matrix.
cas_natorb(mo_coeff=None, ci=None, eris=None, sort=False, casdm1=None, verbose=None, with_meta_lowdin=True)[source]

Transform active orbitals to natrual orbitals, and update the CI wfn

Args:
mc : a CASSCF/CASCI object or RHF object
Kwargs:
sort
: bool
Sort natural orbitals wrt the occupancy.
Returns:
A tuple, the first item is natural orbitals, the second is updated CI coefficients, the third is the natural occupancy associated to the natural orbitals.
cas_natorb_(mo_coeff=None, ci=None, eris=None, sort=False, casdm1=None, verbose=None, with_meta_lowdin=True)[source]

Transform active orbitals to natrual orbitals, and update the CI wfn

Args:
mc : a CASSCF/CASCI object or RHF object
Kwargs:
sort
: bool
Sort natural orbitals wrt the occupancy.
Returns:
A tuple, the first item is natural orbitals, the second is updated CI coefficients, the third is the natural occupancy associated to the natural orbitals.
fix_spin(shift=0.2, ss=None)

Use level shift to control FCI solver spin.

\[(H + shift*S^2) |\Psi\rangle = E |\Psi\rangle\]
Kwargs:
shift
: float
Energy penalty for states which have wrong spin
ss
: number
S^2 expection value == s*(s+1)
fix_spin_(shift=0.2, ss=None)[source]

Use level shift to control FCI solver spin.

\[(H + shift*S^2) |\Psi\rangle = E |\Psi\rangle\]
Kwargs:
shift
: float
Energy penalty for states which have wrong spin
ss
: number
S^2 expection value == s*(s+1)
get_h1cas(casci, mo_coeff=None, ncas=None, ncore=None)

CAS sapce one-electron hamiltonian

Args:
casci : a CASSCF/CASCI object or RHF object
Returns:
A tuple, the first is the effective one-electron hamiltonian defined in CAS space, the second is the electronic energy from core.
get_h1eff(mo_coeff=None, ncas=None, ncore=None)[source]

CAS sapce one-electron hamiltonian

Args:
casci : a CASSCF/CASCI object or RHF object
Returns:
A tuple, the first is the effective one-electron hamiltonian defined in CAS space, the second is the electronic energy from core.
get_h2cas(mo_coeff=None)[source]

Computing active space two-particle Hamiltonian.

Note It is different to get_h2eff when df.approx_hessian is applied, in which get_h2eff function returns the DF integrals while get_h2cas returns the regular 2-electron integrals.

get_h2eff(mo_coeff=None)[source]

Computing active space two-particle Hamiltonian.

Note It is different to get_h2cas when df.approx_hessian is applied. in which get_h2eff function returns the DF integrals while get_h2cas returns the regular 2-electron integrals.

h1e_for_cas(casci, mo_coeff=None, ncas=None, ncore=None)

CAS sapce one-electron hamiltonian

Args:
casci : a CASSCF/CASCI object or RHF object
Returns:
A tuple, the first is the effective one-electron hamiltonian defined in CAS space, the second is the electronic energy from core.
kernel(mo_coeff=None, ci0=None, verbose=None)[source]
Returns:
Five elements, they are total energy, active space CI energy, the active space FCI wavefunction coefficients or DMRG wavefunction ID, the MCSCF canonical orbital coefficients, the MCSCF canonical orbital coefficients.

They are attributes of mcscf object, which can be accessed by .e_tot, .e_cas, .ci, .mo_coeff, .mo_energy

make_rdm1(mo_coeff=None, ci=None, ncas=None, nelecas=None, ncore=None, **kwargs)[source]

One-particle density matrix in AO representation

make_rdm1s(mo_coeff=None, ci=None, ncas=None, nelecas=None, ncore=None, **kwargs)[source]

One-particle density matrices for alpha and beta spin on AO basis

sort_mo(caslst, mo_coeff=None, base=1)[source]

Pick orbitals for CAS space

Args:

casscf : an CASSCF or CASCI object

mo_coeff
: ndarray or a list of ndarray
Orbitals for CASSCF initial guess. In the UHF-CASSCF, it’s a list of two orbitals, for alpha and beta spin.
caslst
: list of int or nested list of int
A list of orbital indices to represent the CAS space. In the UHF-CASSCF, it’s consist of two lists, for alpha and beta spin.
Kwargs:
base
: int
0-based (C-style) or 1-based (Fortran-style) caslst
Returns:
An reoreded mo_coeff, which put the orbitals given by caslst in the CAS space

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> cas_list = [5,6,8,9] # pi orbitals
>>> mo = mc.sort_mo(cas_list)
>>> mc.kernel(mo)[0]
-109.007378939813691
state_average_(weights=(0.5, 0.5))[source]

State average over the energy. The energy funcitonal is E = w1<psi1|H|psi1> + w2<psi2|H|psi2> + ...

Note we may need change the FCI solver to

mc.fcisolver = fci.solver(mol, False)

before calling state_average_(mc), to mix the singlet and triplet states

state_specific_(state=1)[source]

For excited state

Kwargs:
state : int 0 for ground state; 1 for first excited state.
pyscf.mcscf.casci.as_scanner(mc)[source]

Generating a scanner for CASCI PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total CASCI energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters of MCSCF object are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, ...) during calculation.

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mf = scf.RHF(gto.Mole().set(verbose=0))
>>> mc_scanner = mcscf.CASCI(mf, 4, 4).as_scanner()
>>> mc_scanner(gto.M(atom='N 0 0 0; N 0 0 1.1'))
>>> mc_scanner(gto.M(atom='N 0 0 0; N 0 0 1.5'))
pyscf.mcscf.casci.canonicalize(mc, mo_coeff=None, ci=None, eris=None, sort=False, cas_natorb=False, casdm1=None, verbose=3, with_meta_lowdin=True)[source]

Canonicalized CASCI/CASSCF orbitals of effecitve Fock matrix. Effective Fock matrix is built with one-particle density matrix (see also mcscf.casci.get_fock()). For state-average CASCI/CASSCF object, the canonicalized orbitals are based on the state-average density matrix. To obtain canonicalized orbitals for an individual state, you need to pass “casdm1” of the specific state to this function.

Args:
mc: a CASSCF/CASCI object or RHF object
Kwargs:
mo_coeff (ndarray): orbitals that span the core, active and external
space.
ci (ndarray): CI coefficients (or objects to represent the CI
wavefunctions in DMRG/QMC-MCSCF calculations).
eris: Integrals for the MCSCF object. Input this object to reduce the
overhead of computing integrals. It can be generated by mc.ao2mo() method.
sort (bool): Whether the canonicalized orbitals are sorted based on
orbital energy (diagonal part of the effective Fock matrix) within each subspace (core, active, external). If the point group symmetry is not available in the system, the orbitals are always sorted. When the point group symmetry is available, sort=False will keep the symmetry label of input orbitals and only sort the orbitals in each symmetry block while sort=True will reorder all orbitals in each subspace and the symmetry labels may be changed.
cas_natorb (bool): Whether to transform the active orbitals to natual
orbitals
casdm1 (ndarray): 1-particle density matrix in active space. This
density matrix is used to build effective fock matrix. Without input casdm1, the density matrix is computed with the input ci coefficients/object. If neither ci nor casdm1 were given, density matrix is computed by mc.fcisolver.make_rdm1() method. For state-average CASCI/CASCF calculation, this results in a set of canonicalized orbitals of state-average effective Fock matrix. To canonicalize the orbitals for one particular state, you can assign the density matrix of that state to the kwarg casdm1.
Returns:
A tuple, (natural orbitals, CI coefficients, orbital energies) The orbital energies are the diagonal terms of effective Fock matrix.
pyscf.mcscf.casci.cas_natorb(mc, mo_coeff=None, ci=None, eris=None, sort=False, casdm1=None, verbose=None, with_meta_lowdin=True)[source]

Transform active orbitals to natrual orbitals, and update the CI wfn

Args:
mc : a CASSCF/CASCI object or RHF object
Kwargs:
sort
: bool
Sort natural orbitals wrt the occupancy.
Returns:
A tuple, the first item is natural orbitals, the second is updated CI coefficients, the third is the natural occupancy associated to the natural orbitals.
pyscf.mcscf.casci.get_fock(mc, mo_coeff=None, ci=None, eris=None, casdm1=None, verbose=None)[source]

Effective one-electron Fock matrix in AO representation f = sum_{pq} E_{pq} F_{pq} F_{pq} = h_{pq} + sum_{rs} [(pq|rs)-(ps|rq)] DM_{sr}

Ref. Theor. Chim. Acta., 91, 31 Chem. Phys. 48, 157

For state-average CASCI/CASSCF object, the effective fock matrix is based on the state-average density matrix. To obtain Fock matrix of a specific state in the state-average calculations, you can pass “casdm1” of the specific state to this function.

Args:
mc: a CASSCF/CASCI object or RHF object
Kwargs:
mo_coeff (ndarray): orbitals that span the core, active and external
space.
ci (ndarray): CI coefficients (or objects to represent the CI
wavefunctions in DMRG/QMC-MCSCF calculations).
eris: Integrals for the MCSCF object. Input this object to reduce the
overhead of computing integrals. It can be generated by mc.ao2mo() method.
casdm1 (ndarray): 1-particle density matrix in active space. Without
input casdm1, the density matrix is computed with the input ci coefficients/object. If neither ci nor casdm1 were given, density matrix is computed by mc.fcisolver.make_rdm1() method. For state-average CASCI/CASCF calculation, this results in the effective Fock matrix based on the state-average density matrix. To obtain the effective Fock matrix for one particular state, you can assign the density matrix of that state to the kwarg casdm1.
Returns:
Fock matrix
pyscf.mcscf.casci.h1e_for_cas(casci, mo_coeff=None, ncas=None, ncore=None)[source]

CAS sapce one-electron hamiltonian

Args:
casci : a CASSCF/CASCI object or RHF object
Returns:
A tuple, the first is the effective one-electron hamiltonian defined in CAS space, the second is the electronic energy from core.
pyscf.mcscf.casci.kernel(casci, mo_coeff=None, ci0=None, verbose=3)[source]

CASCI solver

UCASCI (CASCI with non-degenerated alpha and beta orbitals, typically UHF orbitals)

pyscf.mcscf.ucasci.h1e_for_cas(casci, mo_coeff=None, ncas=None, ncore=None)[source]

CAS sapce one-electron hamiltonian for UHF-CASCI or UHF-CASSCF

Args:
casci : a U-CASSCF/U-CASCI object or UHF object
pyscf.mcscf.ucasci.kernel(casci, mo_coeff=None, ci0=None, verbose=3)[source]

UHF-CASCI solver

8.5.2. CASSCF

class pyscf.mcscf.mc1step.CASSCF(mf_or_mol, ncas, nelecas, ncore=None, frozen=None)[source]
Args:
mf_or_mol
: SCF object or Mole object
SCF or Mole to define the problem size.
ncas
: int
Number of active orbitals.
nelecas
: int or a pair of int
Number of electrons in active space.
Kwargs:
ncore
: int
Number of doubly occupied core orbitals. If not presented, this parameter can be automatically determined.
Attributes:
verbose
: int
Print level. Default value equals to Mole.verbose.
max_memory
: float or int
Allowed memory in MB. Default value equals to Mole.max_memory.
ncas
: int
Active space size.
nelecas
: tuple of int
Active (nelec_alpha, nelec_beta)
ncore
: int or tuple of int
Core electron number. In UHF-CASSCF, it’s a tuple to indicate the different core eletron numbers.
natorb
: bool
Whether to restore the natural orbital in CAS space. Default is not. Be very careful to set this parameter when CASCI/CASSCF are combined with DMRG solver because this parameter changes the orbital ordering which DMRG relies on.
canonicalization
: bool
Whether to canonicalize orbitals. Default is True.
sorting_mo_energy
: bool
Whether to sort the orbitals based on the diagonal elements of the general Fock matrix. Default is False.
fcisolver
: an instance of FCISolver

The pyscf.fci module provides several FCISolver for different scenario. Generally, fci.direct_spin1.FCISolver can be used for all RHF-CASSCF. However, a proper FCISolver can provide better performance and better numerical stability. One can either use fci.solver() function to pick the FCISolver by the program or manually assigen the FCISolver to this attribute, e.g.

>>> from pyscf import fci
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> mc.fcisolver = fci.solver(mol, singlet=True)
>>> mc.fcisolver = fci.direct_spin1.FCISolver(mol)

You can control FCISolver by setting e.g.:

>>> mc.fcisolver.max_cycle = 30
>>> mc.fcisolver.conv_tol = 1e-7

For more details of the parameter for FCISolver, See fci.

Saved results

e_tot
: float
Total MCSCF energy (electronic energy plus nuclear repulsion)
e_cas
: float
CAS space FCI energy
ci
: ndarray
CAS space FCI coefficients
mo_coeff
: ndarray
When canonicalization is specified, the orbitals are canonical orbitals which make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is specified, the active segment of the mo_coeff is natural orbitls.
mo_energy
: ndarray
Diagonal elements of general Fock matrix (in mo_coeff representation).

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASCI(mf, 6, 6)
>>> mc.kernel()[0]
-108.980200816243354
CASSCF

Extra attributes for CASSCF:

conv_tol
: float
Converge threshold. Default is 1e-7
conv_tol_grad
: float
Converge threshold for CI gradients and orbital rotation gradients. Default is 1e-4
max_stepsize
: float
The step size for orbital rotation. Small step (0.005 - 0.05) is prefered. Default is 0.03.
max_cycle_macro
: int
Max number of macro iterations. Default is 50.
max_cycle_micro
: int
Max number of micro iterations in each macro iteration. Depending on systems, increasing this value might reduce the total macro iterations. Generally, 2 - 5 steps should be enough. Default is 3.
ah_level_shift
: float, for AH solver.
Level shift for the Davidson diagonalization in AH solver. Default is 1e-8.
ah_conv_tol
: float, for AH solver.
converge threshold for AH solver. Default is 1e-12.
ah_max_cycle
: float, for AH solver.
Max number of iterations allowd in AH solver. Default is 30.
ah_lindep
: float, for AH solver.
Linear dependence threshold for AH solver. Default is 1e-14.
ah_start_tol
: flat, for AH solver.
In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 0.2.
ah_start_cycle
: int, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 2.

ah_conv_tol, ah_max_cycle, ah_lindep, ah_start_tol and ah_start_cycle can affect the accuracy and performance of CASSCF solver. Lower ah_conv_tol and ah_lindep might improve the accuracy of CASSCF optimization, but decrease the performance.

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.conv_tol = 1e-10
>>> mc.ah_conv_tol = 1e-5
>>> mc.kernel()[0]
-109.044401898486001
>>> mc.ah_conv_tol = 1e-10
>>> mc.kernel()[0]
-109.044401887945668
chkfile
: str
Checkpoint file to save the intermediate orbitals during the CASSCF optimization. Default is the checkpoint file of mean field object.
ci_response_space
: int
subspace size to solve the CI vector response. Default is 3.
callback
: function(envs_dict) => None
callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current envrionment.

Saved results

e_tot
: float
Total MCSCF energy (electronic energy plus nuclear repulsion)
e_cas
: float
CAS space FCI energy
ci
: ndarray
CAS space FCI coefficients
mo_coeff
: ndarray
Optimized CASSCF orbitals coefficients. When canonicalization is specified, the returned orbitals make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is specified, the active segment of the mo_coeff is natural orbitls.
mo_energy
: ndarray
Diagonal elements of general Fock matrix (in mo_coeff representation).

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.kernel()[0]
-109.044401882238134
as_scanner(mc)

Generating a scanner for CASSCF PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total CASSCF energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters of MCSCF object (conv_tol, max_memory etc) are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, ...) during calculation.

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1.2', verbose=0)
>>> mc_scanner = mcscf.CASSCF(scf.RHF(mol), 4, 4).as_scanner()
>>> e = mc_scanner(gto.M(atom='N 0 0 0; N 0 0 1.1'))
>>> e = mc_scanner(gto.M(atom='N 0 0 0; N 0 0 1.5'))
get_grad(mo_coeff=None, casdm1_casdm2=None, eris=None)[source]

Orbital gradients

get_h2cas(mo_coeff=None)[source]

Computing active space two-particle Hamiltonian.

Note It is different to get_h2eff when df.approx_hessian is applied, in which get_h2eff function returns the DF integrals while get_h2cas returns the regular 2-electron integrals.

get_h2eff(mo_coeff=None)[source]

Computing active space two-particle Hamiltonian.

Note It is different to get_h2cas when df.approx_hessian is applied, in which get_h2eff function returns the DF integrals while get_h2cas returns the regular 2-electron integrals.

kernel(mo_coeff=None, ci0=None, callback=None, _kern=<function kernel>)[source]
Returns:
Five elements, they are total energy, active space CI energy, the active space FCI wavefunction coefficients or DMRG wavefunction ID, the MCSCF canonical orbital coefficients, the MCSCF canonical orbital coefficients.

They are attributes of mcscf object, which can be accessed by .e_tot, .e_cas, .ci, .mo_coeff, .mo_energy

rotate_mo(mo, u, log=None)[source]

Rotate orbitals with the given unitary matrix

solve_approx_ci(h1, h2, ci0, ecore, e_cas, envs)[source]

Solve CI eigenvalue/response problem approximately

pyscf.mcscf.mc1step.as_scanner(mc)[source]

Generating a scanner for CASSCF PES.

The returned solver is a function. This function requires one argument “mol” as input and returns total CASSCF energy.

The solver will automatically use the results of last calculation as the initial guess of the new calculation. All parameters of MCSCF object (conv_tol, max_memory etc) are automatically applied in the solver.

Note scanner has side effects. It may change many underlying objects (_scf, with_df, with_x2c, ...) during calculation.

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1.2', verbose=0)
>>> mc_scanner = mcscf.CASSCF(scf.RHF(mol), 4, 4).as_scanner()
>>> e = mc_scanner(gto.M(atom='N 0 0 0; N 0 0 1.1'))
>>> e = mc_scanner(gto.M(atom='N 0 0 0; N 0 0 1.5'))
pyscf.mcscf.mc1step.kernel(casscf, mo_coeff, tol=1e-07, conv_tol_grad=None, ci0=None, callback=None, verbose=3, dump_chk=True)[source]

quasi-newton CASSCF optimization driver

pyscf.mcscf.mc1step_symm.CASSCF

alias of SymAdaptedCASSCF

class pyscf.mcscf.mc1step_symm.SymAdaptedCASSCF(mf_or_mol, ncas, nelecas, ncore=None, frozen=None)[source]
Args:
mf_or_mol
: SCF object or Mole object
SCF or Mole to define the problem size.
ncas
: int
Number of active orbitals.
nelecas
: int or a pair of int
Number of electrons in active space.
Kwargs:
ncore
: int
Number of doubly occupied core orbitals. If not presented, this parameter can be automatically determined.
Attributes:
verbose
: int
Print level. Default value equals to Mole.verbose.
max_memory
: float or int
Allowed memory in MB. Default value equals to Mole.max_memory.
ncas
: int
Active space size.
nelecas
: tuple of int
Active (nelec_alpha, nelec_beta)
ncore
: int or tuple of int
Core electron number. In UHF-CASSCF, it’s a tuple to indicate the different core eletron numbers.
natorb
: bool
Whether to restore the natural orbital in CAS space. Default is not. Be very careful to set this parameter when CASCI/CASSCF are combined with DMRG solver because this parameter changes the orbital ordering which DMRG relies on.
canonicalization
: bool
Whether to canonicalize orbitals. Default is True.
sorting_mo_energy
: bool
Whether to sort the orbitals based on the diagonal elements of the general Fock matrix. Default is False.
fcisolver
: an instance of FCISolver

The pyscf.fci module provides several FCISolver for different scenario. Generally, fci.direct_spin1.FCISolver can be used for all RHF-CASSCF. However, a proper FCISolver can provide better performance and better numerical stability. One can either use fci.solver() function to pick the FCISolver by the program or manually assigen the FCISolver to this attribute, e.g.

>>> from pyscf import fci
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> mc.fcisolver = fci.solver(mol, singlet=True)
>>> mc.fcisolver = fci.direct_spin1.FCISolver(mol)

You can control FCISolver by setting e.g.:

>>> mc.fcisolver.max_cycle = 30
>>> mc.fcisolver.conv_tol = 1e-7

For more details of the parameter for FCISolver, See fci.

Saved results

e_tot
: float
Total MCSCF energy (electronic energy plus nuclear repulsion)
e_cas
: float
CAS space FCI energy
ci
: ndarray
CAS space FCI coefficients
mo_coeff
: ndarray
When canonicalization is specified, the orbitals are canonical orbitals which make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is specified, the active segment of the mo_coeff is natural orbitls.
mo_energy
: ndarray
Diagonal elements of general Fock matrix (in mo_coeff representation).

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASCI(mf, 6, 6)
>>> mc.kernel()[0]
-108.980200816243354
CASSCF

Extra attributes for CASSCF:

conv_tol
: float
Converge threshold. Default is 1e-7
conv_tol_grad
: float
Converge threshold for CI gradients and orbital rotation gradients. Default is 1e-4
max_stepsize
: float
The step size for orbital rotation. Small step (0.005 - 0.05) is prefered. Default is 0.03.
max_cycle_macro
: int
Max number of macro iterations. Default is 50.
max_cycle_micro
: int
Max number of micro iterations in each macro iteration. Depending on systems, increasing this value might reduce the total macro iterations. Generally, 2 - 5 steps should be enough. Default is 3.
ah_level_shift
: float, for AH solver.
Level shift for the Davidson diagonalization in AH solver. Default is 1e-8.
ah_conv_tol
: float, for AH solver.
converge threshold for AH solver. Default is 1e-12.
ah_max_cycle
: float, for AH solver.
Max number of iterations allowd in AH solver. Default is 30.
ah_lindep
: float, for AH solver.
Linear dependence threshold for AH solver. Default is 1e-14.
ah_start_tol
: flat, for AH solver.
In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 0.2.
ah_start_cycle
: int, for AH solver.

In AH solver, the orbital rotation is started without completely solving the AH problem. This value is to control the start point. Default is 2.

ah_conv_tol, ah_max_cycle, ah_lindep, ah_start_tol and ah_start_cycle can affect the accuracy and performance of CASSCF solver. Lower ah_conv_tol and ah_lindep might improve the accuracy of CASSCF optimization, but decrease the performance.

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.UHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.conv_tol = 1e-10
>>> mc.ah_conv_tol = 1e-5
>>> mc.kernel()[0]
-109.044401898486001
>>> mc.ah_conv_tol = 1e-10
>>> mc.kernel()[0]
-109.044401887945668
chkfile
: str
Checkpoint file to save the intermediate orbitals during the CASSCF optimization. Default is the checkpoint file of mean field object.
ci_response_space
: int
subspace size to solve the CI vector response. Default is 3.
callback
: function(envs_dict) => None
callback function takes one dict as the argument which is generated by the builtin function locals(), so that the callback function can access all local variables in the current envrionment.

Saved results

e_tot
: float
Total MCSCF energy (electronic energy plus nuclear repulsion)
e_cas
: float
CAS space FCI energy
ci
: ndarray
CAS space FCI coefficients
mo_coeff
: ndarray
Optimized CASSCF orbitals coefficients. When canonicalization is specified, the returned orbitals make the general Fock matrix (Fock operator on top of MCSCF 1-particle density matrix) diagonalized within each subspace (core, active, external). If natorb (natural orbitals in active space) is specified, the active segment of the mo_coeff is natural orbitls.
mo_energy
: ndarray
Diagonal elements of general Fock matrix (in mo_coeff representation).

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> mc.kernel()[0]
-109.044401882238134
rotate_mo(mo, u, log=None)[source]

Rotate orbitals with the given unitary matrix

sort_mo_by_irrep(cas_irrep_nocc, cas_irrep_ncore=None, mo_coeff=None, s=None)[source]

Select active space based on symmetry information. See also pyscf.mcscf.addons.sort_mo_by_irrep()

UCASSCF (CASSCF without spin-degeneracy between alpha and beta orbitals) 1-step optimization algorithm

MO integrals for UCASSCF methods

8.5.3. addons

pyscf.mcscf.addons.cas_natorb(casscf, mo_coeff=None, ci=None, sort=False)[source]

Natrual orbitals in CAS space

pyscf.mcscf.addons.caslst_by_irrep(casscf, mo_coeff, cas_irrep_nocc, cas_irrep_ncore=None, s=None, base=1)[source]

Given number of active orbitals for each irrep, return the orbital indices of active space

Args:

casscf : an CASSCF or CASCI object

cas_irrep_nocc
: list or dict
Number of active orbitals for each irrep. It can be a dict, eg {‘A1’: 2, ‘B2’: 4} to indicate the active space size based on irrep names, or {0: 2, 3: 4} for irrep Id, or a list [2, 0, 0, 4] (identical to {0: 2, 3: 4}) in which the list index is served as the irrep Id.
Kwargs:
cas_irrep_ncore
: list or dict
Number of closed shells for each irrep. It can be a dict, eg {‘A1’: 6, ‘B2’: 4} to indicate the closed shells based on irrep names, or {0: 6, 3: 4} for irrep Id, or a list [6, 0, 0, 4] (identical to {0: 6, 3: 4}) in which the list index is served as the irrep Id. If cas_irrep_ncore is not given, the program will generate a guess based on the lowest CASCI.ncore orbitals.
s
: ndarray
overlap matrix
base
: int
0-based (C-like) or 1-based (Fortran-like) caslst
Returns:
A list of orbital indices

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvtz', symmetry=True, verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.kernel()
>>> mc = mcscf.CASSCF(mf, 12, 4)
>>> mcscf.caslst_by_irrep(mc, mf.mo_coeff, {'E1gx':4, 'E1gy':4, 'E1ux':2, 'E1uy':2})
[5, 7, 8, 10, 11, 14, 15, 20, 25, 26, 31, 32]
pyscf.mcscf.addons.get_fock(casscf, mo_coeff=None, ci=None)[source]

Generalized Fock matrix in AO representation

pyscf.mcscf.addons.make_rdm1(casscf, mo_coeff=None, ci=None, **kwargs)[source]

One-particle densit matrix in AO representation

Args:
casscf : an CASSCF or CASCI object
Kwargs:
ci
: ndarray
CAS space FCI coefficients. If not given, take casscf.ci.
mo_coeff
: ndarray
Orbital coefficients. If not given, take casscf.mo_coeff.

Examples:

>>> import scipy.linalg
>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='sto-3g', verbose=0)
>>> mf = scf.RHF(mol)
>>> res = mf.scf()
>>> mc = mcscf.CASSCF(mf, 6, 6)
>>> res = mc.kernel()
>>> natocc = numpy.linalg.eigh(mcscf.make_rdm1(mc), mf.get_ovlp(), type=2)[0]
>>> print(natocc)
[ 0.0121563   0.0494735   0.0494735   1.95040395  1.95040395  1.98808879
  2.          2.          2.          2.        ]
pyscf.mcscf.addons.make_rdm1s(casscf, mo_coeff=None, ci=None, **kwargs)[source]

Alpha and beta one-particle densit matrices in AO representation

pyscf.mcscf.addons.map2hf(casscf, mf_mo=None, base=1, tol=0.4)[source]

The overlap between the CASSCF optimized orbitals and the canonical HF orbitals.

pyscf.mcscf.addons.project_init_guess(casscf, init_mo, prev_mol=None)[source]

Project the given initial guess to the current CASSCF problem. The projected initial guess has two parts. The core orbitals are directly taken from the Hartree-Fock orbitals, and the active orbitals are projected from the given initial guess.

Args:

casscf : an CASSCF or CASCI object

init_mo
: ndarray or list of ndarray
Initial guess orbitals which are not orth-normal for the current molecule. When the casscf is UHF-CASSCF, the init_mo needs to be a list of two ndarrays, for alpha and beta orbitals
Kwargs:
prev_mol
: an instance of Mole
If given, the inital guess orbitals are associated to the geometry and basis of prev_mol. Otherwise, the orbitals are based of the geometry and basis of casscf.mol
Returns:
New orthogonal initial guess orbitals with the core taken from Hartree-Fock orbitals and projected active space from original initial guess orbitals

Examples:

import numpy
from pyscf import gto, scf, mcscf
mol = gto.Mole()
mol.build(atom='H 0 0 0; F 0 0 0.8', basis='ccpvdz', verbose=0)
mf = scf.RHF(mol)
mf.scf()
mc = mcscf.CASSCF(mf, 6, 6)
mo = mcscf.sort_mo(mc, mf.mo_coeff, [3,4,5,6,8,9])
print('E(0.8) = %.12f' % mc.kernel(mo)[0])
init_mo = mc.mo_coeff
for b in numpy.arange(1.0, 3., .2):
    mol.atom = [['H', (0, 0, 0)], ['F', (0, 0, b)]]
    mol.build(0, 0)
    mf = scf.RHF(mol)
    mf.scf()
    mc = mcscf.CASSCF(mf, 6, 6)
    mo = mcscf.project_init_guess(mc, init_mo)
    print('E(%2.1f) = %.12f' % (b, mc.kernel(mo)[0]))
    init_mo = mc.mo_coeff
pyscf.mcscf.addons.sort_mo(casscf, mo_coeff, caslst, base=1)[source]

Pick orbitals for CAS space

Args:

casscf : an CASSCF or CASCI object

mo_coeff
: ndarray or a list of ndarray
Orbitals for CASSCF initial guess. In the UHF-CASSCF, it’s a list of two orbitals, for alpha and beta spin.
caslst
: list of int or nested list of int
A list of orbital indices to represent the CAS space. In the UHF-CASSCF, it’s consist of two lists, for alpha and beta spin.
Kwargs:
base
: int
0-based (C-style) or 1-based (Fortran-style) caslst
Returns:
An reoreded mo_coeff, which put the orbitals given by caslst in the CAS space

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvdz', verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.scf()
>>> mc = mcscf.CASSCF(mf, 4, 4)
>>> cas_list = [5,6,8,9] # pi orbitals
>>> mo = mc.sort_mo(cas_list)
>>> mc.kernel(mo)[0]
-109.007378939813691
pyscf.mcscf.addons.sort_mo_by_irrep(casscf, mo_coeff, cas_irrep_nocc, cas_irrep_ncore=None, s=None)[source]

Given number of active orbitals for each irrep, construct the mo initial guess for CASSCF

Args:

casscf : an CASSCF or CASCI object

cas_irrep_nocc
: list or dict
Number of active orbitals for each irrep. It can be a dict, eg {‘A1’: 2, ‘B2’: 4} to indicate the active space size based on irrep names, or {0: 2, 3: 4} for irrep Id, or a list [2, 0, 0, 4] (identical to {0: 2, 3: 4}) in which the list index is served as the irrep Id.
Kwargs:
cas_irrep_ncore
: list or dict
Number of closed shells for each irrep. It can be a dict, eg {‘A1’: 6, ‘B2’: 4} to indicate the closed shells based on irrep names, or {0: 6, 3: 4} for irrep Id, or a list [6, 0, 0, 4] (identical to {0: 6, 3: 4}) in which the list index is served as the irrep Id. If cas_irrep_ncore is not given, the program will generate a guess based on the lowest CASCI.ncore orbitals.
s
: ndarray
overlap matrix
Returns:
sorted orbitals, ordered as [c,..,c,a,..,a,v,..,v]

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='N 0 0 0; N 0 0 1', basis='ccpvtz', symmetry=True, verbose=0)
>>> mf = scf.RHF(mol)
>>> mf.kernel()
>>> mc = mcscf.CASSCF(mf, 12, 4)
>>> mo = mc.sort_mo_by_irrep({'E1gx':4, 'E1gy':4, 'E1ux':2, 'E1uy':2})
>>> # Same to mo = sort_mo_by_irrep(mc, mf.mo_coeff, {2: 4, 3: 4, 6: 2, 7: 2})
>>> # Same to mo = sort_mo_by_irrep(mc, mf.mo_coeff, [0, 0, 4, 4, 0, 0, 2, 2])
>>> mc.kernel(mo)[0]
-108.162863845084
pyscf.mcscf.addons.spin_square(casscf, mo_coeff=None, ci=None, ovlp=None)[source]

Spin square of the UHF-CASSCF wavefunction

Returns:
A list of two floats. The first is the expectation value of S^2. The second is the corresponding 2S+1

Examples:

>>> from pyscf import gto, scf, mcscf
>>> mol = gto.M(atom='O 0 0 0; O 0 0 1', basis='sto-3g', spin=2, verbose=0)
>>> mf = scf.UHF(mol)
>>> res = mf.scf()
>>> mc = mcscf.CASSCF(mf, 4, 6)
>>> res = mc.kernel()
>>> print('S^2 = %.7f, 2S+1 = %.7f' % mcscf.spin_square(mc))
S^2 = 3.9831589, 2S+1 = 4.1149284
pyscf.mcscf.addons.state_average(casscf, weights=(0.5, 0.5))

State average over the energy. The energy funcitonal is E = w1<psi1|H|psi1> + w2<psi2|H|psi2> + ...

Note we may need change the FCI solver to

mc.fcisolver = fci.solver(mol, False)

before calling state_average_(mc), to mix the singlet and triplet states

pyscf.mcscf.addons.state_average_(casscf, weights=(0.5, 0.5))[source]

State average over the energy. The energy funcitonal is E = w1<psi1|H|psi1> + w2<psi2|H|psi2> + ...

Note we may need change the FCI solver to

mc.fcisolver = fci.solver(mol, False)

before calling state_average_(mc), to mix the singlet and triplet states

pyscf.mcscf.addons.state_average_mix(casscf, fcisolvers, weights=(0.5, 0.5))

State-average CASSCF over multiple FCI solvers.

pyscf.mcscf.addons.state_average_mix_(casscf, fcisolvers, weights=(0.5, 0.5))[source]

State-average CASSCF over multiple FCI solvers.

pyscf.mcscf.addons.state_specific(casscf, state=1)

For excited state

Kwargs:
state : int 0 for ground state; 1 for first excited state.
pyscf.mcscf.addons.state_specific_(casscf, state=1)[source]

For excited state

Kwargs:
state : int 0 for ground state; 1 for first excited state.