Tutorial: Setting up the orbitals for a CAS calculation

This tutorial will describe how to set up the orbitals for a CASSCF calculation, with optional NEVPT2 or MRCI steps to follow.

Note that with the ORCA4 release, a CASSCF tutorial was released that can be found in the Downloads section of the Forum.

This tutorial uses quasi-restricted orbitals (QROs) from a DFT calculation as starting orbitals. They are not the only choice but are convenient starting orbitals and typically lead to fast convergence. ROHF orbitals, localized DFT orbitals, MP2 natural orbitals, previous CASSCF-optimized natural orbitals etc. are other choices (sometimes better).

This procedure requires 4 different ORCA job steps. While some of these steps can be combined by $new_job keywords, the orbitals need to be carefully looked at and selected before starting the CASSCF calculation.

1. Creating initial QRO orbitals from a DFT calculation

First do a single-point DFT calculation on your molecule of interest and then ask ORCA to create quasi-restricted orbitals (QRO) using the UNO keyword (creates both unrestricted natural orbitals (UNOs) and quasi-restricted orbitals (QROs).

job1.inp:

! BP86 def2-SVP def2/J UNO Normalprint

After job completion, the QRO orbitals are found in the file job1.qro.

2.A Inspect QRO orbitals

Next you want to take a look at the QRO orbitals, see their composition or even visualize a surface plot. Keyword MOREAD and %moinp is used to read in the previously calculated QRO orbitals and the NoIter keyword is used to prevent the orbitals from being changed by SCF iterations.

job2.inp:

! BP86 def2-SVP def2/J Normalprint NoIter MOREAD
%moinp "job1.qro"

If visualization of the orbitals in question is required (often the case) then one can use orca_plot to create Cube Files for visualization (see more in Visualization and Printing ) or alternatively, one can print the MO and basis coefficients in the ORCA output file so that a program like Chemcraft can render the orbitals on command:

%output
Print[ P_Basis ] 2
Print[ P_MOs ] 1
end

2.B Rotate QRO orbitals if required

Inspect the orbitals given in the Löwdin reduced orbital analysis or look at the orbitals themselves using a GUI program. The active occupied orbitals that you want to select for the CASSCF active space calculation must be ordered so that they appear last in the list of occupied orbitals and unoccupied orbitals that you want for your CAS space should come after them. If your active orbitals are not what you want, you need to rotate them. In the example below we change the position of orbitals 10 and 20 (as orbital 10 is an orbital that we want in the active space).

job2.inp:

! BP86 def2-SVP def2/J Normalprint NoIter MOREAD
%moinp "job1.qro"
%scf rotate{10,20,90} end end

3. A. Running the CASSCF calculation using rotated DFT-QRO orbitals

Now that we are happy with the starting orbitals and their ordering we can go on to the CASSCF calculation. Again we read in orbitals using MOREAD. Here it should be noted that the basis set used for the CASSCF calculation should ideally be the same as that used in the DFT calculation. Otherwise, ORCA will do a basis set projection step when reading in the DFT-QRO orbitals and the orbitals may end up changing.

See CASSCF calculations on how to set up the CASSCF module. Here setting up a simple CASSCF calculation for a 4-electrons-in-3-orbitals CAS(4,3) calculation.

job3.inp:

! def2-SVP def2-SVP/C Normalprint MOREAD
%moinp "job2.gbw"
%casscf
trafostep ri
nel 4
norb 3
mult 1
end

3. B. Running a CASSCF+NEVPT2 calculation

It is also possible in the same input file to request a dynamic correlation energy calculation in the form of NEVPT2 to be carried out after the CASSCF calculation, by simply adding the NEVPT2 keyword.

job3.inp:

! def2-SVP def2-SVP/C Normalprint MOREAD RI-NEVPT2

%moinp "job2.gbw"

%casscf
trafostep ri
nel 4
norb 3
mult 1
end

4. Optional MRCI calculation:

After checking that the CASSCF calculation ran smoothly one could go on to do a MRCI job, which does an even better (and much more expensive!) job of calculating the dynamical correlation energy.

job4.inp:

! MRCI+Q def2-SVP Normalprint MOREAD Allowrhf NoIter
%moinp "job3.gbw"
%casscf
trafostep ri
nel 4
norb 3
mult 1
end