Selecting an active space#
Selecting the active space is the main step that distinguishes a multireference calculation from a standard quantum chemical calculation. The choice of the orbitals that make up the space is in general guided by intuition and/or by the calculation. Here we will not explain how to choose the active orbitals but how to efficiently select the active space in multipsi.
Starting orbitals#
The choice of the active space is intimately linked to the choice of starting orbitals. Different orbitals may make the selection of the active space more or less easy, but also how quickly the MCSCF calculation will converge. Here are a few typical solutions
From an orbital guess#
MultiPsi implements a form of extended Hückel module to create an orbital guess. This has the advantage of being very quick, a few seconds at most, and to provide reasonably localized starting orbitals. On the other hand, they are not as good quality as Hartree-Fock or DFT orbitals and thus will take more time to converge in the subsequent MCSCF. Let’s try with a furan molecule.
import veloxchem as vlx
import multipsi as mtp
furan_xyz="""9
C -0.86213 -0.90784 0.00007
H -1.63433 -1.64264 -0.00003
C 0.50727 -0.90524 0.00007
C 0.92057 0.47886 -0.00003
C -0.22323 1.23186 -0.00003
O -1.35123 0.40376 -0.00013
H 1.17117 -1.74724 0.00017
H 1.93767 0.81866 0.00007
H -0.46573 2.26986 -0.00013
"""
molecule = vlx.Molecule.from_xyz_string(furan_xyz)
basis = vlx.MolecularBasis.read(molecule,"def2-sv(p)")
#Generate the extended Hückel guess
guessdriver=mtp.OrbitalGuess()
guessorb=guessdriver.compute(molecule,basis)
#Visualize the orbitals using Veloxchel OrbitalViewer
orbviewer=vlx.OrbitalViewer()
orbviewer.plot(molecule,basis,guessorb)
# Run a simple HF to check the quality of the orbitals
# If no active space is given, it assumes RHF/ROHF
space=mtp.OrbSpace(molecule,guessorb)
mcscfdrv=mtp.McscfDriver()
mc_results = mcscfdrv.compute(molecule,basis,space)
Multi-Configurational Self-Consistent Field Driver
====================================================
Active space definition:
------------------------
Number of inactive (occupied) orbitals: 18
Number of active orbitals: 0
Number of virtual orbitals: 60
This is a restricted Hartree-Fock wavefunction
CI expansion:
-------------
Number of determinants: 1
╭────────────────────────────────────╮
│ Driver settings │
╰────────────────────────────────────╯
Number of states : 1
State-specific calculation
- State of interest : 1
Max. iterations : 50
BFGS window : 5
Convergence thresholds:
- Energy change : 1e-08
- Gradient sq. norm : 1e-08
MCSCF Iterations
-------------------
Iter. | Average Energy | E. Change | Grad. Norm | CI Iter. | Trust rad. | Time
---------------------------------------------------------------------------------
-228.0185548989264
1 -228.018554899 0.0e+00 1.6e+00 0 0.40 0:00:00
-228.34604876591195
2 -228.346048766 -3.3e-01 2.9e-01 0 0.40 0:00:00
-228.42353668472376
3 -228.423536685 -7.7e-02 4.5e-02 0 0.40 0:00:00
-228.43246519726623
4 -228.432465197 -8.9e-03 5.0e-03 0 0.40 0:00:00
-228.43352739873924
5 -228.433527399 -1.1e-03 1.2e-04 0 0.40 0:00:00
-228.43356216332572
6 -228.433562163 -3.5e-05 2.4e-05 0 0.40 0:00:00
-228.4335685126977
7 -228.433568513 -6.3e-06 7.4e-07 0 0.40 0:00:00
-228.4335687301369
8 -228.433568730 -2.2e-07 3.3e-08 0 0.40 0:00:00
-228.4335687560168
9 -228.433568756 -2.6e-08 4.2e-09 0 0.48 0:00:00
-228.43356875779236
10 -228.433568758 -1.8e-09 8.2e-10 0 0.48 0:00:00
** Convergence reached in 10 iterations
-228.43356875813006
11 -228.433568758 -3.4e-10 1.4e-10 0 0.48 0:00:00
Final results
-------------
* State 1
- S^2 : 0.00 (multiplicity = 1.0 )
- Energy : -228.43356875813006
Spin Restricted Orbitals
------------------------
Molecular Orbital No. 14:
--------------------------
Occupation: 2.000 Energy: -0.58793 a.u.
( 1 C 1p+1: 0.25) ( 3 C 1p+1: -0.27) ( 4 C 1p+1: 0.30)
( 5 C 1p+1: -0.18) ( 5 C 1p-1: 0.18) ( 7 H 1s : -0.25)
( 8 H 1s : 0.25)
Molecular Orbital No. 15:
--------------------------
Occupation: 2.000 Energy: -0.56268 a.u.
( 3 C 1p-1: -0.31) ( 4 C 1p+1: 0.20) ( 4 C 1p-1: 0.24)
( 6 O 1p+1: 0.37) ( 6 O 2p+1: 0.27) ( 7 H 1s : 0.20)
( 8 H 1s : 0.20)
Molecular Orbital No. 16:
--------------------------
Occupation: 2.000 Energy: -0.54763 a.u.
( 1 C 1p+1: -0.25) ( 3 C 1p+1: 0.23) ( 3 C 1p-1: 0.20)
( 4 C 1p-1: -0.29) ( 5 C 1p+1: -0.24) ( 6 O 1p+1: 0.32)
( 6 O 2p+1: 0.24)
Molecular Orbital No. 17:
--------------------------
Occupation: 2.000 Energy: -0.39986 a.u.
( 3 C 1p0 : -0.32) ( 3 C 2p0 : -0.24) ( 4 C 1p0 : -0.32)
( 4 C 2p0 : -0.24) ( 6 O 1p0 : 0.34) ( 6 O 2p0 : 0.30)
Molecular Orbital No. 18:
--------------------------
Occupation: 2.000 Energy: -0.32333 a.u.
( 1 C 1p0 : -0.35) ( 1 C 2p0 : -0.29) ( 3 C 1p0 : -0.23)
( 3 C 2p0 : -0.21) ( 4 C 1p0 : 0.23) ( 4 C 2p0 : 0.21)
( 5 C 1p0 : 0.35) ( 5 C 2p0 : 0.29)
Molecular Orbital No. 19:
--------------------------
Occupation: 0.000 Energy: 0.14744 a.u.
( 1 C 1p0 : 0.32) ( 1 C 2p0 : 0.56) ( 3 C 1p0 : -0.18)
( 3 C 2p0 : -0.35) ( 4 C 1p0 : -0.18) ( 4 C 2p0 : -0.35)
( 5 C 1p0 : 0.32) ( 5 C 2p0 : 0.56) ( 6 O 1p0 : -0.27)
( 6 O 2p0 : -0.39)
Molecular Orbital No. 20:
--------------------------
Occupation: 0.000 Energy: 0.20445 a.u.
( 1 C 3s : 0.65) ( 1 C 2p+1: -0.30) ( 1 C 2p-1: -0.23)
( 2 H 2s : -0.98) ( 3 C 3s : 0.97) ( 3 C 2p+1: 0.29)
( 3 C 2p-1: -0.39) ( 4 C 3s : 0.97) ( 4 C 2p+1: 0.46)
( 4 C 2p-1: 0.17) ( 5 C 3s : 0.65) ( 5 C 2p-1: 0.36)
( 7 H 2s : -1.31) ( 8 H 2s : -1.31) ( 9 H 2s : -0.98)
Molecular Orbital No. 21:
--------------------------
Occupation: 0.000 Energy: 0.22624 a.u.
( 1 C 1p0 : 0.22) ( 1 C 2p0 : 0.54) ( 3 C 1p0 : -0.32)
( 3 C 2p0 : -0.83) ( 4 C 1p0 : 0.32) ( 4 C 2p0 : 0.83)
( 5 C 1p0 : -0.22) ( 5 C 2p0 : -0.54)
Molecular Orbital No. 22:
--------------------------
Occupation: 0.000 Energy: 0.23036 a.u.
( 1 C 3s : 0.73) ( 1 C 2p+1: -0.85) ( 1 C 2p-1: -0.79)
( 2 H 2s : -1.89) ( 3 C 3s : 1.11) ( 3 C 2p-1: 0.24)
( 4 C 3s : -1.11) ( 4 C 2p-1: 0.20) ( 5 C 3s : -0.73)
( 5 C 1p-1: -0.16) ( 5 C 2p+1: 0.28) ( 5 C 2p-1: -1.13)
( 7 H 2s : -0.36) ( 8 H 2s : 0.36) ( 9 H 2s : 1.89)
Molecular Orbital No. 23:
--------------------------
Occupation: 0.000 Energy: 0.23249 a.u.
( 1 C 3s : -0.93) ( 1 C 2p+1: 0.62) ( 1 C 2p-1: 0.57)
( 2 H 2s : 1.56) ( 3 C 3s : 0.60) ( 3 C 2p+1: 0.33)
( 3 C 2p-1: -0.60) ( 4 C 3s : 0.60) ( 4 C 2p+1: 0.60)
( 4 C 2p-1: 0.32) ( 5 C 3s : -0.93) ( 5 C 2p+1: 0.20)
( 5 C 2p-1: -0.81) ( 7 H 2s : -1.14) ( 8 H 2s : -1.14)
( 9 H 2s : 1.56)
Dipole moment for state 1
---------------------------
X : 0.383480 a.u. 0.974710 Debye
Y : -0.114775 a.u. -0.291728 Debye
Z : 0.000066 a.u. 0.000167 Debye
Total : 0.400288 a.u. 1.017431 Debye
Total MCSCF time: 00:00:02
The guess is nearly instant, and the orbitals are visually easy to analyze. The starting orbitals are not as good quality as those coming from a Hartree-Fock calculation, but the speed and locality of this guess makes it very convenient in a typical workflow.
From HF or DFT canonical orbitals#
Hartree-Fock or Kohn-Sham canonical orbitals are usually a good guess for a MCSCF calculation as they are already reasonably well converged. However, those orbitals are typically very delocalized, which can make it difficult to find the specific orbitals we want. This is particularly true for virtual orbitals. Additionally, in strong correlation cases, convergence of the HF or DFT may be very slow and difficult.
From natural orbitals#
Natural orbitals are often one of the best starting point for a MCSCF calculation. They can be obtained from various means, in particular:
an unrestricted HF or DFT calculation (only for open-shell systems)
a MP2 or other correlated wavefunction methods
another MCSCF calculation
The main advantage is that they provide occupation numbers which are a very intuitive metric of correlation, and we usually want to include the most strongly correlated orbitals in our active space. While unrestricted HF or DFT calculations are not usually counted as “correlated” in the same sense as MP2 or CI, the natural orbital occupation numbers are often a decent indicator of correlation.
Let’s illustrate this with the UHF natural orbitals of NO radical.
NO_xyz="""2
N 0.0 0.0 0.0
O 0.0 0.0 1.157
"""
molecule=vlx.Molecule.from_xyz_string(NO_xyz)
molecule.set_multiplicity(2)
basis = vlx.MolecularBasis.read(molecule,"cc-pvdz")
scfdrv = vlx.ScfUnrestrictedDriver()
hf_results = scfdrv.compute(molecule, basis)
uhf_orbs = scfdrv.natural_orbitals()
Self Consistent Field Driver Setup
====================================
Wave Function Model : Spin-Unrestricted Hartree-Fock
Initial Guess Model : Superposition of Atomic Densities
Convergence Accelerator : Two Level Direct Inversion of Iterative Subspace
Max. Number of Iterations : 50
Max. Number of Error Vectors : 10
Convergence Threshold : 1.0e-06
ERI Screening Threshold : 1.0e-12
Linear Dependence Threshold : 1.0e-06
* Info * Starting Reduced Basis SCF calculation...
* Info * ...done. SCF energy in reduced basis set: -129.185387517916 a.u. Time: 0.26 sec.
Iter. | Hartree-Fock Energy | Energy Change | Gradient Norm | Max. Gradient | Density Change
--------------------------------------------------------------------------------------------
1 -129.256758349985 0.0000000000 0.13580350 0.01935910 0.00000000
2 -129.258795229408 -0.0020368794 0.03970025 0.00682448 0.05175816
3 -129.258974312365 -0.0001790830 0.01961846 0.00354177 0.02008699
4 -129.259027341703 -0.0000530293 0.01015897 0.00264867 0.00683259
5 -129.259066368441 -0.0000390267 0.00678827 0.00140843 0.00569010
6 -129.259081754832 -0.0000153864 0.00667155 0.00141568 0.00205035
7 -129.259096040860 -0.0000142860 0.00665601 0.00148991 0.00199835
8 -129.259109188235 -0.0000131474 0.00671596 0.00156346 0.00191338
9 -129.259073622353 0.0000355659 0.00714468 0.00210129 0.00581311
10 -129.259432418381 -0.0003587960 0.01216600 0.00233020 0.16193419
11 -129.259527889583 -0.0000954712 0.00233284 0.00050382 0.06620341
12 -129.259549124989 -0.0000212354 0.00077024 0.00024174 0.02969819
13 -129.259546755615 0.0000023694 0.00115582 0.00028734 0.01756310
14 -129.259551602414 -0.0000048468 0.00030041 0.00004015 0.01065286
15 -129.259551578633 0.0000000238 0.00017055 0.00003338 0.00068151
16 -129.259551590727 -0.0000000121 0.00009088 0.00002217 0.00162443
17 -129.259551616729 -0.0000000260 0.00001915 0.00000378 0.00086446
18 -129.259551617285 -0.0000000006 0.00000389 0.00000069 0.00009042
19 -129.259551617298 -0.0000000000 0.00000031 0.00000008 0.00001973
*** SCF converged in 19 iterations. Time: 0.97 sec.
Spin-Unrestricted Hartree-Fock:
-------------------------------
Total Energy : -129.2595516173 a.u.
Electronic Energy : -154.8722774691 a.u.
Nuclear Repulsion Energy : 25.6127258518 a.u.
------------------------------------
Gradient Norm : 0.0000003053 a.u.
Ground State Information
------------------------
Charge of Molecule : 0.0
Multiplicity (2S+1) : 2
Magnetic Quantum Number (M_S) : 0.5
Expectation value of S**2 : 0.8133
orbviewer=vlx.OrbitalViewer()
orbviewer.plot_using_k3d(molecule,basis,uhf_orbs)
We can see that the natural orbitals identified the singly-occupied orbital of the NO radical, but also a correlated pair of the \(\pi\) and \(\pi^*\) (occupation respectively 1.973 and 0.027).
Selecting the active space#
The choice of the active space is executed by the OrbSpace class. This class contains several functions to help the user define the active space and the type of CI expansion.
import multipsi as mtp
help(mtp.OrbSpace)
Help on class OrbSpace in module multipsi.orbspace:
class OrbSpace(builtins.object)
| OrbSpace(molecule, molecular_orbitals, comm=<mpi4py.MPI.Intracomm object at 0x10b1e3af0>)
|
| Defines and handles orbital spaces
|
| :param molecule:
| The molecule.
| :param molecular_orbitals:
| The molecular orbitals or the name of a h5 orbital file.
| :param comm:
| The MPI communicator (default: world)
|
| Methods defined here:
|
| __init__(self, molecule, molecular_orbitals, comm=<mpi4py.MPI.Intracomm object at 0x10b1e3af0>)
| Initialize self. See help(type(self)) for accurate signature.
|
| __str__(self)
| :return:
| a string with active space description.
|
| cas(self, n_electrons, n_active)
| Define a complete active space with n_electrons electrons in n_active orbitals
|
| :param n_electrons:
| Number of active electrons.
| :param n_active:
| Number of active orbitals.
|
| cas_orbitals(self, active_list, n_active_electrons=None)
| Choose the orbitals that constitute the active space and define a complete active space model.
|
| :param active_list:
| The list of orbital indices starting at 0.
| :param n_active_electrons:
| The number of active electrons. If not set, the number is chosen
| based on the current occupation of the chosen orbitals.
|
| ci(self, n_excitations, n_frozen=0)
| Define a truncated CI to given order
| for example CI(4) corresponds to CISDTQ
|
| :param n_excitations:
| The maximum excitation order.
| :param n_frozen:
| Number of frozen orbitals (default: 0).
|
| cis(self, n_frozen=0)
| Define a configuration interaction with single excitations
|
| :param n_frozen:
| Number of frozen orbitals (default: 0).
|
| cisd(self, n_frozen=0)
| Define a configuration interaction with single and double excitations
|
| :param n_frozen:
| Number of frozen orbitals (default: 0).
|
| copy(self)
| Create a duplicate of this OrbSpace
|
| fci(self, n_frozen=0)
| Define a full configuration interaction space
| i.e. using all orbitals (except frozen ones)
|
| :param n_frozen:
| Number of frozen orbitals (default: 0).
|
| gas(self, gas_list)
| Define a generalised active space.
|
| :param gas_list:
| a list of each space with cumulated number of orbitals
| and minimum and maximum cumulated number of electrons.
| ex: GAS([[4,6,8],[8,8,8]]) corresponds to a RAS calculation with 4 orbitals in ras1, 4 in ras3
| and up to 2 excitations between ras1 and ras3.
|
| get_ordered_mo_coefs(self)
| Get molecular orbital coefficients sorted as
| inactive - active - virtual
|
| :return
| The mo coefficients matrix.
|
| get_total_density(self, active_1dm)
| Get the total density in AO basis
| given the active density in MO basis
|
| :param active_1dm:
| The active density in MO basis.
| :return
| The AO total density.
|
| guess_act_den(self)
| Guess the active densities
|
| :return
| The 1 and 2 pdm guess
|
| mrci(self, n_excitations, n_frozen=0)
| Define an (uncontracted) multi-reference configuration interaction to a given order
| For example, MRCI(2) corresponds to MRCISD
| This command assumes a previous CASSCF
|
| :param n_excitations:
| The maximum excitation order.
| :param n_frozen:
| Number of frozen orbitals (default: 0).
|
| mrcis(self, n_frozen=0)
| Define an (uncontracted) multi-reference configuration interaction with single excitations
| This command assumes a previous CASSCF
|
| :param n_frozen:
| Number of frozen orbitals (default: 0).
|
| mrcisd(self, n_frozen=0)
| Define an (uncontracted) multi-reference configuration interaction with single and double excitations
| This command assumes a previous CASSCF
|
| :param n_frozen:
| Number of frozen orbitals (default: 0).
|
| orthonormalize(self, overlap)
| Orthonormalize the orbitals.
|
| :param overlap:
| The overlap matrix.
|
| pick_active(self, active_list)
| Choose the orbitals that constitute the active space.
| This comes after the definition of the space (CAS, or GAS).
| Otherwise, by default, the code assumes orbitals ordered as inactive-active-virtual.
|
| :param active_list:
| The list of orbital indices starting at 0.
|
| project_basis(self, basis1, basis2)
| Project the orbitals from basis1 to basis2
|
| :param basis1:
| The current basis
| :param basis2:
| The new basis
|
| ras(self, n_electrons, max_ras1_holes, max_ras3_electrons, n_ras1, n_ras2, n_ras3)
| Define a restricted active space divided into 3 spaces called ras1, ras2 and ras3
| with restriction on the number of holes and electrons in ras1 resp. ras3.
|
| :param n_electrons:
| Number of active electrons.
| :param max_ras1_holes:
| The maximum number of holes allowed in ras1.
| :param max_ras3_holes:
| The maximum number of electrons allowed in ras3.
| :param n_ras1:
| The number of orbitals in ras1.
| :param n_ras2:
| The number of orbitals in ras2.
| :param n_ras3:
| The number of orbitals in ras3.
|
| restrict(self, n_excitations=2)
| Return a copy of the current CAS space but restricting it to a maximum of excitations
|
| :param n_excitations:
| The maximum number of excitations allowed (default 2)
|
| write_hdf5(self, file_name)
| Save the OrbSpace into a h5 file.
|
| :param file_name:
| The output h5 file name.
|
| ----------------------------------------------------------------------
| Data descriptors defined here:
|
| __dict__
| dictionary for instance variables
|
| __weakref__
| list of weak references to the object
This class is initialized using a VeloxChem Molecule object and either a VeloxChem MolecularOrbital or the name of a h5 file containing the orbitals.
Using the number of active electrons and orbitals#
The simplest way to define an active space in MultiPsi is by providing the number of active electrons and orbitals. For example, for a complete active space, one can use the function CAS(\(m,n\)) with \(m\) the number of active electrons and \(n\) the number of active orbitals. In this case, the code will select the orbitals closest in order to the HOMO-LUMO gap. It is not rare that the active orbitals we want are precisely around the HOMO-LUMO gap, especially if one uses natural orbitals ordered by occupation numbers as starting orbitals.
For example, using our NO example, if we want the \(\pi\) orbitals, we just need to specify CAS(5,4).
space=mtp.OrbSpace(molecule,uhf_orbs)
space.cas(5,4)
mcscfdrv=mtp.McscfDriver()
mc_results = mcscfdrv.compute(molecule,basis,space)
Multi-Configurational Self-Consistent Field Driver
====================================================
Active space definition:
------------------------
Number of inactive (occupied) orbitals: 5
Number of active orbitals: 4
Number of virtual orbitals: 19
This is a CASSCF wavefunction: CAS(5,4)
CI expansion:
-------------
Number of determinants: 24
╭────────────────────────────────────╮
│ Driver settings │
╰────────────────────────────────────╯
Number of states : 1
State-specific calculation
- State of interest : 1
Max. iterations : 50
BFGS window : 5
Convergence thresholds:
- Energy change : 1e-08
- Gradient sq. norm : 1e-08
MCSCF Iterations
-------------------
Iter. | Average Energy | E. Change | Grad. Norm | CI Iter. | Trust rad. | Time
---------------------------------------------------------------------------------
-120.73843103721453
1 -129.311457760 0.0e+00 6.2e-03 0 0.40 0:00:00
-120.74715786356627
2 -129.314788870 -3.3e-03 7.5e-04 0 0.40 0:00:00
-120.74275903624346
3 -129.315347735 -5.6e-04 6.3e-05 0 0.48 0:00:00
-120.74246197019936
4 -129.315373471 -2.6e-05 1.8e-06 0 0.58 0:00:00
-120.74232099392256
5 -129.315374340 -8.7e-07 5.0e-08 0 0.58 0:00:00
-120.74245329170331
6 -129.315374368 -2.8e-08 4.3e-09 0 0.58 0:00:00
-120.74248011062907
7 -129.315374371 -2.4e-09 3.3e-10 0 0.69 0:00:00
** Convergence reached in 7 iterations
-120.7424895672272
8 -129.315374371 -3.8e-10 8.3e-11 0 0.80 0:00:00
Final results
-------------
* State 1
- S^2 : 0.75 (multiplicity = 2.0 )
- Energy : -129.31537437094315
- Natural orbitals
1.97487 1.91332 1.02414 0.08767
Spin Restricted Orbitals
------------------------
Molecular Orbital No. 3:
--------------------------
Occupation: 2.000 Energy: -1.58851 a.u.
( 1 N 2s : -0.26) ( 1 N 1p0 : -0.19) ( 2 O 2s : -0.39)
( 2 O 3s : -0.35) ( 2 O 1p0 : 0.17)
Molecular Orbital No. 4:
--------------------------
Occupation: 2.000 Energy: -0.92154 a.u.
( 1 N 2s : -0.36) ( 1 N 3s : -0.34) ( 2 O 2s : 0.30)
( 2 O 3s : 0.43) ( 2 O 1p0 : 0.34) ( 2 O 2p0 : 0.17)
Molecular Orbital No. 5:
--------------------------
Occupation: 2.000 Energy: -0.65920 a.u.
( 1 N 2s : 0.19) ( 1 N 3s : 0.50) ( 1 N 1p0 : -0.46)
( 1 N 2p0 : -0.17) ( 2 O 1p0 : 0.40) ( 2 O 2p0 : 0.25)
Molecular Orbital No. 6:
--------------------------
Occupation: 1.975 Energy: -0.68235 a.u.
( 1 N 1p-1: -0.32) ( 1 N 2p-1: -0.17) ( 2 O 1p+1: -0.16)
( 2 O 1p-1: -0.52) ( 2 O 2p-1: -0.32)
Molecular Orbital No. 7:
--------------------------
Occupation: 1.913 Energy: -0.62471 a.u.
( 1 N 1p+1: -0.35) ( 1 N 2p+1: -0.21) ( 2 O 1p+1: -0.48)
( 2 O 1p-1: 0.15) ( 2 O 2p+1: -0.30)
Molecular Orbital No. 8:
--------------------------
Occupation: 1.024 Energy: -0.11142 a.u.
( 1 N 1p+1: 0.18) ( 1 N 1p-1: 0.58) ( 1 N 2p-1: 0.46)
( 2 O 1p-1: -0.41) ( 2 O 2p-1: -0.36)
Molecular Orbital No. 9:
--------------------------
Occupation: 0.088 Energy: 0.17670 a.u.
( 1 N 1p+1: -0.59) ( 1 N 1p-1: 0.19) ( 1 N 2p+1: -0.33)
( 2 O 1p+1: 0.59) ( 2 O 1p-1: -0.18) ( 2 O 2p+1: 0.19)
Molecular Orbital No. 10:
--------------------------
Occupation: 0.000 Energy: 0.52194 a.u.
( 1 N 3s : 2.05) ( 1 N 1p0 : 0.19) ( 1 N 2p0 : 1.95)
( 2 O 2s : -0.17) ( 2 O 3s : -2.13) ( 2 O 1p0 : 0.22)
( 2 O 2p0 : 1.33)
Molecular Orbital No. 11:
--------------------------
Occupation: 0.000 Energy: 0.92743 a.u.
( 1 N 1p+1: 0.87) ( 1 N 1p-1: -0.27) ( 1 N 2p+1: -1.12)
( 1 N 2p-1: 0.35)
Molecular Orbital No. 12:
--------------------------
Occupation: 0.000 Energy: 0.93060 a.u.
( 1 N 1p+1: -0.28) ( 1 N 1p-1: -0.88) ( 1 N 2p+1: 0.33)
( 1 N 2p-1: 1.06) ( 2 O 1p-1: -0.19)
Dipole moment for state 1
---------------------------
X : 0.000006 a.u. 0.000016 Debye
Y : -0.000002 a.u. -0.000004 Debye
Z : 0.038785 a.u. 0.098580 Debye
Total : 0.038785 a.u. 0.098580 Debye
Total MCSCF time: 00:00:00
However, in some cases, this guess will not be correct and will lead to the wrong orbitals and usually lead to poor MCSCF convergence. Let’s see for example the furan molecule. For some applications, the active space of choice could be the \(\pi\) system of the molecule, or 6 electrons in 5 orbitals. Let’s use the HF canonical orbitals as starting orbitals and use the CAS(6,5) definition.
molecule = vlx.Molecule.from_xyz_string(furan_xyz)
basis = vlx.MolecularBasis.read(molecule,"def2-sv(p)")
scfdrv = vlx.ScfRestrictedDriver()
hf_results = scfdrv.compute(molecule, basis)
space=mtp.OrbSpace(molecule,scfdrv.mol_orbs)
space.cas(6,5)
mcscfdrv=mtp.McscfDriver()
mc_results = mcscfdrv.compute(molecule,basis,space)
Self Consistent Field Driver Setup
====================================
Wave Function Model : Spin-Restricted Hartree-Fock
Initial Guess Model : Superposition of Atomic Densities
Convergence Accelerator : Two Level Direct Inversion of Iterative Subspace
Max. Number of Iterations : 50
Max. Number of Error Vectors : 10
Convergence Threshold : 1.0e-06
ERI Screening Threshold : 1.0e-12
Linear Dependence Threshold : 1.0e-06
* Info * Starting Reduced Basis SCF calculation...
* Info * ...done. SCF energy in reduced basis set: -228.337309732544 a.u. Time: 0.14 sec.
Iter. | Hartree-Fock Energy | Energy Change | Gradient Norm | Max. Gradient | Density Change
--------------------------------------------------------------------------------------------
1 -228.430357327940 0.0000000000 0.16297779 0.01240917 0.00000000
2 -228.433229998540 -0.0028726706 0.05564260 0.00546548 0.06351428
3 -228.433553776746 -0.0003237782 0.00841002 0.00075289 0.01658680
4 -228.433566986233 -0.0000132095 0.00299022 0.00026049 0.00498463
5 -228.433568460938 -0.0000014747 0.00114921 0.00013437 0.00167819
6 -228.433568710663 -0.0000002497 0.00030252 0.00003062 0.00054130
7 -228.433568755137 -0.0000000445 0.00007367 0.00000963 0.00033342
8 -228.433568757939 -0.0000000028 0.00002138 0.00000334 0.00008594
9 -228.433568758030 -0.0000000001 0.00001803 0.00000167 0.00001459
10 -228.433568758068 -0.0000000000 0.00000165 0.00000018 0.00000641
11 -228.433568758069 -0.0000000000 0.00000042 0.00000005 0.00000142
*** SCF converged in 11 iterations. Time: 0.88 sec.
Spin-Restricted Hartree-Fock:
-----------------------------
Total Energy : -228.4335687581 a.u.
Electronic Energy : -387.3224182575 a.u.
Nuclear Repulsion Energy : 158.8888494994 a.u.
------------------------------------
Gradient Norm : 0.0000004225 a.u.
Ground State Information
------------------------
Charge of Molecule : 0.0
Multiplicity (2S+1) : 1
Magnetic Quantum Number (M_S) : 0.0
Multi-Configurational Self-Consistent Field Driver
====================================================
Active space definition:
------------------------
Number of inactive (occupied) orbitals: 15
Number of active orbitals: 5
Number of virtual orbitals: 58
This is a CASSCF wavefunction: CAS(6,5)
CI expansion:
-------------
Number of determinants: 100
╭────────────────────────────────────╮
│ Driver settings │
╰────────────────────────────────────╯
Number of states : 1
State-specific calculation
- State of interest : 1
Max. iterations : 50
BFGS window : 5
Convergence thresholds:
- Energy change : 1e-08
- Gradient sq. norm : 1e-08
MCSCF Iterations
-------------------
Iter. | Average Energy | E. Change | Grad. Norm | CI Iter. | Trust rad. | Time
---------------------------------------------------------------------------------
-221.08098206335433
1 -228.446979427 0.0e+00 4.4e-04 0 0.40 0:00:00
-221.04649405396805
2 -228.450765627 -3.8e-03 2.8e-04 0 0.40 0:00:00
-220.81494856195715
3 -228.453532440 -2.8e-03 2.1e-04 0 0.56 0:00:00
-220.25363314530583
4 -228.457020218 -3.5e-03 6.1e-04 0 0.78 0:00:00
-219.69833942228456
5 -228.456020527 1.0e-03 1.4e-03 0 0.80 0:00:00
-219.76220800660204
6 -228.461486557 -5.5e-03 8.3e-04 0 0.40 0:00:00
-219.61378096976796
7 -228.467331734 -5.8e-03 7.8e-04 0 0.56 0:00:00
-219.19427904264478
8 -228.471988905 -4.7e-03 8.1e-04 0 0.67 0:00:00
-218.95389734228235
9 -228.475632987 -3.6e-03 4.8e-04 0 0.80 0:00:00
-219.04375034567687
10 -228.478452126 -2.8e-03 2.7e-04 0 0.80 0:00:00
-219.0884683415016
11 -228.480503385 -2.1e-03 7.3e-04 0 0.80 0:00:00
-218.9649003902477
12 -228.481646600 -1.1e-03 1.8e-04 0 0.80 0:00:00
-218.72844901106828
13 -228.482624506 -9.8e-04 6.6e-05 0 0.80 0:00:00
-218.5249284343497
14 -228.482976018 -3.5e-04 8.9e-05 0 0.80 0:00:00
-218.55124957465694
15 -228.483107168 -1.3e-04 4.3e-05 0 0.80 0:00:00
-218.56486625660335
16 -228.483230601 -1.2e-04 1.5e-05 0 0.80 0:00:00
-218.56454615804992
17 -228.483364649 -1.3e-04 3.7e-05 0 0.80 0:00:00
-218.49188800941897
18 -228.482614798 7.5e-04 7.0e-03 0 0.80 0:00:00
-218.52984141785004
19 -228.484099389 -1.5e-03 1.2e-03 0 0.40 0:00:00
-218.5359851205917
20 -228.484606673 -5.1e-04 1.5e-03 0 0.40 0:00:00
-218.47789350474702
21 -228.486665255 -2.1e-03 3.9e-04 0 0.56 0:00:00
-218.45317950080158
22 -228.489364387 -2.7e-03 5.9e-04 0 0.78 0:00:00
-218.48941566373094
23 -228.490105357 -7.4e-04 4.6e-04 0 0.80 0:00:00
-218.53602740125694
24 -228.490726224 -6.2e-04 7.8e-05 0 0.80 0:00:00
-218.5362829224382
25 -228.490873747 -1.5e-04 2.2e-05 0 0.80 0:00:00
-218.53304279914053
26 -228.490957029 -8.3e-05 1.3e-05 0 0.80 0:00:00
-218.5359625420554
27 -228.490985956 -2.9e-05 7.2e-06 0 0.80 0:00:00
-218.5406982010631
28 -228.491004941 -1.9e-05 2.7e-06 0 0.80 0:00:00
-218.54278501079574
29 -228.491010718 -5.8e-06 1.3e-06 0 0.80 0:00:00
-218.54458608583047
30 -228.491014542 -3.8e-06 4.2e-07 0 0.80 0:00:00
-218.5460431355291
31 -228.491015873 -1.3e-06 2.4e-07 0 0.80 0:00:00
-218.54728931178005
32 -228.491016349 -4.8e-07 9.7e-08 0 0.80 0:00:00
-218.5479923047976
33 -228.491016466 -1.2e-07 1.8e-08 0 0.80 0:00:00
-218.54820224684417
34 -228.491016487 -2.1e-08 3.3e-09 0 0.80 0:00:00
-218.54823916984265
35 -228.491016491 -3.8e-09 5.9e-10 0 0.80 0:00:00
** Convergence reached in 35 iterations
-218.54820876603407
36 -228.491016492 -9.8e-10 1.4e-10 0 0.80 0:00:00
WARNING! Your active orbitals may have changed significantly
Min_overlap = 0.0553578216270579
Final results
-------------
* State 1
- S^2 : 0.00 (multiplicity = 1.0 )
- Energy : -228.49101649222374
- Natural orbitals
1.98999 1.97270 1.97224 0.03736 0.02771
Spin Restricted Orbitals
------------------------
Molecular Orbital No. 14:
--------------------------
Occupation: 2.000 Energy: -0.54989 a.u.
( 1 C 1p+1: 0.22) ( 2 H 1s : -0.18) ( 3 C 1p+1: -0.17)
( 3 C 1p-1: -0.33) ( 3 C 2p-1: -0.15) ( 4 C 1p-1: 0.38)
( 4 C 2p-1: 0.17) ( 5 C 1p+1: 0.19)
Molecular Orbital No. 15:
--------------------------
Occupation: 2.000 Energy: -0.46672 a.u.
( 3 C 1p0 : -0.17) ( 4 C 1p0 : -0.37) ( 4 C 2p0 : -0.27)
( 5 C 1p0 : -0.36) ( 5 C 2p0 : -0.25)
Molecular Orbital No. 16:
--------------------------
Occupation: 1.990 Energy: -0.34849 a.u.
( 1 C 1p0 : -0.40) ( 1 C 2p0 : -0.30) ( 3 C 1p0 : -0.36)
( 3 C 2p0 : -0.30) ( 5 C 1p0 : 0.23) ( 5 C 2p0 : 0.22)
Molecular Orbital No. 17:
--------------------------
Occupation: 1.973 Energy: -0.78943 a.u.
( 1 C 2s : 0.17) ( 1 C 1p-1: 0.33) ( 6 O 1p+1: 0.21)
( 6 O 1p-1: -0.47) ( 6 O 2p+1: 0.15) ( 6 O 2p-1: -0.27)
Molecular Orbital No. 18:
--------------------------
Occupation: 1.972 Energy: -0.51820 a.u.
( 6 O 1p0 : 0.59) ( 6 O 2p0 : 0.43)
Molecular Orbital No. 19:
--------------------------
Occupation: 0.037 Energy: 0.36068 a.u.
( 1 C 1p0 : 0.47) ( 1 C 2p0 : 0.25) ( 3 C 1p0 : -0.29)
( 3 C 2p0 : -0.28) ( 5 C 1p0 : 0.24) ( 5 C 2p0 : 0.22)
( 6 O 1p0 : -0.66) ( 6 O 2p0 : 0.18)
Molecular Orbital No. 20:
--------------------------
Occupation: 0.028 Energy: 0.71201 a.u.
( 1 C 1s : 0.16) ( 1 C 2s : -0.49) ( 1 C 1p+1: 0.24)
( 1 C 1p-1: -0.62) ( 1 C 2p-1: -0.17) ( 6 O 2s : 0.45)
( 6 O 1p+1: 0.28) ( 6 O 1p-1: -0.73)
Molecular Orbital No. 21:
--------------------------
Occupation: 0.000 Energy: 0.20635 a.u.
( 1 C 3s : 0.74) ( 1 C 2p+1: -0.36) ( 1 C 2p-1: -0.22)
( 2 H 2s : -1.07) ( 3 C 3s : 0.95) ( 3 C 2p+1: 0.25)
( 3 C 2p-1: -0.32) ( 4 C 3s : 0.94) ( 4 C 2p+1: 0.46)
( 4 C 2p-1: 0.19) ( 5 C 3s : 0.72) ( 5 C 2p-1: 0.41)
( 7 H 2s : -1.18) ( 8 H 2s : -1.30) ( 9 H 2s : -1.09)
Molecular Orbital No. 22:
--------------------------
Occupation: 0.000 Energy: 0.22211 a.u.
( 1 C 2p0 : 0.39) ( 3 C 1p0 : -0.28) ( 3 C 2p0 : -0.72)
( 4 C 1p0 : 0.36) ( 4 C 2p0 : 0.89) ( 5 C 1p0 : -0.28)
( 5 C 2p0 : -0.65)
Dipole moment for state 1
---------------------------
X : 0.301536 a.u. 0.766427 Debye
Y : -0.058177 a.u. -0.147871 Debye
Z : 0.000044 a.u. 0.000112 Debye
Total : 0.307097 a.u. 0.780562 Debye
Total MCSCF time: 00:00:10
You see the convergence is very poor. Additionally, we get a final warning that our active space has changed significantly. The reason for it is that we started with the wrong active orbitals and the MCSCF had to make up for this. If we visualize the final result using the orbital viewer, we can see that 2 of our final active orbitals are actually a \(\sigma\) \(\sigma^*\) pair.
If the orbitals are not conveniently situated around the HOMO-LUMO gap, we need a way to visualize them to choose the right ones.
Visualizing the starting orbitals#
Veloxchem and MultiPsi provide a fast and convenient way to visualize the orbitals on a jupyter notebook using the Orbital Viewer. The Veloxchem Orbital Viewer allows to visualize the orbitals, but MultiPsi has its own version which adds several functionalities facilitating the active space selection.
The plot function requires a vlx.Molecule, a vlx.MolecularBasis and either a vlx.MolecularOrbital, a mtp.OrbSpace or a h5 filename.
orbviewer=mtp.OrbitalViewer()
orbviewer.plot(molecule,basis,scfdrv.mol_orbs)
By simply clickling on the “Active” checkbox, one can make an orbital active or inactive. The viewer shows the list of selected active orbitals and conveniently allows to save a h5 file which contains both the orbitals and the active space information that can be use to initialize an OrbSpace object.
Here we find that the orbitals 12, 17, 18, 19 and 21 are the \(\pi\) orbitals we are looking for. There are now 3 equivalent ways to define these orbitals as our active space.
the cas_orbitals function with the list of active orbitals, which simultaneously defines a CAS expansion and selects the active orbitals
the pick_active function with the list of active orbitals after having defined the type of expansion (for example using the CAS(m,n) function)
using the h5 file directly
Note that the orbital lists are python list, and thus start from 0. You need to remove 1 to the orbital index.
# 1) CASOrb function
space = mtp.OrbSpace(molecule, scfdrv.mol_orbs)
space.cas_orbitals([11,16,17,18,20])
mc_results = mcscfdrv.compute(molecule,basis,space)
# 2) activeOrb
space = mtp.OrbSpace(molecule, scfdrv.mol_orbs)
space.cas(6,5)
space.pick_active([11,16,17,18,20])
mc_results = mcscfdrv.compute(molecule,basis,space)
# 3) h5 file generated by OrbitalViewer
space = mtp.OrbSpace(molecule, "furan-cas.h5")
mc_results = mcscfdrv.compute(molecule,basis,space)
Multi-Configurational Self-Consistent Field Driver
====================================================
Active space definition:
------------------------
Number of inactive (occupied) orbitals: 15
Number of active orbitals: 5
Number of virtual orbitals: 58
This is a CASSCF wavefunction: CAS(6,5)
CI expansion:
-------------
Number of determinants: 100
╭────────────────────────────────────╮
│ Driver settings │
╰────────────────────────────────────╯
Number of states : 1
State-specific calculation
- State of interest : 1
Max. iterations : 50
BFGS window : 5
Convergence thresholds:
- Energy change : 1e-08
- Gradient sq. norm : 1e-08
MCSCF Iterations
-------------------
Iter. | Average Energy | E. Change | Grad. Norm | CI Iter. | Trust rad. | Time
---------------------------------------------------------------------------------
-221.02436914119917
1 -228.473029266 0.0e+00 1.7e-03 0 0.40 0:00:00
-221.03262655673217
2 -228.486533870 -1.4e-02 7.8e-04 0 0.40 0:00:00
-221.0104387781901
3 -228.487831306 -1.3e-03 1.2e-04 0 0.48 0:00:00
-221.009282476445
4 -228.488000819 -1.7e-04 1.5e-05 0 0.48 0:00:00
-221.00561918072574
5 -228.488013713 -1.3e-05 4.3e-07 0 0.58 0:00:00
-221.00554171318456
6 -228.488014060 -3.5e-07 5.3e-08 0 0.58 0:00:00
-221.0053940734887
7 -228.488014076 -1.6e-08 1.5e-09 0 0.58 0:00:00
-221.00542610770245
8 -228.488014077 -5.6e-10 6.4e-11 0 0.58 0:00:00
** Convergence reached in 8 iterations
-221.00542892220935
9 -228.488014077 -3.3e-11 8.8e-12 0 0.69 0:00:00
Final results
-------------
* State 1
- S^2 : 0.00 (multiplicity = 1.0 )
- Energy : -228.4880140768445
- Natural orbitals
1.99451 1.94212 1.90287 0.09707 0.06342
Spin Restricted Orbitals
------------------------
Molecular Orbital No. 14:
--------------------------
Occupation: 2.000 Energy: -0.55816 a.u.
( 3 C 1p-1: -0.32) ( 4 C 1p+1: 0.19) ( 4 C 1p-1: 0.26)
( 6 O 1p+1: 0.35) ( 6 O 2p+1: 0.26) ( 7 H 1s : 0.20)
( 8 H 1s : 0.20)
Molecular Orbital No. 15:
--------------------------
Occupation: 2.000 Energy: -0.54203 a.u.
( 1 C 1p+1: 0.25) ( 3 C 1p+1: -0.23) ( 3 C 1p-1: -0.18)
( 4 C 1p-1: 0.27) ( 5 C 1p+1: 0.24) ( 6 O 1p+1: -0.34)
( 6 O 2p+1: -0.27)
Molecular Orbital No. 16:
--------------------------
Occupation: 1.995 Energy: -0.58222 a.u.
( 6 O 1p0 : -0.59) ( 6 O 2p0 : -0.43)
Molecular Orbital No. 17:
--------------------------
Occupation: 1.942 Energy: -0.41896 a.u.
( 1 C 1p0 : 0.19) ( 3 C 1p0 : 0.34) ( 3 C 2p0 : 0.23)
( 4 C 1p0 : 0.34) ( 4 C 2p0 : 0.23) ( 5 C 1p0 : 0.19)
( 6 O 2p0 : -0.15)
Molecular Orbital No. 18:
--------------------------
Occupation: 1.903 Energy: -0.30937 a.u.
( 1 C 1p0 : -0.36) ( 1 C 2p0 : -0.28) ( 3 C 1p0 : -0.24)
( 3 C 2p0 : -0.20) ( 4 C 1p0 : 0.24) ( 4 C 2p0 : 0.20)
( 5 C 1p0 : 0.36) ( 5 C 2p0 : 0.28)
Molecular Orbital No. 19:
--------------------------
Occupation: 0.097 Energy: 0.18407 a.u.
( 1 C 1p0 : 0.48) ( 1 C 2p0 : 0.30) ( 3 C 1p0 : -0.30)
( 3 C 2p0 : -0.18) ( 4 C 1p0 : -0.30) ( 4 C 2p0 : -0.18)
( 5 C 1p0 : 0.48) ( 5 C 2p0 : 0.30) ( 6 O 1p0 : -0.28)
( 6 O 2p0 : -0.27)
Molecular Orbital No. 20:
--------------------------
Occupation: 0.063 Energy: 0.28638 a.u.
( 1 C 1p0 : -0.36) ( 1 C 2p0 : -0.24) ( 3 C 1p0 : 0.51)
( 3 C 2p0 : 0.39) ( 4 C 1p0 : -0.51) ( 4 C 2p0 : -0.39)
( 5 C 1p0 : 0.36) ( 5 C 2p0 : 0.24)
Molecular Orbital No. 21:
--------------------------
Occupation: 0.000 Energy: 0.20569 a.u.
( 1 C 3s : -0.62) ( 1 C 2p+1: 0.28) ( 1 C 2p-1: 0.21)
( 2 H 2s : 0.93) ( 3 C 3s : -0.99) ( 3 C 2p+1: -0.30)
( 3 C 2p-1: 0.41) ( 4 C 3s : -0.99) ( 4 C 2p+1: -0.48)
( 4 C 2p-1: -0.18) ( 5 C 3s : -0.62) ( 5 C 2p-1: -0.33)
( 7 H 2s : 1.35) ( 8 H 2s : 1.35) ( 9 H 2s : 0.93)
Molecular Orbital No. 22:
--------------------------
Occupation: 0.000 Energy: 0.23274 a.u.
( 1 C 3s : 0.75) ( 1 C 2p+1: -0.85) ( 1 C 2p-1: -0.79)
( 2 H 2s : -1.89) ( 3 C 3s : 1.12) ( 3 C 2p-1: 0.22)
( 4 C 3s : -1.13) ( 4 C 2p-1: 0.19) ( 5 C 3s : -0.74)
( 5 C 1p-1: -0.16) ( 5 C 2p+1: 0.28) ( 5 C 2p-1: -1.12)
( 7 H 2s : -0.40) ( 8 H 2s : 0.41) ( 9 H 2s : 1.87)
Dipole moment for state 1
---------------------------
X : 0.435398 a.u. 1.106672 Debye
Y : -0.130480 a.u. -0.331646 Debye
Z : 0.000065 a.u. 0.000165 Debye
Total : 0.454529 a.u. 1.155297 Debye
Total MCSCF time: 00:00:02
Multi-Configurational Self-Consistent Field Driver
====================================================
Active space definition:
------------------------
Number of inactive (occupied) orbitals: 15
Number of active orbitals: 5
Number of virtual orbitals: 58
This is a CASSCF wavefunction: CAS(6,5)
CI expansion:
-------------
Number of determinants: 100
╭────────────────────────────────────╮
│ Driver settings │
╰────────────────────────────────────╯
Number of states : 1
State-specific calculation
- State of interest : 1
Max. iterations : 50
BFGS window : 5
Convergence thresholds:
- Energy change : 1e-08
- Gradient sq. norm : 1e-08
MCSCF Iterations
-------------------
Iter. | Average Energy | E. Change | Grad. Norm | CI Iter. | Trust rad. | Time
---------------------------------------------------------------------------------
-221.02436914119917
1 -228.473029266 0.0e+00 1.7e-03 0 0.40 0:00:00
-221.03262655673217
2 -228.486533870 -1.4e-02 7.8e-04 0 0.40 0:00:00
-221.0104387781904
3 -228.487831306 -1.3e-03 1.2e-04 0 0.48 0:00:00
-221.00928247644535
4 -228.488000819 -1.7e-04 1.5e-05 0 0.48 0:00:00
-221.00561918072597
5 -228.488013713 -1.3e-05 4.3e-07 0 0.58 0:00:00
-221.00554171318504
6 -228.488014060 -3.5e-07 5.3e-08 0 0.58 0:00:00
-221.00539407348896
7 -228.488014076 -1.6e-08 1.5e-09 0 0.58 0:00:00
-221.0054261077027
8 -228.488014077 -5.6e-10 6.4e-11 0 0.58 0:00:00
** Convergence reached in 8 iterations
-221.00542892220952
9 -228.488014077 -3.4e-11 8.8e-12 0 0.69 0:00:00
Final results
-------------
* State 1
- S^2 : -0.00 (multiplicity = 1.0 )
- Energy : -228.48801407684476
- Natural orbitals
1.99451 1.94212 1.90287 0.09707 0.06342
Spin Restricted Orbitals
------------------------
Molecular Orbital No. 14:
--------------------------
Occupation: 2.000 Energy: -0.55816 a.u.
( 3 C 1p-1: -0.32) ( 4 C 1p+1: 0.19) ( 4 C 1p-1: 0.26)
( 6 O 1p+1: 0.35) ( 6 O 2p+1: 0.26) ( 7 H 1s : 0.20)
( 8 H 1s : 0.20)
Molecular Orbital No. 15:
--------------------------
Occupation: 2.000 Energy: -0.54203 a.u.
( 1 C 1p+1: 0.25) ( 3 C 1p+1: -0.23) ( 3 C 1p-1: -0.18)
( 4 C 1p-1: 0.27) ( 5 C 1p+1: 0.24) ( 6 O 1p+1: -0.34)
( 6 O 2p+1: -0.27)
Molecular Orbital No. 16:
--------------------------
Occupation: 1.995 Energy: -0.58222 a.u.
( 6 O 1p0 : -0.59) ( 6 O 2p0 : -0.43)
Molecular Orbital No. 17:
--------------------------
Occupation: 1.942 Energy: -0.41896 a.u.
( 1 C 1p0 : 0.19) ( 3 C 1p0 : 0.34) ( 3 C 2p0 : 0.23)
( 4 C 1p0 : 0.34) ( 4 C 2p0 : 0.23) ( 5 C 1p0 : 0.19)
( 6 O 2p0 : -0.15)
Molecular Orbital No. 18:
--------------------------
Occupation: 1.903 Energy: -0.30937 a.u.
( 1 C 1p0 : -0.36) ( 1 C 2p0 : -0.28) ( 3 C 1p0 : -0.24)
( 3 C 2p0 : -0.20) ( 4 C 1p0 : 0.24) ( 4 C 2p0 : 0.20)
( 5 C 1p0 : 0.36) ( 5 C 2p0 : 0.28)
Molecular Orbital No. 19:
--------------------------
Occupation: 0.097 Energy: 0.18407 a.u.
( 1 C 1p0 : 0.48) ( 1 C 2p0 : 0.30) ( 3 C 1p0 : -0.30)
( 3 C 2p0 : -0.18) ( 4 C 1p0 : -0.30) ( 4 C 2p0 : -0.18)
( 5 C 1p0 : 0.48) ( 5 C 2p0 : 0.30) ( 6 O 1p0 : -0.28)
( 6 O 2p0 : -0.27)
Molecular Orbital No. 20:
--------------------------
Occupation: 0.063 Energy: 0.28638 a.u.
( 1 C 1p0 : -0.36) ( 1 C 2p0 : -0.24) ( 3 C 1p0 : 0.51)
( 3 C 2p0 : 0.39) ( 4 C 1p0 : -0.51) ( 4 C 2p0 : -0.39)
( 5 C 1p0 : 0.36) ( 5 C 2p0 : 0.24)
Molecular Orbital No. 21:
--------------------------
Occupation: 0.000 Energy: 0.20569 a.u.
( 1 C 3s : -0.62) ( 1 C 2p+1: 0.28) ( 1 C 2p-1: 0.21)
( 2 H 2s : 0.93) ( 3 C 3s : -0.99) ( 3 C 2p+1: -0.30)
( 3 C 2p-1: 0.41) ( 4 C 3s : -0.99) ( 4 C 2p+1: -0.48)
( 4 C 2p-1: -0.18) ( 5 C 3s : -0.62) ( 5 C 2p-1: -0.33)
( 7 H 2s : 1.35) ( 8 H 2s : 1.35) ( 9 H 2s : 0.93)
Molecular Orbital No. 22:
--------------------------
Occupation: 0.000 Energy: 0.23274 a.u.
( 1 C 3s : 0.75) ( 1 C 2p+1: -0.85) ( 1 C 2p-1: -0.79)
( 2 H 2s : -1.89) ( 3 C 3s : 1.12) ( 3 C 2p-1: 0.22)
( 4 C 3s : -1.13) ( 4 C 2p-1: 0.19) ( 5 C 3s : -0.74)
( 5 C 1p-1: -0.16) ( 5 C 2p+1: 0.28) ( 5 C 2p-1: -1.12)
( 7 H 2s : -0.40) ( 8 H 2s : 0.41) ( 9 H 2s : 1.87)
Dipole moment for state 1
---------------------------
X : 0.435398 a.u. 1.106672 Debye
Y : -0.130480 a.u. -0.331646 Debye
Z : 0.000065 a.u. 0.000165 Debye
Total : 0.454529 a.u. 1.155297 Debye
Total MCSCF time: 00:00:03
Multi-Configurational Self-Consistent Field Driver
====================================================
Active space definition:
------------------------
Number of inactive (occupied) orbitals: 15
Number of active orbitals: 5
Number of virtual orbitals: 58
This is a CASSCF wavefunction: CAS(6,5)
CI expansion:
-------------
Number of determinants: 100
╭────────────────────────────────────╮
│ Driver settings │
╰────────────────────────────────────╯
Number of states : 1
State-specific calculation
- State of interest : 1
Max. iterations : 50
BFGS window : 5
Convergence thresholds:
- Energy change : 1e-08
- Gradient sq. norm : 1e-08
MCSCF Iterations
-------------------
Iter. | Average Energy | E. Change | Grad. Norm | CI Iter. | Trust rad. | Time
---------------------------------------------------------------------------------
-221.02436914119912
1 -228.473029266 0.0e+00 1.7e-03 0 0.40 0:00:00
-221.0326265567325
2 -228.486533870 -1.4e-02 7.8e-04 0 0.40 0:00:00
-221.0104387781904
3 -228.487831306 -1.3e-03 1.2e-04 0 0.48 0:00:00
-221.0092824764453
4 -228.488000819 -1.7e-04 1.5e-05 0 0.48 0:00:00
-221.00561918072586
5 -228.488013713 -1.3e-05 4.3e-07 0 0.58 0:00:00
-221.005541713185
6 -228.488014060 -3.5e-07 5.3e-08 0 0.58 0:00:00
-221.00539407348896
7 -228.488014076 -1.6e-08 1.5e-09 0 0.58 0:00:00
-221.00542610770287
8 -228.488014077 -5.7e-10 6.4e-11 0 0.58 0:00:00
** Convergence reached in 8 iterations
-221.00542892220972
9 -228.488014077 -3.3e-11 8.8e-12 0 0.69 0:00:00
Final results
-------------
* State 1
- S^2 : -0.00 (multiplicity = 1.0 )
- Energy : -228.48801407684482
- Natural orbitals
1.99451 1.94212 1.90287 0.09707 0.06342
Spin Restricted Orbitals
------------------------
Molecular Orbital No. 14:
--------------------------
Occupation: 2.000 Energy: -0.55816 a.u.
( 3 C 1p-1: -0.32) ( 4 C 1p+1: 0.19) ( 4 C 1p-1: 0.26)
( 6 O 1p+1: 0.35) ( 6 O 2p+1: 0.26) ( 7 H 1s : 0.20)
( 8 H 1s : 0.20)
Molecular Orbital No. 15:
--------------------------
Occupation: 2.000 Energy: -0.54203 a.u.
( 1 C 1p+1: 0.25) ( 3 C 1p+1: -0.23) ( 3 C 1p-1: -0.18)
( 4 C 1p-1: 0.27) ( 5 C 1p+1: 0.24) ( 6 O 1p+1: -0.34)
( 6 O 2p+1: -0.27)
Molecular Orbital No. 16:
--------------------------
Occupation: 1.995 Energy: -0.58222 a.u.
( 6 O 1p0 : 0.59) ( 6 O 2p0 : 0.43)
Molecular Orbital No. 17:
--------------------------
Occupation: 1.942 Energy: -0.41896 a.u.
( 1 C 1p0 : 0.19) ( 3 C 1p0 : 0.34) ( 3 C 2p0 : 0.23)
( 4 C 1p0 : 0.34) ( 4 C 2p0 : 0.23) ( 5 C 1p0 : 0.19)
( 6 O 2p0 : -0.15)
Molecular Orbital No. 18:
--------------------------
Occupation: 1.903 Energy: -0.30937 a.u.
( 1 C 1p0 : 0.36) ( 1 C 2p0 : 0.28) ( 3 C 1p0 : 0.24)
( 3 C 2p0 : 0.20) ( 4 C 1p0 : -0.24) ( 4 C 2p0 : -0.20)
( 5 C 1p0 : -0.36) ( 5 C 2p0 : -0.28)
Molecular Orbital No. 19:
--------------------------
Occupation: 0.097 Energy: 0.18407 a.u.
( 1 C 1p0 : 0.48) ( 1 C 2p0 : 0.30) ( 3 C 1p0 : -0.30)
( 3 C 2p0 : -0.18) ( 4 C 1p0 : -0.30) ( 4 C 2p0 : -0.18)
( 5 C 1p0 : 0.48) ( 5 C 2p0 : 0.30) ( 6 O 1p0 : -0.28)
( 6 O 2p0 : -0.27)
Molecular Orbital No. 20:
--------------------------
Occupation: 0.063 Energy: 0.28638 a.u.
( 1 C 1p0 : -0.36) ( 1 C 2p0 : -0.24) ( 3 C 1p0 : 0.51)
( 3 C 2p0 : 0.39) ( 4 C 1p0 : -0.51) ( 4 C 2p0 : -0.39)
( 5 C 1p0 : 0.36) ( 5 C 2p0 : 0.24)
Molecular Orbital No. 21:
--------------------------
Occupation: 0.000 Energy: 0.20569 a.u.
( 1 C 3s : -0.62) ( 1 C 2p+1: 0.28) ( 1 C 2p-1: 0.21)
( 2 H 2s : 0.93) ( 3 C 3s : -0.99) ( 3 C 2p+1: -0.30)
( 3 C 2p-1: 0.41) ( 4 C 3s : -0.99) ( 4 C 2p+1: -0.48)
( 4 C 2p-1: -0.18) ( 5 C 3s : -0.62) ( 5 C 2p-1: -0.33)
( 7 H 2s : 1.35) ( 8 H 2s : 1.35) ( 9 H 2s : 0.93)
Molecular Orbital No. 22:
--------------------------
Occupation: 0.000 Energy: 0.23274 a.u.
( 1 C 3s : 0.75) ( 1 C 2p+1: -0.85) ( 1 C 2p-1: -0.79)
( 2 H 2s : -1.89) ( 3 C 3s : 1.12) ( 3 C 2p-1: 0.22)
( 4 C 3s : -1.13) ( 4 C 2p-1: 0.19) ( 5 C 3s : -0.74)
( 5 C 1p-1: -0.16) ( 5 C 2p+1: 0.28) ( 5 C 2p-1: -1.12)
( 7 H 2s : -0.40) ( 8 H 2s : 0.41) ( 9 H 2s : 1.87)
Dipole moment for state 1
---------------------------
X : 0.435398 a.u. 1.106672 Debye
Y : -0.130480 a.u. -0.331646 Debye
Z : 0.000065 a.u. 0.000165 Debye
Total : 0.454529 a.u. 1.155297 Debye
Total MCSCF time: 00:00:03
While the OrbitalViewer is probably the most convenient option, it is currently only available in a jupyter notebook which may not be your preferred workflow. Fortunately, there exists alternatives, namely VIAMD and the VisualizationDriver of VeloxChem to produce cube files.
Restricted and generalized active spaces#
MultiPsi supports a variety of CI expansions, including the restricted active space (RAS) or the more general generalized active space (GAS) with an arbitrary number of spaces. RAS expansions have their own keywords:
help(mtp.OrbSpace.ras)
Help on function ras in module multipsi.orbspace:
ras(self, n_electrons, max_ras1_holes, max_ras3_electrons, n_ras1, n_ras2, n_ras3)
Define a restricted active space divided into 3 spaces called ras1, ras2 and ras3
with restriction on the number of holes and electrons in ras1 resp. ras3.
:param n_electrons:
Number of active electrons.
:param max_ras1_holes:
The maximum number of holes allowed in ras1.
:param max_ras3_holes:
The maximum number of electrons allowed in ras3.
:param n_ras1:
The number of orbitals in ras1.
:param n_ras2:
The number of orbitals in ras2.
:param n_ras3:
The number of orbitals in ras3.
However, we can also use the more general GAS keyword. To define a GAS, you define a list of spaces with the cumulated number of orbitals and the cumulated minimal and maximal number of electrons allowed in that space. For example, in our furan example, we can simply define a CISD expansion in our active space by defining 2 spaces, the first one, corresponding to orbitals occupied in the dominant configuration, contains 3 orbitals and a maximum of 6 electrons. As we want up to 2 excitations, the number of electrons in that first space can go down to 4. The first space is thus defined with the list [3,4,6].
The second space, corresponding to the unoccupied orbitals in the dominant configuration, contains the remaining orbitals and electrons. Since we provide the cumulated number, this space is defined with the list [5,6,6], corresponding to our total 6 electrons in 5 orbitals. Note that the last space always has the same number of minimum and maximum cumulated electrons.
The final keyword is thus:
space.gas([[3,4,6],[5,6,6]])
mc_results = mcscfdrv.compute(molecule,basis,space)
Multi-Configurational Self-Consistent Field Driver
====================================================
Active space definition:
------------------------
Number of inactive (occupied) orbitals: 15
Number of active orbitals: 5
Number of virtual orbitals: 58
This is a GASSCF wavefunction
Cumulated Min cumulated Max cumulated
Space orbitals occupation occupation
1 3 4 6
2 5 6 6
CI expansion:
-------------
Number of determinants: 55
╭────────────────────────────────────╮
│ Driver settings │
╰────────────────────────────────────╯
Number of states : 1
State-specific calculation
- State of interest : 1
Max. iterations : 50
BFGS window : 5
Convergence thresholds:
- Energy change : 1e-08
- Gradient sq. norm : 1e-08
Include Active-active rotations
MCSCF Iterations
-------------------
Iter. | Average Energy | E. Change | Grad. Norm | CI Iter. | Trust rad. | Time
---------------------------------------------------------------------------------
-221.00542892220938
1 -228.486019249 0.0e+00 2.3e-05 0 0.40 0:00:00
-221.00654021439288
2 -228.486040306 -2.1e-05 5.8e-06 0 0.40 0:00:00
-221.0084290554948
3 -228.486055586 -1.5e-05 1.1e-05 0 0.48 0:00:00
-221.00987572069678
4 -228.486070209 -1.5e-05 3.0e-05 0 0.48 0:00:00
-221.011312603487
5 -228.486111030 -4.1e-05 6.0e-05 0 0.48 0:00:00
-221.01003518590122
6 -228.486123718 -1.3e-05 3.2e-05 0 0.48 0:00:00
-221.00818561803908
7 -228.486034478 8.9e-05 8.4e-06 0 0.34 0:00:00
-221.00622605882177
8 -228.486045295 -1.1e-05 1.2e-05 0 0.17 0:00:00
-221.0051065176676
9 -228.486043913 1.4e-06 1.8e-05 0 0.17 0:00:00
-220.99975154829534
10 -228.486114329 -7.0e-05 9.8e-05 0 0.08 0:00:00
-221.0019065441725
11 -228.485969225 1.5e-04 6.1e-05 0 0.08 0:00:00
-221.00869735790138
12 -228.486069585 -1.0e-04 1.1e-04 0 0.04 0:00:00
-221.00851152398207
13 -228.486040839 2.9e-05 2.5e-05 0 0.04 0:00:00
-221.00759144393243
14 -228.486140739 -1.0e-04 4.0e-05 0 0.04 0:00:00
-221.00605359010058
15 -228.486133693 7.0e-06 5.4e-05 0 0.04 0:00:00
-221.00551986738841
16 -228.486122872 1.1e-05 7.1e-05 0 0.04 0:00:00
-221.00512279463715
17 -228.486113986 8.9e-06 8.6e-05 0 0.04 0:00:00
-221.00481124953617
18 -228.486106314 7.7e-06 9.8e-05 0 0.04 0:00:00
-221.00453751282046
19 -228.486099540 6.8e-06 1.1e-04 0 0.04 0:00:00
-221.00428597883393
20 -228.486093498 6.0e-06 1.2e-04 0 0.04 0:00:00
-221.00408815308114
21 -228.486088057 5.4e-06 1.3e-04 0 0.04 0:00:00
-221.0039089875504
22 -228.486083335 4.7e-06 1.4e-04 0 0.04 0:00:00
-221.0037606147604
23 -228.486079266 4.1e-06 1.4e-04 0 0.04 0:00:00
-221.0036316006972
24 -228.486075785 3.5e-06 1.5e-04 0 0.04 0:00:00
-221.00351993652592
25 -228.486072776 3.0e-06 1.5e-04 0 0.04 0:00:00
-221.00343264170542
26 -228.486069988 2.8e-06 1.6e-04 0 0.04 0:00:00
-221.00335453826952
27 -228.486067379 2.6e-06 1.6e-04 0 0.04 0:00:00
-221.00328707972523
28 -228.486064907 2.5e-06 1.7e-04 0 0.04 0:00:00
-221.0032281655846
29 -228.486062545 2.4e-06 1.7e-04 0 0.04 0:00:00
-221.00317829076408
30 -228.486060271 2.3e-06 1.7e-04 0 0.04 0:00:00
-221.0031351788487
31 -228.486058073 2.2e-06 1.8e-04 0 0.04 0:00:00
-221.00309707207842
32 -228.486055938 2.1e-06 1.8e-04 0 0.04 0:00:00
-221.00306355114225
33 -228.486053859 2.1e-06 1.8e-04 0 0.04 0:00:00
-221.00303375304847
34 -228.486051832 2.0e-06 1.9e-04 0 0.04 0:00:00
-221.00300720276857
35 -228.486049853 2.0e-06 1.9e-04 0 0.04 0:00:00
-221.0029833807892
36 -228.486047920 1.9e-06 1.9e-04 0 0.04 0:00:00
-221.00296193588605
37 -228.486046032 1.9e-06 1.9e-04 0 0.04 0:00:00
-221.00294266258675
38 -228.486044190 1.8e-06 2.0e-04 0 0.04 0:00:00
-221.00292534270451
39 -228.486042393 1.8e-06 2.0e-04 0 0.04 0:00:00
-221.0029098370974
40 -228.486040643 1.7e-06 2.0e-04 0 0.04 0:00:00
-221.0028960152451
41 -228.486038941 1.7e-06 2.1e-04 0 0.04 0:00:00
-221.00288378053185
42 -228.486037287 1.7e-06 2.1e-04 0 0.04 0:00:00
-221.00287305335712
43 -228.486035682 1.6e-06 2.1e-04 0 0.04 0:00:00
-221.00286375790765
44 -228.486034128 1.6e-06 2.1e-04 0 0.04 0:00:00
-221.0028558291671
45 -228.486032625 1.5e-06 2.1e-04 0 0.04 0:00:00
-221.00284920326982
46 -228.486031173 1.5e-06 2.2e-04 0 0.04 0:00:00
-221.00284381988914
47 -228.486029774 1.4e-06 2.2e-04 0 0.04 0:00:00
-221.00283961899103
48 -228.486028428 1.3e-06 2.2e-04 0 0.04 0:00:00
-221.00283654000043
49 -228.486027134 1.3e-06 2.2e-04 0 0.04 0:00:00
-221.0028345220985
50 -228.486025893 1.2e-06 2.2e-04 0 0.04 0:00:00
** Convergence not reached after 50 iterations
Final results
-------------
* State 1
- S^2 : -0.00 (multiplicity = 1.0 )
- Energy : -228.48602589256444
- Natural orbitals
1.99472 1.94830 1.90957 0.09016 0.05725
Spin Restricted Orbitals
------------------------
Molecular Orbital No. 14:
--------------------------
Occupation: 2.000 Energy: -0.55793 a.u.
( 3 C 1p-1: -0.32) ( 4 C 1p+1: 0.19) ( 4 C 1p-1: 0.26)
( 6 O 1p+1: 0.35) ( 6 O 2p+1: 0.26) ( 7 H 1s : 0.20)
( 8 H 1s : 0.20)
Molecular Orbital No. 15:
--------------------------
Occupation: 2.000 Energy: -0.54198 a.u.
( 1 C 1p+1: -0.25) ( 3 C 1p+1: 0.23) ( 3 C 1p-1: 0.18)
( 4 C 1p-1: -0.28) ( 5 C 1p+1: -0.24) ( 6 O 1p+1: 0.34)
( 6 O 2p+1: 0.26)
Molecular Orbital No. 16:
--------------------------
Occupation: 1.995 Energy: -0.58231 a.u.
( 6 O 1p0 : 0.59) ( 6 O 2p0 : 0.43)
Molecular Orbital No. 17:
--------------------------
Occupation: 1.948 Energy: -0.41965 a.u.
( 1 C 1p0 : 0.18) ( 3 C 1p0 : 0.35) ( 3 C 2p0 : 0.24)
( 4 C 1p0 : 0.35) ( 4 C 2p0 : 0.24) ( 5 C 1p0 : 0.18)
Molecular Orbital No. 18:
--------------------------
Occupation: 1.910 Energy: -0.31015 a.u.
( 1 C 1p0 : 0.37) ( 1 C 2p0 : 0.28) ( 3 C 1p0 : 0.23)
( 3 C 2p0 : 0.19) ( 4 C 1p0 : -0.23) ( 4 C 2p0 : -0.19)
( 5 C 1p0 : -0.37) ( 5 C 2p0 : -0.28)
Molecular Orbital No. 19:
--------------------------
Occupation: 0.090 Energy: 0.18247 a.u.
( 1 C 1p0 : 0.48) ( 1 C 2p0 : 0.31) ( 3 C 1p0 : -0.29)
( 3 C 2p0 : -0.18) ( 4 C 1p0 : -0.29) ( 4 C 2p0 : -0.18)
( 5 C 1p0 : 0.48) ( 5 C 2p0 : 0.31) ( 6 O 1p0 : -0.28)
( 6 O 2p0 : -0.27)
Molecular Orbital No. 20:
--------------------------
Occupation: 0.057 Energy: 0.28369 a.u.
( 1 C 1p0 : -0.35) ( 1 C 2p0 : -0.24) ( 3 C 1p0 : 0.52)
( 3 C 2p0 : 0.41) ( 4 C 1p0 : -0.52) ( 4 C 2p0 : -0.41)
( 5 C 1p0 : 0.35) ( 5 C 2p0 : 0.25)
Molecular Orbital No. 21:
--------------------------
Occupation: 0.000 Energy: 0.20583 a.u.
( 1 C 3s : -0.63) ( 1 C 2p+1: 0.29) ( 1 C 2p-1: 0.22)
( 2 H 2s : 0.94) ( 3 C 3s : -0.99) ( 3 C 2p+1: -0.30)
( 3 C 2p-1: 0.41) ( 4 C 3s : -0.99) ( 4 C 2p+1: -0.47)
( 4 C 2p-1: -0.18) ( 5 C 3s : -0.63) ( 5 C 2p-1: -0.34)
( 7 H 2s : 1.34) ( 8 H 2s : 1.34) ( 9 H 2s : 0.94)
Molecular Orbital No. 22:
--------------------------
Occupation: 0.000 Energy: 0.23262 a.u.
( 1 C 3s : -0.74) ( 1 C 2p+1: 0.86) ( 1 C 2p-1: 0.79)
( 2 H 2s : 1.89) ( 3 C 3s : -1.12) ( 3 C 2p-1: -0.22)
( 4 C 3s : 1.12) ( 4 C 2p-1: -0.20) ( 5 C 3s : 0.74)
( 5 C 1p-1: 0.16) ( 5 C 2p+1: -0.28) ( 5 C 2p-1: 1.12)
( 7 H 2s : 0.38) ( 8 H 2s : -0.39) ( 9 H 2s : -1.88)
Dipole moment for state 1
---------------------------
X : 0.421492 a.u. 1.071326 Debye
Y : -0.126325 a.u. -0.321087 Debye
Z : 0.000056 a.u. 0.000142 Debye
Total : 0.440015 a.u. 1.118408 Debye
Total MCSCF time: 00:00:15
In a RAS or GAS case, the selection of active orbitals is done by the activeOrb keyword after having defined the RAS or GAS expansion. Alternatively, it is possible to use the h5 files from the orbital viewer to initialize the OrbSpace object, and then define the RAS or GAS expansion. In this case, you have to be careful that the number of active orbitals match exactly, and you also need to be careful of the order in which you select the orbitals in the orbital viewer, as the viewer creates a list in the order you select them (even though the active orbital list shown in the viewer is sorted).
RASSCF and GASSCF are tricky calculations to set-up, and are notoriously a little more difficult to converge, but for large active spaces, they can reduce significantly the number of determinants and thus the computational cost.
Saving the active space#
The OrbSpace object is by default updated at every Mcscf calculation, allowing to easily reuse in a script. However, one may in some cases want to create a copy to prevent our OrbSpace from being erased. This can easily be done using the copy() function:
space_backup = space.copy()
The active space can also be saved in a h5 file for later use.
space.write_hdf5("my-casscf.h5")