# 3.2.1. Electronic Energy at Fixed Nuclear Geometry¶

The Molcas suite of Quantum Chemical programs is modular in design, and a desired calculation is achieved by executing a list of Molcas program modules in succession, occasionally manipulating the program information files. If the information files from a previous calculation are saved, then a subsequent calculation need not recompute them. This is dependent on the correct information being preserved in the information files for the subsequent calculations. Each module has keywords to specify the functions to be carried out, and many modules rely on the specification of keywords in previous modules.

In the present examples the calculations will be designed by preparing a single file in which the input for the different programs is presented sequentially. The initial problem will be to compute an electronic energy at a fixed geometry of the nuclei, and this will be performed using different methods and thus requiring different Molcas program modules.

First, the proper Molcas environment has to be set up which requires that following variables must be properly defined, for instance:

export MOLCAS=/home/molcas/molcas
export Project=CH4
export WorkDir=/home/user/tmp


If not defined, Molcas provides default values for the above environment variables:

• The MOLCAS variable will be set to the latest implemented version of the code.

This variable is set directly in the Molcas home directory

• Project and WorkDir have the default values None and $PWD, respectively. It is very important that the molcas driver, called by command molcas, and built during the installation of the code, is included in the$PATH.

The first run involves a calculation of the SCF energy of the methane ($$\ce{CH4}$$) molecule. Three programs should be used: GATEWAY to specify information about the system, SEWARD to compute and store the one- and two-electron integrals, and SCF to obtain the Hartree–Fock SCF wave function and energy.

The three Molcas programs to be used leads to three major entries in the input file: GATEWAY, SEWARD, and SCF. The GATEWAY program contains the nuclear geometry in cartesian coordinates and the label for the one-electron basis set. The keyword coord allows automatic insertion of GATEWAY input from a standard file containing the cartesian coordinates in Angstrom which can be generated by programs like LUSCUS or MOLDEN). No symmetry is being considered so the keyword group=C1 is used to force the program not to look for symmetry in the $$\ce{CH4}$$ molecule, and ,thus, input for SEWARD is not required. In closed-shell cases, like $$\ce{CH4}$$, input for SCF is not required. All the input files discussed here can be found at $MOLCAS/doc/samples/problem_based_tutorials, including the file SCF.energy.CH4 described below. *SCF energy for CH4 at a fixed nuclear geometry. *File: SCF.energy.CH4 * &GATEWAY Title = CH4 molecule coord = CH4.xyz basis = STO-3G group = C1 &SEWARD &SCF  where the content of the CH4.xyz file is: 5 distorted CH4 coordinates in Angstroms C 0.000000 0.000000 0.000000 H 0.000000 0.000000 1.050000 H 1.037090 0.000000 -0.366667 H -0.542115 -0.938971 -0.383333 H -0.565685 0.979796 -0.400000  To run Molcas, simply execute the command molcas SCF.energy.CH4.input > SCF.energy.CH4.log 2 > SCF.energy.CH4.err  where the main output is stored in file SCF.energy.CH4.log or molcas -f SCF.energy.CH4.input  where the main output is stored in SCF.energy.CH4.log, and the default error file in SCF.energy.CH4.err. The most relevant information is contained in the output file, where the GATEWAY program information describing the nuclear geometry, molecular symmetry, and the data regarding the one-electron basis sets and the calculation of one- and two-electron integrals, as described in Section 3.3.4. Next, comes the output of program SCF with information of the electronic energy, wave function, and the Hartree–Fock (HF) molecular orbitals (see Section 3.3.5). Files containing intermediate information, integrals, orbitals, etc, will be kept in the$WorkDir directory for further use. For instance, files $Project.OneInt and$Project.OrdInt contain the one- and two-electron integrals stored in binary format. File $Project.ScfOrb stores the HF molecular orbitals in ASCII format, and$Project.RunFile is a communication file between programs. All these files can be used later for more advanced calculations avoiding a repeat of certain calculations.

There are graphical utilities that can be used for the analysis of the results. By default, Molcas generates files which can be read with the MOLDEN program and are found in the $WorkDir including the file CH4.scf.molden. This file contains information about molecular geometry and molecular orbitals, and requires the use if Density Mode in MOLDEN. However, Molcas has its own graphical tool, program LUSCUS, which is a viewer based on openGL and allows the visualization of molecular geometries, orbitals, densities, and density differences. For example, a graphical display of the $$\ce{CH4}$$ molecule can be obtained from a standard coordinate file by the following command: luscus CH4.xyz  In order to obtain the information for displaying molecular orbitals and densities, it is necessary to run the Molcas program called GRID_IT: *SCF energy for CH4 at a fixed nuclear geometry plus a grid for visualization. *File: SCF.energy_grid.CH4 * &GATEWAY Title = CH4 molecule coord = CH4.xyz basis = STO-3G Group = C1 &SEWARD; &SCF &GRID_IT All  Now, execute the Molcas program: molcas SCF.energy_grid.CH4.input -f  In the$WorkDir and $PWD directories a new file is generated, CH4.lus which contains the information required by the GRID_IT input. The file can be visualized by LUSCUS (Open source program, which can be downloaded and installed to your Linux, Windows, or MacOS workstation or laptop). By typing the command: luscus CH4.lus  a window will be opened displaying the molecule and its charge density. By proper selection of options with the mouse buttons, the shape and size of several molecular orbitals can be visualized. GRID_IT can also be run separately, if an orbital file is specified in the input, and the$WorkDir directory is available.

As an alternative to running a specific project, the short script provided below can be placed in the directory $MOLCAS/doc/samples/problem_based_tutorials with the name project.sh. Simply execute the shell script, project.sh$Project, where $Project is the Molcas input, and output files, error files, and a$WorkDir directory called $Project.work will be obtained. #!/bin/bash export MOLCAS=$PWD
export MOLCAS_DISK=2000
export MOLCAS_MEM=64
export MOLCAS_PRINT=3

export Project=$1 export HomeDir=$MOLCAS/doc/samples/problem_based_tutorials
export WorkDir=$HomeDir/$Project.work
mkdir $WorkDir 2>/dev/null molcas$HomeDir/$1 >$HomeDir/$Project.log 2>$HomeDir/$Project.err exit  In order to run a Kohn–Sham density functional calculation, Molcas uses the same SCF module, and, therefore, the only change needed are the specification of the DFT option and required functional (e.g. B3LYP) in the SCF input: *DFT energy for CH4 at a fixed nuclear geometry plus a grid for visualization. *File: DFT.energy.CH4 * &GATEWAY Title = CH4 molecule coord = CH4.xyz basis = STO-3G group = C1 &SEWARD &SCF KSDFT = B3LYP &GRID_IT All  Similar graphical files can be found in$WorkDir and $PWD. The next step is to obtain the second-order Møller–Plesset perturbation (MP2) energy for methane at the same molecular geometry using the same one-electron basis set. Program MBPT2 is now used, and it is possible to take advantage of having previously computed the proper integrals with SEWARD and the reference closed-shell HF wave function with the SCF program. In such cases, it is possible to keep the same definitions as before and simply prepare a file containing the MBPT2 input and run it using the molcas command. The proper intermediate file will be already in$WorkDir. On the other hand, one has to start from scratch, all required inputs should be placed sequentially in the MP2.energy.CH4 file. If the decision is to start the project from the beginning, it is probably a good idea to remove the entire $WorkDir directory, unless it is known for certain the exact nature of the files contained in this directory. *MP2 energy for CH4 at a fixed nuclear geometry. *File: MP2.energy.CH4 * &GATEWAY Title = CH4 molecule coord = CH4.xyz basis = STO-3G group = C1 &SEWARD &SCF &MBPT2 Frozen = 1  In addition to the calculation of a HF wave function, an MP2 calculation has been performed with a frozen deepest orbital, the carbon 1s, of $$\ce{CH4}$$. Information about the output of the MBPT2 program can be found on Section 3.3.6. The SCF program works by default with closed-shell systems with an even number of electrons at the Restricted Hartee–Fock (RHF) level. If, instead there is a need to use the Unrestricted Hartree–Fock (UHF) method, this can be schieved by invoking the keyword UHF. This is possible for both even and odd electron systems. For instance, in a system with an odd number of electrons such as the $$\ce{CH3}$$ radical, with the following Cartesian coordinates 4 CH3 coordinates in Angstrom C 0.000000 0.000000 0.000000 H 0.000000 0.000000 1.050000 H 1.037090 0.000000 -0.366667 H -0.542115 -0.938971 -0.383333  the input to run an open-shell UHF calculation is easily obtained: *SCF/UHF energy for CH3 at a fixed nuclear geometry *File: SCF.energy_UHF.CH3 * &GATEWAY Title = CH3 molecule coord = CH3.xyz basis = STO-3G group = C1 &SEWARD &SCF UHF  If the system is charged, this must be indicated in the SCF input, for example, by computing the cation of the $$\ce{CH4}$$ molecule at the UHF level: *SCF/UHF energy for CH4+ at a fixed nuclear geometry *File: SCF.energy_UHF.CH4plus * &GATEWAY Title = CH4+ molecule coord = CH4.xyz basis = STO-3G group = c1 &SEWARD &SCF UHF Charge = +1  The Kohn–Sham DFT calculation can be also run using the UHF algorithm: *DFT/UHF energy for CH4+ at a fixed nuclear geometry *File: DFT.energy.CH4plus * &GATEWAY Title = CH4+ molecule coord = CH4.xyz basis = STO-3G group = C1 &SEWARD &SCF KSDFT = B3LYP UHF Charge = +1  For the UHF and UHF/DFT methods it is also possible to specify $$\alpha$$ and $$\beta$$ orbital occupations in two ways. 1. First, the keyword ZSPIn can be invoked in the SCF program, which represents the difference between the number of $$\alpha$$ and $$\beta$$ electrons. For example, setting the keyword to 2 forces the program to converge to a result with two more $$\alpha$$ than $$\beta$$ electrons. *DFT/UHF energy for different electronic occupation in CH4 at a fixed nuclear geometry *File: DFT.energy_zspin.CH4 * &GATEWAY Title = CH4 molecule coord = CH4.xyz basis = STO-3G group = c1 &SEWARD &SCF Title = CH4 molecule zspin 2 UHF; ZSPIN = 2 KSDFT = B3LYP  The final occupations in the output will show six $$\alpha$$ and four $$\beta$$ orbitals. 2. Alternatively, instead of ZSPIn, it is possible to specify occupation numbers with keyword Occupation at the beginning of the SCF calculation. This requires an additional input line containing the occupied $$\alpha$$ orbitals (e.g. 6 in this case), and a second line with the $$\beta$$ orbitals (e.g. 4 in this case). Sometimes, SCF convergence may be improved by using this option. Different sets of methods use other Molcas modules. For example, to perform a Complete Active Space (CAS) SCF calculation, the RASSCF program has to be used. This module requires starting trial orbitals, which can be obtained from a previous SCF calculation or, automatically, from the SEWARD program which provides trial orbitals by using a model Fock operator. Recommended keywords are • Nactel defines the total number of active electrons, holes in Ras1, and particles in Ras3, respectively. The last two values are only for RASSCF-type calculations. • Inactive indicates the number of inactive orbitals where the occupation is always 2 in the CASSCF reference, and • Ras2 defines the number of active orbitals. By default, the wave function for the lowest state corresponds to the symmetry with spin multiplicity of 1. Most of the input may not be necessary, if one has prepared and linked an INPORB file with the different orbital types defined by a program like LUSCUS. *CASSCF energy for CH4 at a fixed nuclear geometry *File: CASSCF.energy.CH4 * &GATEWAY coord = CH4.xyz basis = STO-3G group = C1 &SEWARD &RASSCF Title = CH4 molecule Spin = 1; Nactel = 8 0 0; Inactive = 1; Ras2 = 8 &GRID_IT All  In this case, the lowest singlet state (i.e. the ground dstate) is computed, since this is a closed-shell situation with an active space of eight electrons in eight orbitals and with an inactive C 1s orbital, the lowest orbital of the $$CH4$$ molecule. This is a CASSCF example in which all the valence orbitals and electrons (C 2s, C 2p and 4 × H 1s) are included in the active space and allows complete dissociation into atoms. If this is not the goal, then the three almost degenerate highest energy occupied orbitals and the corresponding antibonding unoccupied orbitalsmust be active, leading to a 6 in 6 active space. Using the CASSCF wave function as a reference, it is possible to perform a second-order perturbative, CASPT2, correction to the electronic energy by employing the CASPT2 program. If all previously calculated files are retained in the$WorkDir directory, in particular, integral files (CH4.OneInt, CH4.OrdInt), the CASSCF wave function information file (CH4.JobIph), and communication file CH4.RunFile), it will not be necessary to re-run programs SEWARD, and RASSCF. In this case case, it is enough to prepare a file containing input only for the CASPT2 program followed be execution. Here, however, for the sake of completness, input to all Molcas modules is provided:

*CASPT2 energy for CH4 at a fixed nuclear geometry
*File: CASPT2.energy.CH4
*
&GATEWAY
coord = CH4.xyz; basis = STO-3G; group = C1
&SEWARD
&RASSCF
LumOrb
Spin = 1; Nactel = 8 0 0; Inactive = 1; Ras2 = 8
&CASPT2
Multistate = 1 1


In most of cases, the Hartree–Fock orbitals will be a better choice as starting orbitals. In that case, the RASSCF input has to include keyword LumOrb to read from any external source of orbitals other than those generated by the SEWARD program. By modifying input to the SCF program, it is possible to generate alternative trial orbitals for the RASSCF program. Since a new set of trial orbitals is used, the input to the RASSCF program is also changed. Now, the number of active orbitals, as well as the number of active electrons, are 6.

The two lowest orbitals (Inactive 2) are excluded from the active space and one other orbital is placed in the secondary space. If the previous (8,8) full valence space was used, the CASPT2 program would not be able to include more electronic correlation energy, considering that the calculation involves a minimal basis set. The input for the CASPT2 program includes a frozen C 1s orbital, the lowest orbital in the $$\ce{CH4}$$ molecule.

The charge and multiplicity of our wave function can be changed by computing the $$\ce{CH4^+}$$ cation with the same methods. The RASSCF program defines the character of the problem by specifying the number of electrons, the spin multiplicity, and the spatial symmetry. In the example below, there is one less electron giving rise to doublet multiplicity:

*CASSCF energy for CH4+ at a fixed nuclear geometry
*File: CASSCF.energy.CH4plus
*
&GATEWAY
Title = CH4+ molecule
coord = CH4.xyz; basis = STO-3G; Group = C1
&SEWARD
&RASSCF
Spin = 2; Nactel = 7 0 0; Inactive = 1; Ras2 = 8


No further modification is needed to the CASPT2 input:

*CASPT2 energy for CH4+ at a fixed nuclear geometry
*File: CASPT2.energy.CH4plus
*
&GATEWAY
coord = CH4.xyz; basis = STO-3G; group = C1
&SEWARD
&RASSCF
Title = CH4+ molecule
Spin = 2; Nactel = 1 0 0; Inactive = 4; Ras2 = 1
&CASPT2


A somewhat more sophisticated calculation can be performed at the Restricted Active Space (RAS) SCF level. In such a situation, the level of excitation in the CI expansion can be controlled by restricting the number of holes and particles present in certain orbitals.

*RASSCF energy for CH4 at a fixed nuclear geometry
*File: RASSCF.energy.CH4
*
&GATEWAY
coord = CH4.xyz; basis = STO-3G; group = C1
&SEWARD
&RASSCF
Title = CH4 molecule
Spin = 1; Nactel = 8 1 1
Inactive = 1; Ras1 = 1; Ras2 = 6; Ras3 = 1


In particular, the previous calculation includes one orbital within the Ras1 space and one orbital within the Ras3 space. One hole (single excitation) at maximum is allowed from Ras1 to Ras2 or Ras3, while a maximum of one particle is allowed in Ras3, derived from either Ras1 or Ras2. Within Ras2, all types of orbital occupations are allowed. The RASSCF wave functions can be used as reference for multiconfigurational perturbation theory (RASPT2), but this approach has not been as extensively tested as CASPT2, and, so experience is still somewhat limited.

Molcas also has the possibility of computing electronic energies at different CI levels by using the MRCI program. The input provided below involves a Singles and Doubles Configuration Interaction (SDCI) calculation on the $$\ce{CH4}$$ molecule. To set up the calculations, program MOTRA which transforms the integrals to molecular basis, and program GUGA which computes the coupling coefficients, must be run before the MRCI program. In MOTRA the reference orbitals are specifiedi, and those employed here are from an HF SCF calculation including frozen orbitals. In GUGA the reference for the CI calculation is described by the number of correlated electrons, the spatial and spin symmetry, the inactive orbitals always occupation 2 in the reference space, and the type of CI expansion.

*SDCI energy for CH4 at a fixed nuclear geometry
*File: SDCI.energy.CH4
*
&GATEWAY
coord = CH4.xyz; basis = STO-3G; group = c1
&SEWARD
&SCF
&MOTRA
Lumorb
Frozen= 1
&GUGA
Electrons = 8
Spin = 1
Inactive= 4
Active= 0
Ciall= 1
&MRCI
SDCI


To use reference orbitals from a previous CASSCF calculation, the RASSCF program will have to be run before the MOTRA module. Also, if the spatial or spin symmetry are changed for the CI calculation, the modifications will be introduced in the input to GUGA program. Many alternatives are possible for performing an MRCI calculation as shown in the next example below, in which the reference space to perform the CI is multiconfigurational:

*MRCI energy for CH4 at a fixed nuclear geometry
*File: MRCI.energy.CH4
*
&GATEWAY
Title = CH4 molecule
coord = CH4.xyz; basis = STO-3G; group = c1
&SEWARD
&SCF
&RASSCF
LumOrb
Spin= 1; Nactel= 6 0 0; Inactive= 2; Ras2= 6
&MOTRA
Lumorb
Frozen= 1
&GUGA
Electrons= 8
Spin= 1
Inactive= 2
Active= 3
Ciall= 1
&MRCI
SDCI


The MRCI program also allows the calculation of electronic energies using the ACPF method. Another Molcas program, CPF, offers the possibility to use the CPF, MCPF, and ACPF methods with a single reference function. The required input is quite similar to that for the MRCI program:

*CPF energy for CH4 at a fixed nuclear geometry
*File: CPF.energy.CH4
*
&GATEWAY
Title= CH4 molecule
coord = CH4.xyz; basis = STO-3G; group = c1
&SEWARD
&SCF
&MOTRA
Lumorb
Frozen= 1
&GUGA
Electrons= 8
Spin = 1
Inactive = 4
Active = 0
Ciall= 1
&CPF
CPF
End Of Input


Finally, Molcas can also perform closed- and open-shell coupled cluster calculations at the CCSD and CCSD(T) levels. These calculations are controlled by the CCSDT program, whose main requirement is that the reference function has to be generated with the RASSCF program. The following input is required to obtain a CCSD(T) energy for the $$\ce{CH4}$$ molecule:

*CCSDT energy for CH4 at a fixed nuclear geometry
*File: CCSDT.energy.CH4
*
&GATEWAY
Title= CH4 molecule
coord = CH4.xyz; basis = STO-3G; group = c1
&SEWARD
&RASSCF
Spin= 1; Nactel= 0 0 0; Inactive= 5; Ras2= 0
OutOrbitals
Canonical
&MOTRA
JobIph
Frozen= 1
&CCSDT
CCT


Since this is a closed-shell calculation, the RASSCF input computes a simple RHF wave function with zero active electrons and orbitals using keywords OutOrbitals and Canonical. The MOTRA must include the keyword JobIph to extract the wave function information from file JOBIPH which is automatically generated by RASSCF. Finally, the keywork CCT in program CCSDT leads to the calculation of the CCSD(T) energy using the default algorithms.

The CCSDT program in Molcas is especially suited to compute open-shell cases. The input required to obtain the electronic energy of the $$\ce{CH4^+}$$ cation with the CCSD(T) method is:

*CCSDT energy for CH4+ at a fixed nuclear geometry
*File: CCSDT.energy.CH4plus
*
&GATEWAY
Title= CH4+ molecule
coord = CH4.xyz; basis = STO-3G; group = c1
&SEWARD
&RASSCF
Spin= 2; Nactel= 1 0 0; Inactive= 4; Ras2= 1
OutOrbitals
Canonical
&MOTRA
JobIph
Frozen= 1
&CCSDT
CCT


where the RASSCF program generated the proper Restricted Open-Shell Hartree–Fock (ROHF) reference. Different levels of spin adaptation are also available.

If solvent effects are desired, Molcas includes two models: Kirkwood and PCM. Adding a solvent effect to a ground state at HF, DFT, or CASSCF levels, simply requires the inclusion of the keyword RF-input within the input for the SEWARD which calculates a self-consistend reaction field.

*DFT energy for CH4 in water at a fixed nuclear geometry
*File: DFT.energy_solvent.CH4
*
&GATEWAY
Title= CH4 molecule
coord = CH4.xyz; basis = STO-3G; group = c1
RF-input
PCM-model; solvent= water
End of RF-input
&SEWARD
&SCF
KSDFT= B3LYP


Other programs such as CASPT2, RASSI, and MOTRA require that the reaction field is included as a perturbation with keyword RFPErturbation. In the next example the correction is added at both the CASSCF and CASPT2 levels.

*CASPT2 energy for CH4 in acetone at a fixed nuclear geometry
*File: CASPT2.energy_solvent.CH4
*
&GATEWAY
Title= CH4 molecule
coord = CH4.xyz; basis = STO-3G; group = c1
RF-input
PCM-model; solvent= acetone; AAre= 0.2
End of RF-input
&SEWARD
&RASSCF
Spin= 1; Nactel= 6 0 0; Inactive= 2; Ras2= 6
&CASPT2
Frozen= 1
Multistate= 1 1
RFPert


Notice that the tesserae of the average area in the PCM model (keyword has been changed to the value required for acetone by the keyword Aare, while the default is 0.4 Å:math:^2 for water (see Section 4.2.20.1.5). More detailed examples can be found in Section 5.1.6.