# 12 Keywords

In this section the keywords of the Mrcc input file are listed in alphabetical order.

active

The active orbitals for multi-reference (active-space) CI/CC calculations can be specified using this keyword. Note that this keyword overwrites the effect of keywords nacto and nactv.

Options:
none

All orbitals are inactive (i.e., single-reference calculation).

serialno

Using this option one can select the active orbitals specifying their serial numbers. The latter should be given in the subsequent line as $$,$$,…,$$-$$,…, where $n_{i}$’s are the serial numbers of the correlated orbitals. Serial numbers separated by dash mean that $$ through $$ are active. Note that the numbering of the orbitals is relative to the first correlated orbital, that is, frozen orbitals are excluded.

vector

Using this option one can set the active/inactive feature for each correlated orbital. In the subsequent line an integer vector should be supplied with as many elements as the number of correlated orbitals. The integers must be separated by spaces. Type 1 for active orbitals and 0 for inactive ones.

Default:

active=none

Examples:
1.

We have 20 correlated orbitals. Orbitals 1, 4, 5, 6, 9, 10, 11, 12, and 14 are active. Using the serialno option the input should include the following two lines:
active=serialno
1,4-6,9-12,14

2.

The same using the vector option:
active=vector
1 0 0 1 1 1 0 0 1 1 1 1 0 1 0 0 0 0 0 0

agrid

Specifies the angular integration grid for DFT calculations. The grid construction follows the design principles of Becke [7], the smoothing function for the Voronoi polyhedra are adopted from Ref. 103 with $m_{\mu}$ = 10. Angular grids are taken from the Grid file which is located in the BASIS directory created at the installation. By default, the 6-, 14-, 26-, 38-, 50-, 74-, 86-, 110-, 146-, 170-, 194-, 230-, 266-, 302-, 350-, 434-, 590-, 770-, 974-, 1202-, 1454-, and 1730-point Lebedev quadratures [7] are included in the file, which are labeled, respectively, by LD0006, LD0014, etc. In addition to the above grids, any angular integration grid can be used by adding it to the BASIS/Grid file or alternatively to the GENBAS file to be placed in the directory where Mrcc is executed. The format is as follows. On the first line give the label of the grid as XXNNNN, where XX is any character and NNNN is the number of the grid points (see the above examples). The subsequent NNNN lines must contain the Cartesian coordinates and the weights for the grid points.

For the selection of the angular grids, by default, an adaptive scheme motivated by Ref. 72 is used. The angular grids are selected for each radial point so that the error in the angular integrals will not be larger than $10^{-{\tt grtol}}$. The important difference is that the grids are optimized for each atom separately to avoid discontinuous potential energy surfaces. For the construction of the radial integration grid see the description of keyword rgrid. See also the description of keyword grtol.

Options:
$<$name of the grid$>$

the name of the quadrature as it is specified in the BASIS/Grid (or GENBAS) file. This angular quadrature will be used in each radial point.

LDMMMM-LDNNNN

An adaptive integration grid will be used. For each radial point, depending on its distance from the nucleus, a different Lebedev grid will be selected. The minimal and maximal number of points is MMMM and NNNN, respectively.

Default:

agrid=LD0006-LD0590

Examples:
1.

for a 974-point Lebedev grid set agrid=LD0974

2.

to use an adaptive grid with at least 110 and at most 974 angular points set agrid=LD0110-LD0974

3.

for a very fine grid use
agrid=LD0110-LD0974
grtol=12

basis

Specifies the basis set used in all calculations. By default the basis sets are taken from the files named by the chemical symbol of the elements, which can be found in the BASIS directory created at the installation. The basis sets are stored in the format used by the Cfour package (see Sect. 6.9). In addition to the basis sets provided by default, any basis set can be used by adding it to the corresponding files in the BASIS directory. Alternatively, you can also specify your own basis sets in the file GENBAS which must be copied to the directory where Mrcc is executed.

Options:
$<$basis set label$>$

If the same basis set is used for all atoms, the label of the basis set must be given.

atomtype

If different basis set are used, but the basis sets are identical for atoms of the same type, basis=atomtype should be given, and the user must specify the basis sets for each atomtype in the subsequent lines as $<$atomic symbol$>$:$<$basis set$>$ .

mixed

Mixed basis sets will be used, that is, different basis sets will be used for different groups of atoms specified by their serial number. The number of groups, the basis sets, and corresponding atoms must be specified in the subsequent lines as
$<$number of groups$>$
$<$basis set label 1$>$ $$,$$,…,$$-$$,…
$<$basis set label 2$>$ $$,$$,…,$$-$$,…

where $n_{i}$’s, $m_{i}$’s, …are the serial numbers of the atoms. Serial numbers separated by dash mean that atoms $$ through $$ are included.

embed

A mixed basis set composed of two AO bases will be used in the case of an embedding calculation. It only works if keyword embed is also specified. The two basis sets must be given in the following two lines. The first basis will be used for the environment, while the second one is the AO basis for the embedded subsystem (see also the description of keyword embed).

corembed

It is the same as embed, but the partitioning defined by keyword corembed will be used.

special

In the general case, if different basis set are used for each atom, then one should give basis=special and specify the basis sets for each atom in the subsequent lines by giving the label of the corresponding basis sets in the order the atoms appear at the specification of the geometry.

Notes:
• 1.

By default the following basis sets are available for elements H to Kr in Mrcc:

• Dunning’s correlation consistent basis sets [29, 66, 160, 161, 124]: cc-pV$X$Z, cc-pCV$X$Z, aug-cc-pV$X$Z, aug-cc-pCV$X$Z ($X$ = D, T, Q, 5, 6)

• Gaussian basis sets of Pople and co-workers [50, 73, 54, 26, 37, 12, 43, 94, 19]: STO-3G, 3-21G, 6-31G, 6-311G, 6-31G*, 6-311G*, 6-31G**, 6-311G**, 6-31+G*, 6-31+G**, 6-31++G**, 6-311+G*, 6-311+G**, 6-311++G**

• the def2 Gaussian basis sets of Weigend and Ahlrichs [156]: def2-SV(P), def2-SVP, def2-TZVP, def2-TZVPP, def2-QZVP, def2-QZVPP

• the augmented def2 Gaussian basis sets of Rappoport and Furche [132]: def2-SVPD, def2-TZVPD, def2-TZVPPD, def2-QZVPD, def2-QZVPPD

• F12 basis sets for explicitly correlated wave functions developed by Peterson et al. [123]: cc-pV$X$Z-F12 ($X$ = D, T, Q)

• the Gaussian basis sets of Dunning and Hay (LANL2DZ) [28]

• the auxiliary basis sets of Weigend et al. for correlation calculations using the density fitting/resolution of the identity approximation [157, 158]: cc-pV$X$Z-RI, aug-cc-pV$X$Z-RI ($X$ = D, T, Q, 5, 6); def2-SV(P)-RI, def2-SVP-RI, def2-TZVP-RI, def2-TZVPP-RI, def2-QZVP-RI, def2-QZVPP-RI

• the auxiliary basis sets of Hellweg and Rappoport for the augmented def2 Gaussian basis sets [56]: def2-SVPD-RI, def2-TZVPD-RI, def2-TZVPPD-RI, def2-QZVPD-RI, def2-QZVPPD-RI

• Weigend’s Coulomb/exchange auxiliary basis sets for density fitting/resolution of the identity SCF calculations [159]: cc-pV$X$Z-RI-JK, aug-cc-pV$X$Z-RI-JK ($X$ = D, T, Q, 5), def2-QZVPP-RI-JK

From Na to La and from Hf to Rn the following basis sets are available, which must be used together with the corresponding ECP (see also the description of keyword ECP):

• the LANL2DZ basis sets of Hay and Wadt [52, 155, 51]

• the def2 Gaussian basis sets of Weigend and Ahlrichs [156]: def2-SV(P), def2-SVP, def2-TZVP, def2-TZVPP, def2-QZVP, def2-QZVPP

• the augmented def2 Gaussian basis sets of Rappoport and Furche [132]: def2-SVPD, def2-TZVPD, def2-TZVPPD, def2-QZVPD, def2-QZVPPD

• the correlation consistent PP basis sets of Peterson and co-workers [127, 126, 122, 125, 33]: cc-pV$X$Z-PP and aug-cc-pV$X$Z-PP ($X$ = D, T, Q, 5)

• the auxiliary basis sets of Hellweg and Rappoport for the augmented def2 Gaussian basis sets [56]: def2-SVPD-RI, def2-TZVPD-RI, def2-TZVPPD-RI, def2-QZVPD-RI, def2-QZVPPD-RI

• the auxiliary basis sets of Hättig for correlation calculations with the PP basis sets: cc-pV$X$Z-PP-RI and aug-cc-pV$X$Z-PP-RI ($X$ = D, T, Q, 5)

Please note that some of the above basis sets are not available for all elements.

• 2.

• 3.

If you use your own basis sets, these must be copied to the end of the corresponding file in the BASIS directory. Alternatively, you can also create a file called GENBAS in the directory where Mrcc is executed, and then you should copy your basis sets to that file.

• 4.

The labels of the basis sets must be identical to those used in the BASIS/* files (or the GENBAS file). For the default basis sets just type the usual name of the basis set as given above, e.g., cc-pVDZ, 6-311++G**, etc. If you employ non-default basis sets, you can use any label.

• 5.

For Dunnings’s aug-cc-p(C)V$X$Z basis sets one, two, or three additional diffuse function sets can be automatically added by attaching the prefix d-, t-, or q-, respectively, to the name of the basis set. To generate a d-aug basis set one even tempered diffuse function is added to each primitive set. Its exponent is calculated by multiplying the exponent of the most diffuse function by the ratio of the exponents of the most diffuse and the second most diffuse functions in the primitive set. If there is only one function in the set, the exponent of the most diffuse function is divided by 2.5. To generate t-aug and q-aug sets this procedure is repeated.

• 6.

For Dunnings’s basis sets, to use the aug-cc-p(C)V$X$Z set for the non-hydrogen atoms and the corresponding cc-p(C)V$X$Z set for the hydrogens give aug'-cc-p(C)V$X$Z. Then the diffuse functions will be automatically removed from the hydrogen atoms.

• 7.

Only the conventional AO basis set can be specified with this keyword. For the fitting basis sets used in density fitting approximations see the description of keywords dfbasis_*.

• 8.

The cc-pVDZ-RI-JK basis set has been generated from cc-pVTZ-RI-JK by dropping the functions of highest angular momentum. The aug-cc-pV$X$Z-RI-JK (def2-QZVPPD-RI-JK) basis sets are constructed automatically from the corresponding cc-pV$X$Z-RI-JK (def2-QZVPP-RI-JK) sets by adding diffuse functions as described above for the d-aug-cc-p(C)V$X$Z basis sets.

• 9.

For Dunnings’s and Pople’s basis sets add the -min postfix to the basis set name to generate a minimal basis set dropping all the polarization (correlation) functions.

• 10.

If the (aug-)cc-pV$X$Z-PP basis set does not exist for an element with $Z\leq 28$, the program will automatically attempt to use the corresponding aug-cc-pV$X$Z basis instead.

Default:

none, that is, the basis set must be specified (excepting the case when Mrcc is used together with another code, that is, iface $\neq$ none).

Examples:
1.

Consider any molecule and suppose that the cc-pVDZ basis set is used for all atoms. The input must include the following line:
basis=cc-pVDZ

2.

To use Dunning’s doubly augmented cc-pVDZ basis set (d-aug-cc-pVDZ) for all atoms the input must include the following line:
basis=d-aug-cc-pVDZ

3.

Consider the water molecule and use the cc-pVDZ basis set for the hydrogens and cc-pVTZ for the oxygen. The input must include the following lines:
basis=atomtype
O:cc-pVTZ
H:cc-pVDZ

4.

Consider water again and use the cc-pVQZ, cc-pVTZ, and cc-pVDZ basis sets for the oxygen atom, for the first hydrogen, and for the second hydrogen, respectively. Note that the order of the basis set labels after the basis=special statement must be identical to the order of the corresponding atoms in the Z-matrix/Cartesian coordinates:
geom
O
H 1 R
H 1 R 2 A

R=0.9575
A=104.51

basis=special
cc-pVQZ
cc-pVTZ
cc-pVDZ

5.

Consider the water molecule and use the cc-pVTZ basis set for the hydrogens and aug-cc-pVTZ for the oxygen. The following two inputs are identical:
basis=atomtype
O:aug-cc-pVTZ
H:cc-pVTZ

or
basis=aug'-cc-pVTZ

6.

Consider the water molecule. If you specify
basis=cc-pVTZ-min
minimal basis sets generated from cc-pVTZ will be used for the atoms, that is, only one $s$ function (two $s$ and one $p$ shells) will be retained from the $s$$p$ kernel of the H (O) cc-pVTZ basis set.

7.

Consider the PbO molecule. If you want to use the cc-pVDZ basis set for O and the cc-pVDZ-PP basis with the corresponding ECP for Pb, you only need to set
basis=cc-pVDZ-PP
in the MINP file.

8.

Mixed basis approach with two basis sets, the cc-pVTZ basis is used for atoms 1, 2, 3, and 5, while cc-pVDZ is employed for atoms 4, 6, 7, 8:
basis=mixed
2
cc-pVTZ 1-3,5
cc-pVDZ 4,6-8

basis_sm

Specifies the small basis set used in dual basis-set calculations as well as for generating SCF initial guess (scfiguess=small).

Options:

the options are the same as for keyword basis, but there is an additional one, none, which means that no small basis is defined.

Default:

basis_sm=none

Examples:
1.

To restart an SCF calculation with the cc-pVQZ basis set from the densities obtained with the cc-pVDZ basis give
basis=cc-pVQZ
basis_sm=cc-pVDZ
scfiguess=small

2.

To perform a dual basis set DF-HF calculation with the 6-311G** and 6-31G** basis sets you need:
basis=6-311G**
basis_sm=6-31G**
dual=on
calc=DF-HF

basopt

Use this keyword to turn on/off basis set optimization. Besides setting this keyword a user supplied GENBAS file is also required for basis set optimization jobs. It is also possible to set the value of basopt to be equal to an appropriate energy. In this case the basis set parameters are optimized so that the absolute value of the difference between this value and the actual energy is minimized. This option comes handy when optimizing a density fitting basis set. In this case the difference between the actual and non-density-fitting energy (obtained from a previous calculation) will be minimized. See also Sect. 6.9.

Options:

on, off, or $<$any real number$>$

Default:

basopt=off

Examples:
1.

To optimize a basis set variationally set basopt=on

2.

To optimize a basis set minimizing the difference of the calculated energy and -76.287041 E${}_{h}$ set basopt=-76.287041

bfbasis

Specifies the bond function (BF) basis (see Ref. 95 for details).

Options:
none

No BFs are used.

$<$BF basis name$>$

name of the BF basis to be used.

Notes:
• 1.

The format of the name of the BF basis, $<$BF basis name$>$, is $<$AO basis name$>$-$<$BF type$>$. E.g., 6-31G-1s1p is a BF basis optimized for the 6-31G AO basis and one s and one p function set are placed on the corresponding bonds.

• 2.

The BF basis sets are stored in the BASIS/Bond file but the BF basis can also be specified in the GENBAS file similar to the AO basis sets (see the description of keyword basis). The format of the label of the BF basis in the file is B$<$bond name$>$:$<$BF basis name$>$. E.g., BCH:6-31G-1s1p is 6-31G-1s1p BF basis optimized for the C–H bond.

• 3.

If BF bases are used, the geometry must be given in mol format (see the description of keyword geom)

Default:

bfbasis=none

Example:

hydrogen-fluoride molecule, the 6-31G basis and the 6-31G-1s1p bond function basis are used:
basis=6-31G
bfbasis=6-31G-1s1p
geom=mol
2 1
0.00000000 0.00000000 0.00000000 F
0.00000000 0.00000000 0.91690000 H
1 2 1

bpcompo

Boughton–Pulay completeness criterion [16] for occupied orbitals. In various local correlation approaches the Boughton–Pulay procedure is used to identify the atoms on which an LMO is localized. The least-squares residual of the parent LMO and the LMO truncated to the selected atoms is required to be less than one minus this criterion.

Options:
$<$any real number in the [0,1] interval$>$

This number will be used as the completeness criterion.

Default:

bpcompo=0.985

Note:

Atom domains determined by bpcompo are also utilized to construct local fitting domains in the case of localcc=2016 or 2018 according to Ref. 105.

Example:

to set a threshold of 0.99 type bpcompo=0.99

bpcompv

Boughton–Pulay completeness criterion [16] for virtual orbitals (projected atomic orbitals). See also keyword bpcompo.

Options:
$<$any real number in the [0,1] interval$>$

This number will be used as the completeness criterion.

Default:

bpcompv=0.98

Example:

to set a threshold of 0.95 type bpcompv=0.95

bpedo

Boughton–Pulay completeness criterion [16] for the occupied orbitals of an extended domain. See also keyword bpcompo.

Options:
$<$any real number in the [0,1] interval$>$

This number will be used as the completeness criterion.

Default:

bpedo=0.9999 is set as default in the case of localcc=2018 and lcorthr=normal for both LMP2 and LNO-CC computations, according to Refs. 105 and 106. See the description of lcorthr for further details on the predefined values of bpedo for other cases

Note:

bpedo=bpcompo is set if bpedo is not specified and not employed in the local correlation calculation

Example:

to set a threshold of 0.9998 type bpedo=0.9998

bpedv

Boughton–Pulay completeness criterion [16] for the virtual orbitals (projected atomic orbitals) of an extended domain. See also keyword bpcompo.

Options:
$<$any real number in the [0,1] interval$>$

This number will be used as the completeness criterion.

Default:

bpedv=0.995 is set as default in the case of localcc=2016 or 2018 according to Ref. 105.

Note:

bpedv=bpcompv is set if bpedv is not specified and not employed in the local correlation calculation

Example:

to set a threshold of 0.99 type bpedv=0.99

bppdo

Boughton–Pulay completeness criterion [16] for the occupied orbital of a primary domain. See also keyword bpcompo.

Options:
$<$any real number in the [0,1] interval$>$

This number will be used as the completeness criterion.

Default:

bppdo=0.999 is set as default in the case of localcc=2016 or 2018 according to Ref. 105.

Note:

bppdo=bpcompo is set if bppdo is not specified and not employed in the local correlation calculation

Example:

to set a threshold of 0.99 type bppdo=0.99

bppdv

Boughton–Pulay completeness criterion [16] for virtual orbitals (projected atomic orbitals) of a primary domain. See also keywords bppdo and bpcompo.

Options:
$<$any real number in the [0,1] interval$>$

This number will be used as the completeness criterion.

Default:

bppdv=bpcompv

Example:

to set a threshold of 0.99 type bppdv=0.99

calc

Specifies the type of the calculation.

Options:
SCF or HF

Hartree–Fock SCF calculation, the type of the Hartree–Fock wave function can be controlled by keyword scftype (see also keyword scftype).

RHF, UHF, ROHF

Restricted, unrestricted, or restricted open-shell Hartree–Fock SCF calculation, respectively. The type of the Hartree–Fock wave function is also defined at the same time if these options are chosen, and it is not necessary to set scftype. That is, calc=RHF is equivalent to calc=SCF plus scftype=RHF, etc.

B3LYP, PBE0, B3PW91, B3LYP-D3, B2PLYP-D3, …

Kohn–Sham SCF calculation with the specified density functional. The type of the Kohn–Sham procedure (i.e., RKS of UKS) can be controlled by keyword scftype (see also keyword scftype). The options are identical to those of keyword dft (except for off, user, and userd), see the description of keyword dft. Note that for a correlated calculation with KS orbitals you can only select the functional with keyword dft, the value of keyword calc must be set to the desired correlation method. Note also that for DFT calculations the density fitting approximation is used by default, i.e., dfbasis_scf is set to auto. To run a conventional KS calculation set dfbasis_scf=none.

TDHF

Time-dependent HF (TD-HF, also known as random-phase approximation). If calc=SCF and number of the states is greater than one (set by keywords nsing, ntrip, or nstate), also TD-HF calculations are performed for the excited states. It is only available with density fitting.

TDDFT

Full time-dependent DFT (TD-DFT). The density functional must be set using keyword dft. Alternatively, if calc is set to the name of the functional, and the number of the states is greater than one (set by keywords nsing, ntrip, or nstate), also TD-DFT calculations are performed for the excited states using the given functional. For HF reference it is equivalent to TD-HF. It is only available with density fitting.

TDA

TD-DFT in the Tamm–Dancoff approximation (TDA). For HF reference it is equivalent to CIS. It is only available with density fitting.

MP2

Second-order Møller–Plesset (MP2) calculation, the spin-component scaled MP2 (SCS-MP2) [46] and the scaled opposite-spin MP2 (SOS-MP2) [63] energy will also be computed (see also keywords scsps and scspt). Note that efficient MP2 calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-MP2 ($\equiv$ RI-MP2) calculation is performed (that is, options MP2, DF-MP2, and RI-MP2 are synonyms). If you are still interested in the MP2 energy without DF, you can, e.g., run a CCSD calculation (without DF), where the MP2 energy is also calculated.

SOS-MP2

Scaled opposite-spin second-order Møller–Plesset (SOS-MP2) calculation [63] using an $N^{4}$-scaling algorithm based on the Cholesky decomposition/Laplace transform of energy denominators (in practice one dRPA iteration is performed, see below). Note that it is only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-SOS-MP2 ($\equiv$ RI-SOS-MP2) calculation is performed (that is, options SOS-MP2, DF-SOS-MP2, and RI-SOS-MP2 are synonyms).

SCS-MP2

For canonical calculations it is equivalent to option MP2. If a local correlation calculation is executed, only the spin-component scaled MP2 (SCS-MP2) energy will be computed.

dRPA

Direct random-phase approximation (dRPA) calculation (see Eqs. 7 and 8 in Ref. 152). Note that dRPA calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-dRPA ($\equiv$ RI-dRPA) calculation is performed (that is, options dRPA, DF-dRPA, and RI-dRPA are synonyms).

RPA

Random-phase approximation (RPA) calculation (see Eqs. 10 and 13 in Ref. 152, where it is referred to as RPAx-SO2). Note that RPA calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-RPA ($\equiv$ RI-RPA) calculation is performed (that is, options RPA, DF-RPA, and RI-RPA are synonyms).

SOSEX

Second-order screened exchange (SOSEX) [48] calculation (see Eqs. 7 and 9 in Ref. 152), the dRPA energy is also computed. Note that SOSEX calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-SOSEX ($\equiv$ RI-SOSEX) calculation is performed (that is, options SOSEX, DF-SOSEX, and RI-SOSEX are synonyms).

RPAX2

RPAX2 calculation (see Eqs. 17 to 19 in Ref. 58). Note that RPAX2 calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-RPAX2 ($\equiv$ RI-RPAX2) calculation is performed (that is, options RPAX2, DF-RPAX2, and RI-RPAX2 are synonyms).

CIS

Configuration interaction singles (CIS) calculation [96]. Efficient CIS calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-CIS ($\equiv$ RI-CIS) calculation is performed (that is, options CIS, DF-CIS, and RI-CIS are synonyms). If you are still interested in the CIS energy without DF, set ccprog=mrcc, dfbasis_scf=none, and dfbasis_cor=none.

CIS(Di)

Iterative doubles correction to configuration interaction singles [CIS(D${}_{\infty}$)] calculation [53, 96]. Note that CIS(D${}_{\infty}$) calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-CIS(D${}_{\infty}$) [$\equiv$ RI-CIS(D${}_{\infty}$)] calculation is performed [that is, options CIS(Di), DF-CIS(Di), and RI-CIS(Di) are synonyms].

Second-order algebraic diagrammatic construction [ADC(2)] calculation [141, 96]. Note that ADC(2) calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-ADC(2) [$\equiv$ RI-ADC(2)] calculation is performed [that is, options ADC(2), DF-ADC(2), and RI-ADC(2) are synonyms].

CC2

Second-order coupled-cluster singles and doubles (CC2) calculation [18, 96]. Efficient CC2 calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-CC2 ($\equiv$ RI-CC2) calculation is performed (that is, options CC2, DF-CC2, and RI-CC2 are synonyms). If you are still interested in the CC2 energy without DF, set ccprog=mrcc, dfbasis_scf=none, and dfbasis_cor=none.

CCS, CCSD, CCSDT, CCSDTQ, CCSDTQP, CC($<$n$>$)

The corresponding single-reference CC calculation if the number of active orbitals is zero (see Ref. 82); the corresponding SRMRCCSD, SRMRCCSDT, etc. calculation otherwise (see Ref. 83).

CCSD[T], CCSDT[Q], CCSDTQ[P], CC($<$n-1$>$)[$<$n$>$]

The corresponding single-reference CC calculation with perturbative corrections (see Ref. 77).

CCSD(T), CCSDT(Q), CCSDTQ(P), CC($<$n-1$>$)($<$n$>$)

The corresponding single-reference CC calculation with perturbative corrections (see Ref. 77).

CCSD(T)_L, CCSDT(Q)_L, CCSDTQ(P)_L, CC($<$n-1$>$)($<$n$>$)_L

The corresponding CCSD(T)${}_{\Lambda}$, CCSDT(Q)${}_{\Lambda}$, etc. calculation (see Ref. 77).

CCSDT-1a, CCSDTQ-1a, CCSDTQP-1a, CC($<$n$>$)-1a

The corresponding iterative approximate single-reference CC calculation (see Ref. 77).

CCSDT-1b, CCSDTQ-1b, CCSDTQP-1b, CC($<$n$>$)-1b

The corresponding iterative approximate single-reference CC calculation (see Ref. 77).

CC2, CC3, CC4, CC5, CC$<$n$>$

The corresponding iterative approximate single-reference CC calculation (see Ref. 77).

CCSDT-3, CCSDTQ-3, CCSDTQP-3, CC($<$n$>$)-3

The corresponding iterative approximate single-reference CC calculation (see Ref. 77).

CCSDT[Q]/A, CCSDTQ[P]/A, CC($<$n-1$>$)[$<$n$>$]/A

The corresponding single-reference CC calculation with perturbative corrections using ansatz A (see Ref. 79).

CCSDT[Q]/B, CCSDTQ[P]/B, CC($<$n-1$>$)[$<$n$>$]/B

The corresponding single-reference CC calculation with perturbative corrections using ansatz B (see Ref. 79).

CCSDT(Q)/A, CCSDTQ(P)/A, CC($<$n-1$>$)($<$n$>$)/A

The corresponding single-reference CC calculation with perturbative corrections using ansatz A (see Ref. 79).

CCSDT(Q)/B, CCSDTQ(P)/B, CC($<$n-1$>$)($<$n$>$)/B

The corresponding single-reference CC calculation with perturbative corrections using ansatz B (see Ref. 79).

CIS, CISD, CISDT, CISDTQ, CISDTQP, CI($<$n$>$), FCI

The corresponding single-reference CI calculation if the number of active orbitals is zero (see Ref. 82), the corresponding MRCISD, MRCISDT, etc. calculation otherwise (see Ref. 83).

Notes:
• 1.

In the above options $n$ is a positive integer, which is the excitation level of the highest excitation. $n$ is supposed to be equal to or greater than 6 since for smaller $n$’s the CC($<$n$>$) and similar options are equivalent to one of the other options, e.g., CC(5) is equivalent to CCSDTQP or CC(3)(4) is identical with CCSDT(Q).

• 2.

For excited-state calculations with the TD-HF, TDA, TD-DFT, CIS, CIS(D${}_{\infty}$), ADC(2), CC2 and various CC and CI methods the number of states should be greater than one (keywords nsing, ntrip, or nstate). If more than one state is requested for CC calculations, the corresponding linear-response (LR) CC (for excitation energies it is equivalent to equation-of-motion CC, EOM-CC) calculation is performed automatically for the excited states. If more than one state is requested and calc=SCF, TD-HF (dft=off) or TD-DFT (dft$\neq$off) calculations will be carried out for the excited states.

• 3.

The active orbitals can be selected and the MRCI/CC calculations can be controlled by keywords nacto, nactv, active, maxex, and maxact.

• 4.

In principle, all methods can be used with the density fitting (resolution-of-identity) approximation. It is possible in two ways. You can attach the prefix DF- or RI- to the corresponding option from the above list. Then, for a HF calculation keyword dfbasis_scf will be set to auto, while for a correlated calculation both dfbasis_scf and dfbasis_cor will be given the value auto. Alternatively, you can also set the values for keywords dfbasis_scf and dfbasis_cor, see their description.

• 5.

Local correlation methods can be run if the prefix “L” is added to the corresponding option of the keyword, e.g., as LMP2, LdRPA, LCCSD(T), etc. Additionally, the prefix “LNO-” can also be used as a synonym in the case of local coupled-cluster approaches, e.g., as LNO-CCSD, LNO-CCSD(T), LNO-CCSDT, etc. Both options are equivalent to setting localcc=on.

• 6.

For the dRPA, RPA, and SOSEX methods the use of PBE orbitals is recommended.

• 7.

For the RPAX2 method the use of PBEx orbitals is recommended.

Default:

calc=SCF

Examples:
1.

To run a CCSD(T) calculation the user should set calc=CCSD(T)

2.

For DF-HF (RI-HF) calculations type:
calc=DF-HF
which is equivalent to the following input:
calc=SCF
dfbasis_scf=auto

3.

For a local CCSD(T) calculation using the local natural orbital approximation set calc=LCCSD(T) or calc=LNO-CCSD(T)

4.

For a RI-MP2 calculation set calc=MP2

5.

For a DFT calculation with the B3LYP functional set calc=B3LYP

6.

Direct RPA calculation with Kohn–Sham orbitals calculated with the PBE functional:
calc=dRPA
dft=PBE

7.

TD-DFT calculation for the 3 lowest singlet excited states of a molecule using the PBE functional:
calc=TDDFT
dft=PBE
nsing=4

A somewhat less complicated input for the same purpose:
calc=PBE
nsing=4

ccmaxit

Maximum number of iteration steps in correlated calculations (CC, CI, RPA, …).

Options:

$<$any positive integer$>$

Default:

ccmaxit=50

Example:

to increase the maximum number of CC iterations to 100 give ccmaxit=100

ccprog

Specifies the CC program to be used.

Options:
mrcc

The automated, string-based CC program mrcc will be called.

ccsd

The very fast, hand-coded CCSD(T) codes, ccsd or uccsd, will be executed (currently the spatial symmetry cannot be utilized).

cis

The very fast, hand-coded, integral direct DF-CIS code cis will be executed (currently the spatial symmetry cannot be utilized).

Note:

Please note that the mrcc code was optimized for high-order CC calculations, such as CCSDT(Q) and CCSDTQ, which require different algorithms than CCSD(T). Thus it is slow for CCSD(T), but optimal for high-order CC models.

Default:

ccprog=ccsd for CCSD and CCSD(T) calculations, ccprog=cis for CIS, CIS(D${}_{\infty}$), ADC(2), and CC2 calculations, ccprog=mrcc otherwise.

Example:

to use the mrcc code for CCSD or CCSD(T) calculations give ccprog=mrcc

cctol

Convergence threshold for the energy in correlated calculations (CC, CI, dRPA, RPA, etc.). The energy will be accurate to $10^{-{\tt cctol}}$ E${}_{h}$.

Options:

$<$any integer$>$

Default:

cctol=8 for property calculations, cctol=[-$\log_{10}$(optetol)]+2 for geometry optimizations, cctol=5 for localcc=2016 and
localcc=2018, cctol=6 otherwise

Example:

for an accuracy of $10^{-8}$ E${}_{h}$ one must give cctol=8

charge

Charge of the system.

Options:

$<$any integer$>$

Default:

charge=0

Example:

for the Cl${}^{-}$ ion one should give charge=-1

cialg

Specifies what type of algorithm is to be used in CIS, TDA, TD-HF, and TD-DFT calculations.

Options:
disk

Conventional algorithm, two-electron integrals are stored on disk

direct

Completely I/O-free, integral-direct algorithm, two-electron integrals are recalculated in each iteration step.

direct2

Partially I/O-free, integral-direct algorithm; recommended if the I/O is fast and/or few states are required.

direct3

Variant of direct2, but usually slower.

auto

Based on the size of the molecule the program will automatically select the most efficient one from the above options.

Default:

cialg=auto

Example:

to use disk-based algorithm set cialg=disk

ciguess

The initial guess vectors for CI and LR-CC calculations can be specified using this keyword.

Options:
on

The initial trial vectors are supplied by the user and should be given in the subsequent lines as follows. For each state the corresponding initial guess vector must given by the number of non-zero elements of the vector on the first line, followed by as many lines as the number of non-zero elements. In each line the corresponding excitation operator and the value for this element of the vector must be provided in the following format:
$\dots\dots\dots<% i_{n}>$
where $$ is the level of excitation, and the electrons are promoted from occupied orbitals $\dots$ to virtual orbitals $\dots$ with spins $\dots$ ($$ is 1 for alpha and 0 for beta), respectively. $$ is the corresponding coefficient.

off

Initial trial vectors are not specified, the program applies simple unit vectors as initial guess. The unit vectors are determined on the basis of the diagonal elements of the Hamiltonian: if $n$ roots are requested, $n$ unit vectors corresponding to the $n$ lowest diagonals will be used.

Default:

ciguess=off

Example:

Suppose that we have two excited states in a LR-CC calculation. Then the initial guess can be given as follows.
ciguess=on
1
1 1 6 4 1.0
3
1 1 7 3 0.1
2 1 0 7 7 5 5 1.0
2 1 1 7 6 3 4 0.1

For the first state there is only one entry, a single excitation of the alpha electron from orbital 4 to orbital 6 with a coefficient of 1.0. For the second root the initial guess vector contains three entries. A single excitation from orbital 3 to orbital 7 with alpha spin and a relative weight of 0.1, a double excitation from orbital 5 to orbital 7 with a weight of 1.0, and another double excitation of the alpha electrons from orbitals 3 and 4 to orbitals 6 and 7 with a weight of 0.1.

Notes:
1.

For $M_{S}=0$ states the vector is automatically spin-adapted, and you do not need to specify the coefficients for the corresponding spin-reversed excitations. E.g., in the above example, for root 1 the 1 0 6 4 1.0 entry is unnecessary.

2.

The guess vector is not required to be normalized, it is done automatically.

3.

In the case of four-component relativistic calculations (Dirac interface) the serial numbers of the spinors should be specified. In addition, the second number in the above strings must be 1 (that is, all excitations are formally considered as excitations of alpha electrons).

cmpgrp

Specifies the computational point group. All calculations will use the specified Abelian group. See Sect. 13 for more details.

Options:
auto

The molecular symmetry is automatically recognized.

$<$point group symbol$>$

Schönflies symbol of the Abelian point group such as C1, Ci, Cs, C2, C2v, C2h, D2, D2h

Note:

cmpgrp=C1 is equivalent to symm=off

Default:

cmpgrp=auto

Example:

to use $C_{2v}$ point group for benzene set cmpgrp=C2v

core

Specifies whether the core electrons are correlated.

Options:
frozen

Frozen core approximation

corr

All core electrons are correlated

$<$any non-negative integer n$>$

The lowest (according to orbital energy order) $n$ pieces of spatial orbitals (the lowest $n$ pieces of alpha and $n$ pieces of beta spin orbitals for UHF/semicanonical ROHF reference) will be dropped.

Default:

core=frozen

Example:

to correlate all core electrons set core=corr or core=0

corembed

This keyword controls the models and subsystems selected for multi-level local correlation methods. Currently it is only available for closed-shell systems using density-fitting.

Options:
off

Conventional case, a single model defined by calc is used for the entire system.

on

Multi-level calculation is performed with different local correlation methods for the active (high-level) and the environmental (low-level) subsystems. The three input lines following corembed define the list of active atoms, the computational model for the environment level, and the number of embedded orbitals (if it is specified). The syntax for these three lines is analogous with that for keyword embed. (See the description of keyword embed.) The high-level method for the active region should be specified by the keyword calc.

Default:

corembed=off

Notes:
1.

Local correlation methods available with localcc=2015,
localcc=2016, and localcc=2018 (e.g., MP2 or arbitrary single-reference CC) can be chosen for both the active and the environmental subsystem. Additionally, HF or HF+LRC are also available choices for the low-level model. If the latter is set, the environment is treated at the HF level but the long-range correlation (LRC) between the active subsystem and its environment is also taken into account (see Ref. 62). Note that models with KS-DFT reference, such as dRPA, SOSEX, etc., are not available for multi-level local correlation calculations.

2.

The threshold settings of the local correlation method chosen for the high-level model can be given (as in the case of corembed=off) by the keywords controlling the local correlation methods (see their list in Sect. 6.8). Default settings according to lcorthr=normal and localcc=2018 (or for previous versions according to lcorthr=loose and localcc=2015 or localcc=2016) are employed for the low-level model of the environment.

Examples:
1.

LNO-CCSD(T)-in-LMP2 scheme, where LNO-CCSD(T) is performed for the active orbitals with tight thresholds, atoms 1 and 2 are included in the high-level region, and the number of the active orbitals is determined automatically:
calc=LNO-CCSD(T)
lcorthr=tight
corembed=on
1-2
LMP2
0

2.

LNO-CCSDT-in-LNO-CCSD scheme, where the local CCSDT calculation is performed with the mrcc program for the active orbitals and the local CCSD is calculation performed with the ccsd program for the environment:
calc=LCCSDT
corembed=on
1-2
LCCSD
0

3.

LNO-CCSD(T)-in-HF+LRC embedding where only HF is used for the environment but the additional LRC term accounts for the interaction of the active and environmental parts. Atoms 1, 2, 3, and 5 define the active subsystem, and 10 orbitals are included in the active region:
calc=LNO-CCSD(T)
corembed=on
1-3,5
HF+LRC
10

dboc

Diagonal Born–Oppenheimer correction (DBOC) (available only with Cfour).

Options:

on or off

Default:

dboc=off

Example:

for a DBOC calculation set dboc=on

dendec

Selects the algorithm for the decomposition of energy denominators, Cholesky-decomposition or Laplace transform, for canonical SOS-MP2 and dRPA (also required for SOSEX) as well as for local MP2 and dRPA calculations. The dRPA calculation is performed using the modified algorithm of Heßelmann [58] based on the decomposition of energy denominators. For the calculation of the SOS-MP2 energy, in practice one dRPA iteration is performed with the aforementioned algorithm. In the case of local MP2 and dRPA calculations the correlation energy contributions are also evaluated with the aid of the decomposition of energy denominators (see Ref. 85). The algorithm for the decomposition can be set using this keyword in all of the above cases. The number of retained Cholesky vectors/quadrature points can be controlled by keyword nchol.

Options:
Cholesky

Cholesky decomposition will be used

Laplace

Laplace transform will be used

Default:

dendec=Laplace for SOS-MP2, dendec=Cholesky otherwise

Notes:
1.

The algorithms based on the Laplace-transformed technique use minimax quadratures obtained from Ref. 150.

2.

The default quadratures are taken from the Quad file which is located in the BASIS directory created at the installation. In addition to the default quadratures, any further quadrature can be used by adding it to the BASIS/Quad file or alternatively to the GENBAS file to be placed in the directory where Mrcc is executed. The format is as follows. On the first line give the label of the quadrature as KNNRXXX, where NN is the number of the quadrature points and XXX is the upper limit of the interval in which the Laplace transform is approximated (variable $R$ in Ref. 150). The subsequent NN lines must contain, respectively, the weights and quadrature points.

Example:

to use Laplace transform give dendec=Laplace

dens

Construction of density, derivative density, and transition density matrices for property calculations. If mod(dens,2)=1, only one-particle, if mod(dens,2)=0, both one- and two-particle density matrices will be calculated and contracted with the available property integrals. See Refs. 74, 75, 76, 78, 112, 113 for more details.

Options:
1, 2

Density-matrix calculation (for geometry optimizations, first-order properties, etc.)

3, 4

Density-matrix first derivatives (for second-order property calculations, available only with Cfour)

5, 6

Transition density matrices (for transition moment calculations)

7, 8

Second and third derivatives of the density-matrix (for third-order property calculations, available only with Cfour)

Default:

dens=2 for geometry optimizations and QM/MM calculations, dens=0 otherwise

Notes:
1.

Transition moment as well as excited-state gradient calculations can be performed for only one excited state at a time, that is, nsing, ntrip, or nstate cannot exceed 2. To compute the transition moment or gradient for a higher excited state you need to converge the equations to that root. The best practice is to run a calculation with the desired number of excited states, and then restart the calculation selecting a higher solution (see the description of keyword rest). You can also try to start the calculation from a good initial guess (see the description of keyword ciguess).

2.

If dens $\neq$ 0, a population analysis is also performed, and Mulliken and Löwdin atomic charges as well as Mayer bond orders are computed.

Example:

for the calculation of both one- and two-particle density matrices set dens=2

dfalg

Specifies how the inverse of the two-center Coulomb integral matrix is decomposed in density fitting direct SCF calculations.

Options:
LinEq

The fitting coefficients are computed by solving the corresponding system of linear equations. It is efficient and numerically stable. It is the best choice for very large auxiliary basis sets for which the diagonalization of the two-center integral matrix is prohibitive.

InvSqrt

Inverse square root of the two-center integral matrix is used. It is relatively stable numerically, but the diagonalization is slow and requires much memory.

Cholesky

Cholesky decomposition of the inverse of the two-center integral matrix is used. It is an efficient algorithm but numerically unstable if the two-center matrix tends to be singular.

Default:

dfalg=InvSqrt for property calculations, dfalg=LinEq otherwise

Example:

to use Cholesky decomposition set dfalg=Cholesky

dfbasis_cor

Specifies whether the density fitting approximation will be used in the correlated calculations and also specifies the fitting basis set.

Options:
none

The density fitting approximation is not used for the correlated calculation.

$<$basis set label$>$, atomtype, special

The density fitting approximation is invoked, and the specified basis set is used as fitting basis set. For the specification of the basis the same rules apply as for keyword basis, see the description of keyword basis.

auto

This option can only be used if Dunning’s (aug-)cc-pV$X$Z, Weigend and Ahlrichs’ def2, the augmented def2 basis sets of Rappoport and Furche, Peterson’s cc-pV$X$Z-F12 or (aug-)cc-pV$X$Z-PP, or Pople’s basis sets are used as the normal basis set. In this case, if dfbasis_cor=auto, the density fitting approximation is invoked. For the (aug-)cc-pV$X$Z(-PP) basis sets the corresponding (aug-)cc-pV$X$Z(-PP)-RI basis sets will be used automatically as the fitting basis sets, while for a cc-pV$X$Z-F12 basis set the corresponding aug-cc-pV$X$Z-RI basis will be taken. For the (augmented) def2 basis sets also the corresponding RI basis sets will be used, e.g., def2-TZVPP-RI for def2-TZVPP, def2-QZVPP-RI for def2-QZVPP, def2-TZVPPD-RI for def2-TZVPPD, etc. For Pople-type minimal and double-$\zeta$ basis sets (i.e., STO-3G, 3-21G, 6-31G**, etc.) the cc-pVDZ-RI basis set, while for triple-$\zeta$ basis sets (i.e., 6-311G, 6-311G**, etc.) the cc-pVTZ-RI basis set will be used as the auxiliary basis; if the basis also includes diffuse functions (i.e., 6-31+G**, 6-311++G**, etc.) the aug-cc-pVDZ-RI and aug-cc-pVTZ-RI basis sets are employed by default.

Notes:
1.

For the available fitting basis sets see the notes for keyword basis on page 1..

2.

The density fitting approximation can also be invoked by attaching the prefix DF- or RI- to the corresponding option of keyword calc, see the description of calc.

Default:

dfbasis_cor=auto for all the correlation methods that use the density fitting approximation by default as well as for local correlation calculations (i.e., localcc $\neq$ off), dfbasis_cor=none otherwise.

Examples:
1.

To use the cc-pVTZ-RI fitting basis in the correlated calculation for all atoms the input must include dfbasis_cor= cc-pVTZ-RI

2.

Consider the water molecule and use the cc-pVTZ-RI fitting basis set for the hydrogens and aug-cc-pVTZ-RI for the oxygen. The following inputs are equivalent:
dfbasis_cor=atomtype
O:aug-cc-pVTZ-RI
H:cc-pVTZ-RI

or
dfbasis_cor=aug'-cc-pVTZ-RI

3.

Consider the water molecule and use the cc-pVTZ (cc-pVTZ-RI) basis set (fitting basis set) for the hydrogens and aug-cc-pVTZ (aug-cc-pVTZ-RI.) for the oxygen in a local correlation calculation. The following inputs are equivalent:
calc=CCSD(T)
localcc=on
basis=aug'-cc-pVTZ
dfbasis_scf=aug'-cc-pVTZ-RI
dfbasis_cor=aug'-cc-pVTZ-RI

or
calc=LCCSD(T)
basis=aug'-cc-pVTZ

4.

To run a DF-HF calculation with the cc-pVTZ-F12 basis set and the aug-cc-pVTZ-RI auxiliary basis the input should only include the following lines:
basis=cc-pVTZ-F12
calc=DF-HF

dfbasis_scf

Specifies whether the density fitting approximation will be used in the HF- or KS-SCF calculation and also specifies the fitting basis set. For the syntax see the description of keyword dfbasis_cor. The important difference is that, if dfbasis_scf=auto, the (aug-)cc-pV$X$Z-RI-JK basis sets will be used as auxiliary basis sets for Dunning’s, Peterson’s, and Pople’s basis sets, while for the def2 basis sets the def2-QZVPP-RI-JK auxiliary basis is taken. For the augmented def2 as well as for the aug-cc-pV$X$Z-PP basis sets the def2-QZVPPD-RI-JK auxiliary basis will be used.

Default:

dfbasis_scf=auto if dfbasis_cor$\neq$none and for DFT calculations, dfbasis_scf=none otherwise.

dfintran

Specifies the integral transformation program to be used for the transformation of three-center Coulomb integrals.

Options:
drpa

the drpa program will be called

ovirt

the ovirt program will be called

Default:

dfintran=ovirt if ovirt$\neq$off, dfintran=drpa otherwise.

Example:

to use the ovirt code set dfintran=ovirt

dft

Use this keyword to perform DFT calculations and to specify the functional.

Options:
off

No DFT calculation is carried out.

$<$functional name$>$

The name of the functional, see Table LABEL:FuncTable for the available functionals.

$<$Libxc identifier$>$

The identifier of a functional implemented in the Libxc library (if installed), such as LDA_X, LDA_C_VWN_1, GGA_X_B88, etc. (see the homepage of the Libxc project [59]).

user

User-defined functional. Any combination of the following contributions can be defined:

• the available standalone functionals, see column “User” in Table LABEL:FuncTable.

• the functionals available in the Libxc library (if installed), use simply the Libxc identifier of the functionals (see the homepage of the Libxc project [59]).

• the HF exchange, denoted by HFx

• the MP2, dRPA, and SOSEX correlation, denoted, respectively, by MP2, dRPA, and SOSEX;

• the antiparallel- and parallel-spin components of the latter correlation corrections, add the s and t postfix to the above labels, respectively, e.g., instead of the MP2 label, the MP2s and MP2t labels should be used.

Note that for hybrid functionals, such as B97, the HF exchange will be neglected. The combination should be specified in the subsequent lines as follows (see also the examples below):
$<$number of entries$>$
$<$coefficient 1$>$ $<$functional name 1$>$
$<$coefficient 2$>$ $<$functional name 2$>$
$<$coefficient 3$>$ $<$functional name 3$>$

userd

User-defined functional, but different functionals are used for the calculation of the density and the energy. It is useful for defining special double-hybrid functionals. The combination should be specified in the subsequent lines as follows (see also the examples below):
$<$number of entries for density$>$
$<$coefficient 1$>$ $<$functional name 1$>$
$<$coefficient 2$>$ $<$functional name 2$>$
$<$coefficient 3$>$ $<$functional name 3$>$

$<$number of entries for energy$>$
$<$coefficient 1’$>$ $<$functional name 1’$>$
$<$coefficient 2’$>$ $<$functional name 2’$>$
$<$coefficient 3’$>$ $<$functional name 3’$>$

See option user for the possible values of $<$functional name n$>$ and $<$functional name n’$>$. The weight of the HF exchange (HFx), if any, can be different for the density and the energy, and, in contrast to previous versions of Mrcc, must be specified also in the second block.

Table 4: Options for keyword redcost_exc. CS-NAF – complete MO space NAFs are used, NO – frozen natural orbitals are used, Can. – the NO space is augmented with canonical virtual orbitals, RS-NAF – restricted NO space NAFs are used. See Ref. for more details.
 Functional Description User LDA exchange functionals LDA Slater–Dirac exchange (local density approximation) [27, 144, 69] Yes LDA correlation functionals VWN1 functional I of Vosko, Wilk, and Nusair [153] Yes VWN2 functional II of Vosko, Wilk, and Nusair [153] Yes VWN3 functional III of Vosko, Wilk, and Nusair [153] Yes VWN4 functional IV of Vosko, Wilk, and Nusair [153] Yes VWN5 functional V of Vosko, Wilk, and Nusair [153] Yes PZ Perdew–Zunger 1981 correlation functional [117] Yes PW Perdew–Wang 1992 correlation functional [116] Yes GGA exchange functionals B88 Becke’s 1988 exchange functional [8] Yes PBEx functional of Perdew, Burke, and Ernzerhof [118] Yes PBEh 1988 revision of PBEx by Ernzerhof and Perdew [31] Yes PW91x Perdew–Wang 1991 exchange functional [114] Yes G96 exchange functional of Gill [118] Yes mPW91x modified PW91x functional of Adamo and Barone [2] Yes GGA correlation functionals LYP correlation functional of Lee, Yang, and Parr [86] Yes P86 Perdew’s 1986 correlation functional [120] Yes PBEc functional of Perdew, Burke, and Ernzerhof [118] Yes PW91c Perdew–Wang 1991 correlation functional [114] Yes GGA exchange-correlation functionals BLYP Becke’s 1988 exchange functional [8] and the correlation functional of Lee, Yang, and Parr (B88 + LYP) [86] No BP86 BP86 exchange-correlation functional (B88 + P86) [8, 120] No PBE exchange-correlation functional of Perdew, Burke, and Ernzerhof (PBEx + PBEc) [118] No PW91 Perdew and Wang 1991 exchange-correlation functional (PW91x + PW91c) [114] No HCTH120 HCTH120 exchange-correlation functional of Boese and co-workers [13] Yes HCTH147 HCTH147 exchange-correlation functional of Boese and co-workers [13] Yes HCTH407 HCTH407 exchange-correlation functional of Boese and Handy [14] Yes XLYP exchange-correlation functional of Xu and Goddard [162] Yes mPWLYP1w exchange-correlation functional of Dahlke and Truhlar optimized for water [22] Yes Hybrid GGA exchange-correlation functionals BHLYP Becke’s half-and-half exchange in combination with the LYP correlation functional (0.5 B88 + 0.5 HF exchange + LYP) [8, 86, 9] No B3LYP Becke’s three-parameter hybrid functional including the correlation functional of Lee, Yang, and Parr (0.08 LDA + 0.72 B88 + 0.2 HF exchange + 0.19 VWN5 + 0.81 LYP) [8, 27, 144, 153, 10, 86] No B3LYP3 Becke’s three-parameter hybrid functional including the correlation functional of Lee, Yang, and Parr (0.8 LDA + 0.72 B88 + 0.2 HF exchange + 0.19 VWN3 + 0.81 LYP) [8, 27, 144, 153, 10, 86, 146]. Note that this is equivalent to the B3LYP functional of the Gaussian package. Yes B3PW91 Becke’s three-parameter hybrid functional including the 1991 correlation functional of Perdew and Wang (0.08 LDA + 0.72 B88 + 0.2 HF exchange + 0.19 VWN5 + 0.81 PW91c) [8, 27, 144, 153, 10, 114] No B1LYP modified B3LYP functional of Adamo and Barone [1] Yes O3LYP modified B3LYP functional of Cohen and Handy [20] Yes B97 Becke’s 1997 exchange-correlation functional (including 0.1943 HF exchange) [11] Yes PBE0 hybrid functional of Perdew, Burke, and Ernzerhof (0.75 PBEx + 0.25 HF exchange + PBEc) [118, 119] No X3LYP hybrid functional of Xu and Goddard [162] Yes Meta-GGA exchange functionals TPSSx exchange functional of Tao, Perdew, Staroverov, and Scuseria [151] Yes revTPSSx revised TPSS exchange of Perdew et al. [115] Yes SCANx exchange functional of Sun, Ruzsinszky, and Perdew [149] Yes Meta-GGA correlation functionals B95 Becke’s 1995 correlation functional [6] Yes TPSSc correlation functional of Tao, Perdew, Staroverov, and Scuseria [151] Yes revTPSSc revised TPSS correlation of Perdew et al. [115] Yes SCANc correlation functional of Sun, Ruzsinszky, and Perdew [149] Yes Meta-GGA exchange-correlation functionals TPSS exchange-correlation functional of Tao, Perdew, Staroverov, and Scuseria [151] No revTPSS revised TPSS functional of Perdew et al. [115] No M06-L 2006 exchange-correlation functional of Zhao and Truhlar [167, 166] No B97M-V exchange-correlation functional of Mardirossian and Head-Gordon [91] Yes SCAN exchange-correlation functional of Sun, Ruzsinszky, and Perdew [149] No Hybrid meta-GGA exchange-correlation functionals M06-2X 29-parameter exchange-correlation functional of Zhao and Truhlar including 0.54 HF exchange [167] No M08-HX 47-parameter exchange-correlation functional of Zhao and Truhlar including 0.5223 HF exchange [168] Yes M08-SO 44-parameter exchange-correlation functional of Zhao and Truhlar including 0.5679 HF exchange [168] Yes TPSSh hybrid version of TPSS including 0.1 HF exchange [145] Yes revTPSSh revised TPSSh of Csonka, Perdew, and Ruzsinszky including 0.1 HF exchange [145, 21] Yes mPW1B95 mixture of mPW91x and B95 by Zhao and Truhlar [164] Yes PW6B95 mixture of PW91x and B95 by Zhao and Truhlar [165] Yes SCAN0 hybrid version of SCAN including 0.25 HF exchange [149, 60] No Double hybrid functionals B2PLYP Grimme’s two-parameter double hybrid functional including MP2 correction (0.47 B88 + 0.53 HF exchange + 0.73 LYP + 0.27 MP2 correlation) [47] No B2GPPLYP two-parameter double hybrid functional including MP2 correction of Martin and co-workers (0.35 B88 + 0.65 HF exchange + 0.64 LYP + 0.36 MP2 correlation) [64] No DSDPBEP86 dispersion corrected, spin-component scaled double hybrid functional of Kozuch and Martin (0.30 PBEx + 0.70 HF exchange + 0.43 P86 + 0.53 MP2 antiparallel-spin correlation + 0.25 MP2 parallel-spin correlation) [70, 71]. Note that the dispersion correction is only included if the -D3 postfix is added (see the note below). No DSDPBEhB95 dispersion corrected, spin-component scaled double hybrid functional of Kozuch and Martin (0.34 PBEh + 0.66 HF exchange + 0.55 B95 + 0.47 MP2 antiparallel-spin correlation + 0.09 MP2 parallel-spin correlation) [71]. Note that the dispersion correction is only included if the -D3 postfix is added (see the note below). No XYG3 double hybrid functional of Zhang, Xu, and Goddard (0.2107 B88 - 0.014 LDA + 0.8033 HF exchange + 0.6789 LYP + 0.3211 MP2 correlation evaluated with B2LYP orbitals) [163, 42] No SCAN0-2 SCAN-based double-hybrid of Hui and Chai (0.793701 HF exchange + 0.206299 SCANx + 0.5 SCANc + 0.5 MP2 correlation [149, 60] No dRPA75 the dual-hybrid random phase approximation (dRPA75) method of Mezei et al. [99]. The KS orbitals are obtained with the “0.25 PBEx + 0.75 HF exchange + PBEc” functional, while the energy is calculated using the “0.25 PBEx + 0.75 HF exchange + dRPA correlation” expression. Dispersion correction [17] can be included if the -D3 postfix is added. No SCS-dRPA75 the spin-component scaled dual-hybrid random phase approximation (SCS-dRPA75) method of Mezei et al. [99, 100]. The KS orbitals are obtained with the “0.25 PBEx + 0.75 HF exchange + PBEc” functional, while the energy is calculated using the “0.25 PBEx + 0.75 HF exchange + 1.5 dRPA antiparallel-spin correlation + 0.5 dRPA parallel-spin correlation” expression. No van der Waals density functionals VV10NL the nonlocal part (the $\beta N$ term is ignored) of the 2010 van der Waals density functional of Vydrov and Van Voorhis [154], both self-consistent and non-self-consistent implementations are available, see also the the comment below Yes
Default:

dft=off

Notes:
1.

The built-in functionals implemented in Mrcc were obtained from the Density Functional Repository [38, 148]. Other functionals are available via the Libxc interface [92, 59] and require the Libxc library, see Sect. 7.2 for the installation of Libxc.

2.

Empirical dispersion corrections can be calculated for particular functionals and also for the HF energy using the DFT-D3 approach of Grimme and co-workers [44, 45] by attaching the -D3 postfix to the corresponding options: BLYP-D3, BHLYP-D3, B3LYP-D3, B3PW91-D3, BP86-D3, PBE-D3, PBE0-D3, HCTH120-D3, B2PLYP-D3, mPW1B95-D3, TPSS-D3, TPSSh-D3, B2GPPLYP-D3, DSDPBEP86-D3, DSDPBEhB95-D3, dRPA75-D3,

3.

For a simple DFT calculation (i.e., without subsequent correlation calculations) the value of keyword calc can be SCF, HF, RHF, or UHF. Note that you do not need to set its value since it is set to SCF by default. Alternatively, you can select the DFT functional using keyword calc, and in this case you do not have to set keyword dft (see the description of calc).

4.

For a correlated calculation with KS orbitals you should select the functional with this keyword, and the value of keyword calc must be set to the desired correlation method. Note that you can also accelerate the post-KS calculation using local correlation schemes (e.g., local dRPA). See the examples below.

5.

For a correlated calculation with KS orbitals (excluding calculations with double hybrid functionals) the HF energy computed with KS orbitals is used as reference energy.

6.

For the B2PLYP, B2GPPLYP, DSDPBEP86, DSDPBEhB95, dRPA75, etc. double hybrid functionals as well as for user-defined double hybrid functionals including MP2 (SCS-MP2), dRPA, etc. correlation calc is automatically set to MP2, dRPA, etc. Note that you can accelerate the MP2, dRPA, …part of a double hybrid DFT calculation for large molecules using local correlation approaches. For the built-in double hybrid functionals just add the “L” prefix, while for the user-defined functionals set localcc=on. See the examples below.

7.

The DSDPBEP86, DSDPBEhB95, and dRPA75 functionals use special parameters for the calculation of the D3 correction which are read by the DFT-D3 program from the .dftd3par.\$HOST file located in your home directory. This file will be created by the program, but you must be sure that the program is able to access your home directory. Also note that, if you already have this file in your home, it will be overwritten, so please do not forget to save it before executing Mrcc.

8.

For the VV10 van der Waals functional you can modify parameters $b$ and $C$ (see Ref. 154) if it is used with the user or userd options. For that purpose the two parameters should be specified after the VV10NL flag separated by spaces, see the example below. If the parameters are not set, those of Ref. 154 will be used.

Examples:
1.

To perform a DFT calculation with the B3LYP functional give dft=B3LYP or calc=B3LYP

2.

The B3LYP functional can also be defined using the user option as
calc=scf
dft=user
5
0.08 LDA
0.72 B88
0.20 HFx
0.19 VWN5
0.81 LYP

3.

The B2PLYP double-hybrid functional can also be defined using the user option as
calc=scf
dft=user
4
0.47 B88
0.73 LYP
0.53 HFx
0.27 MP2

4.

The DSDPBEP86 double-hybrid functional can also be defined using the user option as
calc=SCF
dft=user
5
0.30 PBEx
0.43 P86
0.70 HFx
0.53 MP2s
0.25 MP2t

5.

SOSEX calculation with Kohn–Sham orbitals calculated with the LDA exchange functional:
calc=SOSEX
dft=LDA

6.

To perform a DFT calculation with the B2PLYP double-hybrid functional and add the D3 dispersion correction set dft=B2PLYP-D3 or calc=B2PLYP-D3

7.

B2PLYP calculation, the MP2 contribution is evaluated using local MP2 approximation:
calc=LB2PLYP

8.

User-defined functional, different functionals are used for the calculation of the density (0.25 PBEx + 0.75 HF exchange + PBEc) and the energy (0.50 PBEx + 0.50 HF exchange + MP2 correlation).
dft=userd
3
0.75 HFx
0.25 PBEx
1.00 PBEc
3
0.50 HFx
0.50 PBEx
1.00 MP2

9.

The dRPA75 dual-hybrid functional can also be defined using the userd option as
dft=userd
3
0.75 HFx
0.25 PBEx
1.00 PBEc
3
0.75 HFx
0.25 PBEx
1.00 dRPA

10.

Local dRPA calculation with Kohn–Sham orbitals calculated with the PBE functional:
calc=LdRPA
dft=PBE

11.

To perform a DFT calculation with the B3LYP functional using its Libxc implementation set calc=HYB_GGA_XC_B3LYP5

12.

The B3LYP functional can also be defined using the user option and the functionals implemented in the Libxc library as
dft=user
5
0.08 LDA_X
0.72 GGA_X_B88
0.20 HFx
0.19 LDA_C_VWN
0.81 GGA_C_LYP

13.

DFT calculation with a user-defined PBE0-VV10 functional. Parameter $b$ of VV10 is modified, while for $C$ its default value, 0.0093, is used. If you do not want to modify either $b$ or $C$, simply drop the two numbers for the VV10NL entry.
dft=userd
3
0.75 PBEx
0.25 HFx
1.00 PBEc
4
0.75 PBEx
0.25 HFx
1.00 PBEc
1.00 VV10NL 8.0 0.0093

diag

Type of diagonalization algorithm used for the CI and LR-CC calculations.

Options:
david

Standard Davidson diagonalization

olsen

Another algorithm proposed by Olsen using only two expansion vectors (see Refs. 111, 81, and 82), useful for very large CI/LR-CC vectors

follow

Davidson diagonalization with root-following, recommended for excited-state calculations if the initial guess is given manually or the calculation is restarted

Default:

diag=david

Example:

for root-following type diag=follow

Radius of atom domains for the local correlation method of Ref. 138 (localcc=2013). For each localized MO (LMO), using the Boughton–Pulay procedure [16], we assign those atoms to the LMO on which it is localized. Then, for each LMO an atom domain is constructed in two steps, the LMO is called the central LMO of the domain. In the first step, those atoms are included in the domain whose distance from the atoms assigned to the central LMO is smaller than domrad. In the second step, those LMOs are identified which are localized on the atoms selected in the first step, and the domain is extended to include all atoms assigned to these LMOs.

Options:
$<$any positive real number$>$

In the first step of the construction of atom domains all atoms whose distance from the atoms assigned to the central LMO is smaller than this number (in bohr) will be included in the domain.

inf

Infinite radius will be applied, i.e., there is only one atom domain including all atoms.

Default:

Note:

To use the local CC methods as defined in Ref. 135 set domrad=inf, that is, use only one atom domain including all atoms.

Example:

to set a threshold of 12.0 bohr type domrad=12.0

drpaalg

Specifies the type of the algorithm for the solution of the dRPA equations or the calculation of SOS-MP2 energies. See Ref. 85 for more details.

Options:
fit

The algorithm of Ref. 58 will be used, the fitting of integral lists will be performed before the dRPA iterations (SOS-MP2 calculation).

nofit

The algorithm of Ref. 85 will be executed, the fitting of the integrals is not performed. This algorithm is efficient for large molecules.

plasmon

The dRPA correlation energy is calculated using the plasmon formula.

auto

The algorithm is automatically selected on the basis of the size of the molecule (canonical dRPA) or the HOMO-LUMO gap (local dRPA).

Notes:
1.

For SOSEX calculations drpaalg=fit is the only option, which is forced by the program.

2.

For canonical dRPA the algorithm using the plasmon formula scales as $N^{6}$, it is only competitive for smaller molecules but inefficient for bigger ones. It avoids, however, the problems of the other algorithms, that is, convergence problems and unphysical solutions. Thus, it is useful for testing.

3.

For local dRPA drpaalg=plasmon is also linear scaling but typically 2- to 4-times slower than drpaalg=fit. It is advantageous for the aforementioned reasons. If drpaalg=auto, the plasmon formula-based algorithm is executed if the HOMO-LUMO gap is lower than 0.05 E${}_{h}$.

Default:

drpaalg=fit and drpaalg=auto for canonical and local dRPA, respectively.

Example:

to set the second option give drpaalg=nofit

dual

Activates dual basis set calculations. For these calculations two basis sets must be specified: a smaller one by keyword basis_sm (see the description of the keyword) and a bigger basis defined by keyword basis. The energy evaluated with the bigger basis set is estimated from a small-basis calculation. See Ref. for more details.

Options:

off No dual basis calculation.

on Dual basis set calculation for conventional SCF and correlated methods. First, an SCF calculation will be performed using the small basis set. Second, one iteration of a SCF calculation is carried out with the large basis, and the energy is extrapolated using a first-order formula. If a correlation calculation is requested, the orbitals obtained in the second SCF step will be used for that purpose.

e1 Dual basis set embedding, Ansatz 1 of Ref. . A Huzinaga-embedding calculation is performed with the small basis set. The steps of the large-basis Huzinaga-embedding calculation are non-iterative. See also the description of keyword embed.

e2 Dual basis set embedding, Ansatz 2 of Ref. . Similar to e1, but there is also an iterative step with the large basis.

Default:

dual=off

Examples:
1.

To perform a dual basis set PBE calculation with the cc-pVTZ and cc-pVDZ basis sets you need:
basis=cc-pVTZ
basis_sm=cc-pVDZ
dual=on
calc=PBE

2.

Dual basis set PBE-in-LDA calculation with Ansatz 1 using the cc-pVTZ and cc-pVDZ basis sets as large and small basis set, respectively; atoms 1 to 5 are included in the embedded subsystem:
basis=cc-pVTZ
basis_sm=cc-pVDZ
calc=PBE
dual=e1
embed=Huzinaga
1-5
LDA
0

3.

Dual basis set PBE-in-LDA calculation with Ansatz 1 using a mixed large basis (cc-pVTZ for atoms 1 to 5, cc-pVDZ otherwise) and the cc-pVDZ basis sets as the small basis set; atoms 1 to 5 are included in the embedded subsystem:
basis=embed
cc-pVTZ
cc-pVDZ
basis_sm=cc-pVDZ
calc=PBE
dual=e1
embed=Huzinaga
1-5
LDA
0

ecp

Specifies the effective core potential (ECP) used in all calculations. By default the ECPs are taken from the files named by the chemical symbol of the elements, which can be found in the BASIS directory created at the installation. The ECPs are stored in the format used by the Cfour package. In addition to the ECPs provided by default, any ECP can be used by adding it to the corresponding files in the BASIS directory. Alternatively, you can also specify your own ECP in the file GENBAS which must be copied to the directory where Mrcc is executed.

Options:
none

No ECPs will be used.

auto

The ECPs will be automatically selected: no ECP will be used for atoms with all-electron basis sets, while the ECP adequate for the basis set of the atom will be selected otherwise.

$<$ECP label$>$

If the same ECP is used for all atoms, the label of the ECP can be given here.

atomtype

If different ECPs are used or no ECP is used for particular atoms, but the atoms of the same type are treated in the same way, ecp=atomtype should be given, and the user must specify the ECP for each atomtype (for which an ECP is used) in the subsequent lines as $<$atomic symbol$>$:$<$ECP label$>$ .

special

In the general case, if different ECPs are used for each atom, then one should give ecp=special and specify the ECP for each atom in the subsequent lines by giving the label of the corresponding ECP (or none if no ECP is used for that atom) in the order the atoms appear at the specification of the geometry.

Notes:
• 1.

By default the following ECP are available for elements Na to Rn in Mrcc:

• the LANL2DZ ECP’s of Hay and Wadt [52, 155, 51]: LANL2DZ-ECP-10, LANL2DZ-ECP-18, LANL2DZ-ECP-28, LANL2DZ-ECP-36, LANL2DZ-ECP-46, LANL2DZ-ECP-60, LANL2DZ-ECP-68, LANL2DZ-ECP-78

• the Stuttgart–Köln ECPs for the def2 basis sets [4, 65, 87]: def2-ECP-28, def2-ECP-46, def2-ECP-60

• the Stuttgart–Köln multiconfiguration Dirac–Hartree–Fock-adjusted ECPs [97, 98, 126, 125, 121, 33, 34]: MCDHF-ECP-10, MCDHF-ECP-28, MCDHF-ECP-60

Please note that some of the above ECPs are not available for all elements.

• 2.

• 3.

If you use your own ECPs, these must be copied to the end of the corresponding file in the BASIS directory. Alternatively, you can also create a file called GENBAS in the directory where Mrcc is executed, and then you should copy your ECPs to that file.

• 4.

The labels of the ECPs must be identical to those used in the BASIS/* files (or the GENBAS file). For the default ECPs just type the name of the ECPs as given above, e.g., LANL2DZ-ECP-10, def2-ECP-28, etc. If you employ non-default ECPs, you can use any label.

Default:

ecp=auto

Examples:
1.

To use the MCDHF-ECP-10 pseudopotential for all atoms the input must include ecp=MCDHF-ECP-10

2.

Consider the PbO molecule and use the def2-SVP basis set for both elements as well as the def2-ECP-60 pseudopotential for Pb. The following inputs are equivalent.
Input 1:
basis=def2-SVP
geom
Pb
O 1 R

R=1.921813

Input 2:
basis=def2-SVP
ecp=atomtype
Pb:def2-ECP-60
geom
Pb
O 1 R

R=1.921813

Input 3:
basis=def2-SVP
ecp=special
def2-ECP-60
none
geom
Pb
O 1 R

R=1.921813

edisp

This keyword controls the calculation of empirical dispersion corrections for DFT and HF calculations using the DFT-D3 approach of Grimme and co-workers [44, 45]. The corrections are evaluated by the dftd3 program of the latter authors, which is available at http://www.thch.uni-bonn.de/tc/ and interfaced to Mrcc. You need to separately install this code and add the directory where the dftd3 executable is located to your PATH environmental variable.

Options:
off

No dispersion correction will be computed.

auto

The dispersion correction will be automatically evaluated to the KS or HF energy. Note that it is only possible for particular functionals listed in the description of keyword dft (and the HF method). For these methods, however, you can also turn on the calculations of the dispersion corrections by attaching the -D3 postfix to the corresponding options, e.g., as BLYP-D3, B3LYP-D3, B2PLYP-D3, etc. (see the description of keyword dft).

$<$any options of the dftd3 program$>$

You can directly give any options of the dftd3 code. The options will be passed over to dftd3 without any consistency check, the user should take care of the compatibility of these options with the calculation performed by Mrcc. Note that the coordinate file name must not be specified here, the coordinates will be taken from the COORD.xyz file generated by Mrcc.

Note:

If edisp=auto or the -D3 postfix is added to the corresponding options, the empirical dispersion correction is by default evaluated with the Becke and Johnson (BJ) damping function [45].

Default:

edisp=off

Example:
1.

to calculate the D3 dispersion correction including BJ damping to the B3LYP energy give calc=B3LYP-D3

2.

to calculate the D3 dispersion correction to the B3LYP energy without the BJ damping the input should include:
calc=B3LYP
edisp=-func b3-lyp -zero

embed

This keyword controls DFT embedding calculations. Currently it is only available for closed-shell systems using density-fitting.

Options:
off

No embedding.

Huzinaga

The Huzinaga-equation-based embedding approach [62] will be used. The embedded atoms and the low-level DFT approach (or HF) used for the embedding must be specified and the number of embedded orbitals can be given in the subsequent lines as follows. The embedded atoms must be given by their serial numbers in the first line as $$,$$,…,$$-$$,…, where $n_{i}$’s are the serial numbers of the atoms. Serial numbers separated by dash mean that $$ through $$ are embedded atoms. Note that the numbering of the atoms must be identical to that used in the Z-matrix or Cartesian coordinate specification, but dummy atoms must be excluded. The low-level DFT (LDA, GGA, or hybrid) or HF approach must be specified in the second line using the corresponding option of keyword dft (for HF simply HF). The high-level method (any DFT, HF, or any correlation method) for the active region should be specified by the keyword calc (or keywords calc and dft if a KS reference is used in a correlation calculation). In the third line an integer should be given which is the number of the embedded orbitals or zero if the latter is determined automatically. In addition, the algorithm for the selection of the orbitals can also be given in the third line after the integer. The options are aMul and bopu, which mean the Mulliken population- and the Boughton–Pulay algorithm-based schemes of Refs. 62 and , respectively. For aMul you can also specify the Mulliken population threshold (see the appendix of Ref. ). The threshold for the bopu algorithm can be controlled by the bpcompo keyword. The default is aMul with a threshold of 0.3. See also examples below.

project

The projector-based embedding approach of Manby and co-workers [90] will be used. The embedded atoms and the low-level DFT approach can be specified as described above.

Default:

embed=off

Examples:
1.

CCSD(T)-in-PBE embedding with the Huzinaga-equation-based approach, atoms 1 and 2 are included in the embedded region, the number of the embedded orbitals is determined automatically, the aMul algorithm is used for the selection of the orbitals with the default threshold:
calc=DF-CCSD(T)
embed=huzinaga
1-2
PBE
0

2.

The same as example 1, but the aMul algorithm is used with a threshold of 0.25:
calc=DF-CCSD(T)
embed=huzinaga
1-2
PBE
0 amul 0.25

3.

The same as example 1, but the bopu algorithm is used:
calc=DF-CCSD(T)
embed=huzinaga
1-2
PBE
0 bopu

4.

Same as example 1, but a mixed basis set is used, cc-pVDZ for the environment and cc-pVTZ for the embedded subsystem:
calc=DF-CCSD(T)
basis=embed
cc-pVDZ
cc-pVTZ
embed=huzinaga
1-2
PBE
0

5.

dRPA@PBE-in-PBE embedding with the Huzinaga-equation-based approach, atoms 1, 2, 3, and 5 as well as 10 orbitals are included in the active region:
calc=dRPA
dft=PBE
embed=huzinaga
1-3,5
PBE
10

epert

Use this option to add an external perturbation to the Hamiltonian, e.g., an external electric dipole field.

Options:
none

$<$any integer in the [0,9] interval$>$

the number of the operators added to the Hamiltonian. The operators and the corresponding coefficients (in a.u.) should be specified in the subsequent lines as follows:
$<$operator 1$>$ $<$coefficient 1$>$
$<$operator 2$>$ $<$coefficient 2$>$
$<$operator 3$>$ $<$coefficient 3$>$

where the operator can be x, y, z, xx, yy, zz, xy, xz, yz, xxx, xxy, xxz, xyy, xyz, xzz, yyy, yyz, yzz, zzz.

Note:

The symmetry of the perturbation is not taken care of automatically. If the perturbation lowers the symmetry of the system, you must change the computational point group (keyword cmpgrp) or turn off symmetry (symm=off).

Default:

epert=none

Example:

to add the $\hat{y}$ and $\hat{z}$ dipole length operators to the Hamiltonian with coefficients 0.01 and 0.001 a.u., respectively, the input should include the following lines
epert=2
y 0.01
z 0.001

eps

Threshold for the cumulative populations of MP2 natural orbitals (NOs) or optimized virtual orbitals (OVOs), to be used together with keyword ovirt. The cumulative population for an MO is calculated by summing up the occupation number of that particular MO and all the MOs with larger occupation numbers, and then this number is divided by the number of electrons. See Ref. 136 for more details.

Options:
$<$any real number in the [0,1] interval$>$

Virtual orbitals with cumulative populations of higher than this number will be dropped.

Default:

eps=0.975

Example:

to set a threshold of 0.95 type eps=0.95

Sets the radius of local fitting domains for the exchange contribution in direct density-fitting SCF calculations [130, 105]. In direct DF-SCF calculations, in each iteration step, the MOs are localized. For each localized MO Löwdin atomic charges are computed, and all atoms are selected which have a charge greater than 0.05. All further atoms will be included in the fitting domain of the MO for which the electron repulsion integrals including the corresponding AOs and the basis functions residing on the atoms selected in the first step are estimated to be greater than a threshold.

Options:
$<$any positive real number$>$

Threshold for the integrals (in E${}_{h}$).

0

A threshold of zero will be applied, i.e., a conventional direct DF-SCF calculation will be executed.

Notes:
1.

Local fitting domains are currently available only for RHF and RKS wave functions.

2.

For average organic molecules with localized electronic structure excrad=1.0 is a good choice. For more complicated systems other thresholds may be necessary. For excrad=1.0 excrad_fin is set to $10^{-3}$, which is fine with basis sets excluding diffuse functions. For basis sets with diffuse functions excrad_fin=1e-5 or tighter is recommended.

Default:

Example:

to set a threshold of 1.0 E${}_{h}$ type excrad=1.0

In density-fitting SCF calculations, if excrad and excrad_fin differ, an extra iteration is performed to get an accurate SCF energy. excrad_fin specifies the radius of local fitting domains for the exchange contribution in this iteration step. See also notes for keyword excrad.

Options:

See the description of keyword excrad.

Default:

Example:

to avoid the use of local fitting domains in the extra iteration step give excrad_fin=0.0

freq

Requests harmonic vibrational frequency calculation (numerical). Ideal gas thermodynamic properties will also be evaluated in the rigid-rotor harmonic-oscillator approximation. If geometry optimization is also carried out, i.e., gopt $\neq$ off, the frequencies are calculated at the optimized geometry.

Options:

on or off

Note:

You should also set this keyword if you are interested thermodynamic properties of atoms.

Default:

freq=off

Example:

for a frequency calculation set freq=on

gauss

Specifies whether spherical harmonic or Cartesian Gaussian basis functions will be used.

Options:
spher

Spherical harmonic Gaussians will be used

cart

Cartesian Gaussians will be used

Notes:
1.

For calculations using the density fitting (DF) approximation, if intalg=os or intalg=auto, the Coulomb integrals are evaluated by algorithms [3, 139] which only enable the use of spherical harmonic Gaussians. Consequently, Cartesian Gaussians are only available with intalg=rys in DF calculations (see the description of keyword intalg).

2.

The derivative integrals are evaluated by the solid-harmonic Hermite scheme [133] (see the description of keyword intalg), consequently, differentiated integrals, and thus energy derivatives cannot be evaluated with Cartesian Gaussian basis sets.

Default:

gauss=spher

Example:

for Cartesian Gaussians the user should set gauss=cart

geom

Specifies the format of molecular geometry. The geometry must be given in the corresponding format in the subsequent lines.

Options:
zmat

Usual Z-matrix format. In the Z-matrix the geometrical parameters can only be specified as variables, and the variables must be defined after the matrix, following a blank line. Another blank line is required after the variables. This Z-matrix format is compatible to that of Cfour and nearly compatible to that for Gaussian and Molpro. Z-matrices can be generated by Molden (see also Sect. 14.1), then the Gaussian-style Z-matrix format must be chosen. The symbol for dummy atoms is “X”.

xyz

Cartesian coordinates in xyz format, that is, the number of atoms, a blank line, then for each atom the atomic symbol or atomic number and the $x$, $y$, and $z$ components of Cartesian coordinates. Cartesian coordinates in xyz format can also be generated by Molden (see also Sect. 14.1).

tmol

Cartesian coordinates in a format similar to that used by the Turbomole package, that is, the number of atoms, a blank line, then for each atom the $x$, $y$, and $z$ components of Cartesian coordinates and the atomic symbol or atomic number.

mol

Cartesian coordinates and connectivity in .mol format, that is, the number of atoms and number of bonds in the first line, then for each atom the $x$, $y$, and $z$ components of Cartesian coordinates and the atomic symbol, then for each bond the serial number of the atoms connected by the bond and the type of the bond (1 for single bond, etc.). This geometry specifications is needed if the specified method requires the connectivity.

Note:

For the use of ghost atoms see the description of keyword ghost.

Default:

geom=zmat, which is equivalent to geom, i.e., if it is not specified whether the geometry is supplied in Z-matrix format or in other formats, Z-matrix format is supposed. Nevertheless, the coordinates must be given in the subsequent lines in any case.

Examples:

the following five geometry inputs for H${}_{2}$O${}_{2}$ are equivalent

1.

Z-matrix format, bond lengths in Å:
geom
H
O 1 R1
O 2 R2 1 A
H 3 R1 2 A 1 D

R1=0.967
R2=1.456
A=102.32
D=115.89

2.

xyz format, coordinates in bohr, atoms are specified by atomic symbols:
unit=bohr
geom=xyz
4

H 0.00000000 0.00000000 0.00000000
O 1.82736517 0.00000000 0.00000000
O 2.41444411 2.68807873 0.00000000
H 3.25922198 2.90267673 1.60610134

3.

xyz format, coordinates in bohr, atoms are specified by atomic numbers:
unit=bohr
geom=xyz
4

1 0.00000000 0.00000000 0.00000000
8 1.82736517 0.00000000 0.00000000
8 2.41444411 2.68807873 0.00000000
1 3.25922198 2.90267673 1.60610134

4.

Turbomole format, coordinates in bohr, atoms are specified by atomic symbols:
unit=bohr
geom=tmol
4

0.00000000 0.00000000 0.00000000 H
1.82736517 0.00000000 0.00000000 O
2.41444411 2.68807873 0.00000000 O
3.25922198 2.90267673 1.60610134 H

4.

.mol format, coordinates in bohr:
unit=bohr
geom=mol
4 3
0.00000000 0.00000000 0.00000000 H
1.82736517 0.00000000 0.00000000 O
2.41444411 2.68807873 0.00000000 O
3.25922198 2.90267673 1.60610134 H
1 2 1
2 3 1
3 4 1

ghost

Ghost atoms can be specified using this keyword, e.g., for the purpose of basis set superposition error (BSSE) calculations with the counterpoise method.

Options:
none

There are no ghost atoms.

serialno

Using this option one can select the ghost atoms specifying their serial numbers. The latter should be given in the subsequent line as $$,$$,…,$$-$$,…, where $n_{i}$’s are the serial numbers of the atoms. Serial numbers separated by dash mean that $$ through $$ are ghost atoms. Note that the numbering of the atoms must be identical to that used in the Z-matrix or Cartesian coordinate specification, but dummy atoms must be excluded.

Default:

ghost=none

Examples:
1.

Rectangular HF dimer, the atoms of the second HF molecule are ghost atoms:
geom
H
F 1 R1
H 2 R2 1 A
F 3 R1 2 A 1 D

R1=0.98000000
R2=2.00000000
A=90.00000000
D=0.000000000

ghost=serialno
3-4

2.

Ammonia, the third hydrogen is a ghost atom (note that the serial number of the hydrogen is 4 instead of 5 because of the dummy atom:
geom
X
N 1 R
H 2 NH 1 AL
H 2 NH 1 AL 3 A
H 2 NH 1 AL 3 B

R=1.00000000
NH=1.01000000
AL=115.40000000
A=120.00000000
B=-120.00000000

ghost=serialno
4

gopt

Requests geometry optimization. Currently only the full geometry optimization is supported, geometrical parameters cannot be frozen.

Options:
off

no geometry optimization.

full

full geometry optimization.

Notes:
1.

The coordinates in the MINP file are replaced by the converged ones at the end of the geometry optimization. The initial MINP file is saved as MINP.init.

2.

The Abelian symmetry of the molecule is utilized at the calculation of gradients and update of the coordinates, thus, the computational point group is preserved during the optimization.

Default:

gopt=off

Example:

to carry out a full geometry optimization set gopt=full

gtol

Threshold for automatic point group recognition. Two atoms will be considered symmetry-equivalent if the difference in any component of their Cartesian coordinates after the symmetry operation is less than $10^{-{\tt gtol}}$ bohr.

Options:

$<$any integer$>$

Default:

gtol=7

Example:

for a tolerance of $10^{-4}$ bohr give gtol=4

grdens

This keyword is useful for the analysis of HF, KS, or correlated (MP2, CI, CC, …) one-electron density and its derivatives. The one-electron density, its gradient, and Laplacian will be calculated on a grid used for DFT calculations (see keywords agrid and rgrid) and saved for external use. In the case of correlated calculations the densities are evaluated using the relaxed density matrices.

Options:
off

Densities are not evaluated.

on

The one-electron density and its derivatives are calculated in the grid points. These values together with the grid are written to the unformatted Fortran file DENSITY. If a correlation calculation is performed the densities calculated with the correlation method are stored in the DENSITY file, while the SCF densities are saved to the file DENSITY.SCF. For restricted orbitals the files use the following format:

$N$
$\mathbf{r}_{1}$ $w_{1}$ $\rho(\mathbf{r}_{1})$ $\nabla\rho(\mathbf{r}_{1})$ $\nabla^{2}\rho(\mathbf{r}_{1})$
$\mathbf{r}_{2}$ $w_{2}$ $\rho(\mathbf{r}_{2})$ $\nabla\rho(\mathbf{r}_{2})$ $\nabla^{2}\rho(\mathbf{r}_{2})$
⋮
$\mathbf{r}_{N}$ $w_{N}$ $\rho(\mathbf{r}_{N})$ $\nabla\rho(\mathbf{r}_{N})$ $\nabla^{2}\rho(\mathbf{r}_{N})$

where $N$ is the number of grid points, $\mathbf{r}_{i}$ = ($x_{i}$, $y_{i}$, $z_{i}$) is the Cartesian coordinate of grid point $i$ with $w_{i}$ as the corresponding weight, and $\rho(\mathbf{r}_{i})$ is the density in point $i$. For unrestricted calculations the corresponding $\alpha$ and $\beta$ quantities are stored separately, and the lines of the files change as
$\mathbf{r}_{i}$ $w_{i}$ $\rho_{\alpha}(\mathbf{r}_{i})$ $\rho_{\beta}(\mathbf{r}_{i})$ $\nabla\rho_{\alpha}(\mathbf{r}_{i})$ $\nabla\rho_{\beta}(\mathbf{r}_{i})$ $\nabla^{2}\rho_{\alpha}(\mathbf{r}_{i})$ $\nabla^{2}\rho_{\beta}(\mathbf{r}_{i})$

Default:

grdens=off

Example:

to save the densities give grdens=on

grtol

The keyword controls the fineness of the angular and radial integration grids employed in DFT calculations. The tolerance for the accuracy of angular integrals will be $10^{-{\tt grtol}}$, while the number of radial grid points increases linearly with grtol. See also the description of keywords agrid and rgrid.

Options:

$<$any positive integer$>$

Note:

meta-GGA functionals or molecules with special bonding characteristics may require larger integration grids, and it is recommended to run test calculations to decide if the default value of grtol is sufficient.

Default:

grtol=10

Example:

for a fine integration grid give grtol=12

hamilton

Specifies what type of Hamiltonian is used in relativistic calculations. This keyword has only effect if iface=Dirac.

Options:
X2Cmmf

exact 2-component molecular-mean-field Hamiltonian [143]

DC

other types of relativistic Hamiltonians such as the full Dirac–Coulomb Hamiltonian or the exact 2-component Hamiltonian

Default:

hamilton=DC

Example:

if you use the exact 2-component molecular-mean-field Hamiltonian, set hamilton=X2Cmmf

iface

Specifies whether Mrcc is used together with another program system. In this case the transformed MO integrals are calculated by that program and not by Mrcc. See Sect. 5 for the description of various interfaces.

Options:
none

Transformed MO integrals are calculated by Mrcc.

Cfour

Mrcc is interfaced to Cfour.

Columbus

Mrcc is interfaced to Columbus.

Dirac

Mrcc is interfaced to Dirac.

Molpro

Mrcc is interfaced to Molpro.

Notes:
1.

If you use Mrcc together with Cfour or Molpro, you do not need to use this keyword. The Mrcc input file is automatically written and Mrcc is automatically called by these program systems. The user is not required to write the Mrcc input file, most of the features of Mrcc can be controlled from the input files of these programs. With Cfour the user has the option to turn off the automatic construction of the Mrcc input file by giving INPUT_MRCC=OFF in the Cfour input file ZMAT. In the latter case one should use this keyword.

2.

If you use Mrcc together with Columbus or Dirac, this keyword must be always given.

Default:

iface=none, that is, all calculations will be performed by Mrcc.

Example:

to carry out four-component relativistic calculations using the Dirac interface give iface=Dirac

intalg

Specifies the algorithm used for the evaluation of two-electron integrals over primitive Gaussian-type orbitals.

Options:
os

The $(\mathbf{e}\mathbf{0}|\mathbf{f}\mathbf{0})$ integrals are evaluated by the Obara–Saika procedure using the vertical and transfer recurrence relations [110, 89].

rys

The $(\mathbf{e}\mathbf{0}|\mathbf{f}\mathbf{0})$ integrals are evaluated by the Rys quadrature scheme [67, 89, 35].

auto

Depending on the angular momenta the program automatically determines which of the two algorithms is executed. For integrals of low angular momentum functions the Rys procedure is used, while the Obara–Saika algorithm is executed otherwise.

herm

The integrals over contracted Gaussians are evaluated by the solid-harmonic Hermite scheme of Reine et al. [133].

Notes:
1.

For calculations using the density fitting (DF) approximation intalg=auto is equivalent intalg=os since the Obara–Saika algorithm is more efficient for any integrals.

2.

For DF methods option herm is not available.

3.

For DF methods, if intalg=os or intalg=auto the Coulomb integrals are evaluated by the algorithm of Ref. 139, which only enables the use of spherical harmonic Gaussians. Consequently, Cartesian Gaussians are only available with intalg=rys in DF calculations (see the description of keyword gauss).

4.

The derivative integrals are evaluated by the solid-harmonic Hermite scheme even if another option is used for the undifferentiated integrals. Consequently, differentiated integrals, and thus energy derivatives cannot be evaluated with Cartesian Gaussian basis sets.

Default:

intalg=auto

Example:

to use the Obara–Saika scheme for all angular momenta add intalg=os

itol

Threshold for integral calculation. Integrals less than $10^{-{\tt itol}}$ E${}_{h}$ will be neglected.

Options:

$<$any integer$>$

Default:

itol=max(10, scftol+4, scfdtol), but itol is changed to itol+1 if basis functions are dropped because of linear dependence (see keyword ovltol)

Example:

for an accuracy of $10^{-15}$ E${}_{h}$ one must give itol=15

laptol

Specifies the accuracy of the numerical Laplace transform in the (T) correlation energy term of local CCSD(T) calculations with localcc=2016 or 2018. See also the description of keyword talg and Ref. 104.

Options:
$<$any positive real number$>$

The (T) energy denominator will be approximated using its numerical Laplace transform. The number of quadrature points ($n_{\mathrm{q}}$), and hence the accuracy is determined by this number. The minimum value of $n_{\mathrm{q}}$ is also set by this number as $n_{\mathrm{q}}>|\mathrm{log}_{10}(\texttt{laptol})|$.

Default:

laptol=$10^{-2}$, but laptol=$10^{-3}$ is set if lcorthr=tight. See also the description of lcorthr for further details.

Example:

to use a threshold of $10^{-3}$ type laptol=1e-3

lccrest

Use this keyword to restart local CC calculations in the case of localcc=2015, localcc=2016, and localcc=2018, e.g., after power failure. For the restart with lccrest=on the LMP2 calculation and integral transformations have to be completed, and for the remaining domains the DFINT_AI.*, DFINT_AB.*, DFINT_IJ.*, ajb.*, 55.*, 56, localcc.restart, and VARS files are required. In the case of localcc=2016 or localcc=2018 and talg=lapl or topr the laplbas.* files are also needed.

If the loop for the extended domains has started but was not finish in the LMP2 calculation or in the LMP2 part of and LNO-CC calculation with localcc=2016 or localcc=2018, the calculation can be restarted with the lccrest=domain option. In this case all intermediate files located in the folder of execution have to be kept intact, and only the updated MINP file should be overwritten before the calculation is restarted. Note that the beginning of the computation up to the pair energy evaluation is repeated in the restarted run, and the loop over the extended domains will continue from the index of the first unfinished domain.

Options:

on, off, or domain.

Default:

lccrest=off

Example:

to restart a local CCSD(T) calculation set lccrest=on or lccrest=domain depending on the point of interruption

lcorthr

Controls the accuracy of local correlation calculations by setting the relevant thresholds: bp*, lnoepso, lnoepsv, laptol, naf_cor, osveps, wpairtol, spairtol (see also Refs. 85, 105, and 106 for details).

Options:
Loose

Relatively loose thresholds will be used, see Tables 2 and 3.

Normal

Default threshold set if localcc=2018. lcorthr=Normal is a synonym for lcorthr=Loose if localcc=2015 or 2016. See Tables 2 and 3.

Tight

Tight thresholds will be used, see Tables 2 and 3.

0

The truncation thresholds will be set so that the canonical energy be reproduced, it is only useful for testing.

Notes:
1.

The values of the thresholds controlled by lcorthr are summarized in Table 2 for the default localcc=2018 case and in Table 3 for the previous schemes.

2.

Expected accuracy. Using the Normal settings for local MP2 and CCSD(T), if localcc=2018, 1 kJ/mol (1 kcal/mol) average (maximum) errors are expected in energy differences even for relatively complicated or sizable systems (see Refs. 105 and 106). In the case of localcc=2015 considering local dRPA and dRPA related methods the expected average (maximum) errors for energy differences are 2 kJ/mol (2 kcal/mol) with Loose and 1 kJ/mol (1 kcal/mol) with Tight settings (see Ref. 85).

3.

For local MP2 naf_cor=off is set if localcc=2015. Keywords lnoepso, lnoepsv, and laptol are irrelevant for local MP2 calculations.

Table 2: Keyword values set by lcorthr if localcc=2018.
calc MP2 CC all
localcc 2018 2018
lcorthr Loose Normal Tight Loose Normal Tight 0
bpedo 0.9999 0.9999 0.99995 0.9999 0.9999 0.99995 1.0
wpairtol 3e-5 1e-5 3e-6 3e-5 1e-5 3e-6 0.0
lnoepso - - - 3e-5 1e-5 3e-6 0.0
lnoepsv - - - 3e-6 1e-6 3e-7 0.0
laptol - - - 1e-1 1e-2 1e-3 -
naf_cor 2e-3 2e-3 2e-3 1e-2 1e-2 1e-2 off
bpcompo 0.985 0.985 0.985 0.985 0.985 0.985 1.0
bpcompv 0.98 0.98 0.98 0.98 0.98 0.98 1.0
bppdo 0.999 0.999 0.999 0.999 0.999 0.999 1.0
bppdv 0.98 0.98 0.98 0.98 0.98 0.98 1.0
bpedv 0.995 0.995 0.995 0.995 0.995 0.995 1.0
osveps off off off off off off off
spairtol off off off off off off off
Table 3: Keyword values set by lcorthr for localcc=2015 or 2016. Note that identical settings are used for both lcorthr=Normal and lcorthr=Loose if localcc=2015 or 2016.
calc MP2, dRPA, CC MP2 CC all
localcc 2015 2016 2016
lcorthr Normal Tight Normal Tight Normal Tight 0
bpedo - - 0.9998 0.9999 0.9999 0.99995 1.0
wpairtol 1e-6 1e-7 1.5e-5 1e-5 1e-5 3e-6 0.0
spairtol 1e-4 1e-5 off off off off 0.0
osveps 1e-3 1e-4 off off off off 0.0
lnoepso 3e-5 1e-5 - - 2e-5 1e-5 0.0
lnoepsv 1e-6 3e-7 - - 1e-6 5e-7 0.0
laptol - - - - 1e-2 1e-3 -
naf_cor 1e-2 8e-3 2e-3 2e-3 1e-2 1e-2 off
bpcompo 0.985 0.985 0.985 0.985 0.985 0.985 1.0
bpcompv 0.98 0.98 - - - - 1.0
bppdo - - 0.999 0.999 0.999 0.999 1.0
bppdv - - 0.98 0.98 0.98 0.98 1.0
bpedv - - 0.995 0.995 0.995 0.995 1.0
Default:

lcorthr=Normal

Example:

to use tight thresholds set lcorthr=Tight

lmp2dens

Determines whether the MP2 density matrix fragments are calculated using the “correct” expressions derived for the general type of orbitals, or using the expressions derived for the canonical case (as described in Ref. 135).

Options:
on

The MP2 density matrix fragments are calculated using the correct, non-canonical expressions.

off

The MP2 density matrix fragments are calculated using the approximate canonical expressions (as defined in Ref. 135).

Notes:
1.

To reproduce the method described in Ref. 135 use lmp2dens=off.

2.

The use of lmp2dens=on is recommended since in this case the local CC energy can be corrected by the difference of the local MP2 energy and the approximate local MP2 energy calculated in the local interacting subspaces (see Total CC… energy + correction in the output). This correction usually improves the local CC energy.

Default:

lmp2dens=on

Example:

to use the canonical expressions give lmp2dens=off

lnoepso

Threshold for the occupation numbers of occupied local natural orbitals (LNOs) in the case of local correlation calculations, or for state-averaged MP2/CIS(D) occupied natural orbitals for reduced-cost excited-state calculations, see also keyword lnoepsv. See Ref. 138 as well as Refs. 96 and for more details.

Options:
$<$any real number in the [0,1] interval$>$

Orbitals with occupation numbers greater than 1-lnoepso will be frozen.

Note:

For default settings with options other than localcc=2018 and lcorthr=Normal, see the description of lcorthr.

Default:

lnoepso=1e-5 for local correlation, lnoepso=0 for excited-states

Example:

to set a threshold of $5\cdot 10^{-6}$ type lnoepso=5e-6

lnoepsv

Threshold for the occupation numbers of virtual local natural orbitals (LNOs) in the case of local correlation calculations, or for state-averaged MP2/CIS(D) virtual natural orbitals for reduced-cost excited-state calculations, see also keyword lnoepso. See Ref. 138 as well as Refs. 96 and for more details.

Options:
$<$any real number in the [0,1] interval$>$

Orbitals with occupation numbers smaller than this number will be dropped.

Note:

For default settings with options other than localcc=2018 and lcorthr=Normal, see the description of lcorthr.

Default:

lnoepsv=1e-6 for local correlation, lnoepsv=7.5e-5 for excited-states

Example:

to set a threshold of $5\cdot 10^{-7}$ type lnoepsv=5e-7

localcc

Specifies if local correlation calculation is performed. See Refs. 135, 138, 85, 105, 104, and 106 for more details.

Options:
off

No local correlation calculation is performed.

2013

The algorithm of Ref. 138 is used.

2015

An algorithm based on Ref. 85 is used.

2016

The local MP2 algorithm of Ref. 105 or the local CC algorithm of Ref. 104 is used.

2018

The local MP2 algorithm of Ref. 105 (with tighter settings) or the local CC algorithm of Ref. 106 is used.

on

Equivalent to localcc=2018 in the case of MP2 or CC, and equivalent to localcc=2015 otherwise

Note:

Local correlation methods can also be run if the prefix “L” (or “LNO-” in the case of CC methods) is added to the corresponding option of keyword calc, see the description of calc. Local SCS-MP2 (SOS-MP2) calculations can be run with calc=LSCS-MP2 (calc=LSOS-MP2).

Default:

localcc=off

Example:

for local correlation calculations give localcc=on

maxact

Maximum number of inactive labels. One can impose restrictions on the cluster operator using this keyword. The maximum number of virtual/occupied inactive labels on the singly, doubly, … excited clusters can be specified.

Options:

on or off. If maxact=on, the maximum number of virtual and occupied inactive labels must be specified in the subsequent line as an integer vector. The integers must be separated by spaces. The vector should contain as many elements as the excitation rank of the highest excitation in the cluster operator. The integers are maximum number of virtual/occupied inactive labels allowed on amplitudes of single, double, … excitations, respectively.

Default:

maxact=off

Example:

Suppose that we have up to quadruple excitations, and the single, double, triple, and quadruple excitations are allowed to have maximum of 1, 2, 2, and 1 inactive virtual and occupied labels, respectively. Then the input file should include the following lines:
maxact=on
1 2 2 1

maxdim

Maximum number of trial vectors in the Davidson procedure in the case of CIS, TDA, TD-HF, and TD-DFT calculations. The procedure will be restarted if the dimension of the reduced subspace reaches maxdim.

Options:

$<$any positive integer$>$

Default:

maxdim=max(100, 3*max(ntrip, nsing))

Example:

for a maximum of 50 expansion vectors set maxdim=50

maxex

Level of highest excitation included in the cluster operator in the case of MRCI/CC calculations. In an MR calculation all single, double (or higher) excitations out of the reference determinants are included in the cluster operator (see the description of keyword nacto), however, the very high excitations are frequently irrelevant. Using this option the latter can be dropped. If maxex is set to a positive integer $n$, only up to $n$-fold excitations will be included in the cluster operator. The excitation manifold can be further selected by imposing constrains on the number of active/inactive labels of the excitations (see keyword maxact). See Refs. 83 and 74 for more details.

Options:
0

The excitation manifold is not truncated.

$<$any positive integer$>$

The excitation manifold is truncated at $n$-fold excitations, see above.

Default:

maxex=0

Example:

to truncate the excitation manifold at triple excitations set maxex=3

mem

Specifies the core memory used.

Options:
$<$any positive integer$>$MB

The amount of memory to allocate is specified in megabytes

$<$any positive integer$>$GB

The amount of memory to allocate is specified in gigabytes

Default:

mem=256MB

Example:

to allocate 8 GB core memory the user should set mem=8GB

molden

Specifies whether input file for the Molden program and an xyz-file containing the Cartesian coordinates are written (see also Sect. 14).

Options:
on

Cartesian coordinates, basis set information, and MO coefficients are saved to file MOLDEN. This file can be opened by Molden and used to visualize the structure of the molecule and the MOs. In addition, Cartesian coordinates are also written to file COORD.xyz in xyz (XMol) format, which can be processed by many molecular visualization programs.

off

The construction of the MOLDEN input and the COORD.xyz file is turned off.

Default:

molden=on

Example:

if you do not need Molden input and the COORD.xyz file, add molden=off

mulmet

Specifies the multipole approximation used for the evaluation of pair energies of distant pairs.

Options:
0

The simplified dipole-dipole estimate of Ref. 134 will be used.

1

Full dipole-dipole estimate [128].

2

All the terms are included in the multipole expansion up to the contributions of the quadrupole moment [105].

3

All the terms are included in the multipole expansion up to the contributions of the octapole moment [105].

Default:

mulmet=3 if localcc=2016 or 2018,
mulmet=0 otherwise

Example:

to use the octapole approximation set mulmet=3

mult

Spin multiplicity ($2S+1$) of the Hartree–Fock or Kohn–Sham wave function. If a CI or CC calculation is also performed, the same multiplicity is supposed for the ground-state wave function. For excited states the multiplicity will be arbitrary, only $M_{S}$ is conserved. For closed-shell reference determinants the multiplicity (strictly speaking the parity of $S$) can be controlled by keywords nsing and ntrip, see below.

Options:

$<$any positive integer$>$

Default:

for atoms the corresponding experimental multiplicity is set, for molecules mult=1 (singlet) for an even number of electrons, mult=2 (doublet) otherwise

Example:

for a triplet state one should give mult=3

nacto

Number of active occupied spinorbitals. By default, nacto pieces of spinorbitals under the Fermi level are supposed to be active. This can be overwritten using keyword active, which enables the user to select the active orbitals manually (see the description of keyword active). In a MRCI/CC calculation a complete active space (CAS) is supposed defined by keywords nacto and nactv (or alternatively by active) and up to $n$-fold excitations from the reference determinants of this space are included in the excitation manifold, where $n$ is determined by keyword calc (2 for CCSD, 3 for CCSDT, …). See Ref. 83 for more details. See also keywords nactv, maxex, and maxact.

Options:

$<$any positive integer$>$ or 0

Default:

nacto=0

Example:

for two active occupied spin-orbitals give nacto=2

nactv

Number of active virtual spinorbitals. By default, nactv pieces of spinorbitals above the Fermi level are supposed to be active, which can be overwritten using keyword active. For a detailed description see keyword nacto.

Options:

$<$any positive integer$>$ or 0

Default:

nactv=0

Example:

for two active virtual spin-orbitals give nactv=2

nafalg

Specifies how natural auxiliary functions (NAFs) will be constructed in the case spin-unrestricted MOs. NAFs can be calculated by diagonalizing ($\mathbf{W}^{\alpha}$ $+$ $\mathbf{W}^{\beta}$)/2 or $\mathbf{W}^{\alpha}$ (see Ref. 84 for the definitions). The latter option is somewhat more efficient but can be dangerous for processes involving atoms.

Options:
albe

NAFs are constructed from ($\mathbf{W}^{\alpha}$ $+$ $\mathbf{W}^{\beta}$)/2.

alpha

NAFs are constructed from $\mathbf{W}^{\alpha}$.

Default:

nafalg=albe

Example:

to use $\mathbf{W}^{\alpha}$ set nafalg=alpha

naf_cor

Specifies whether natural auxiliary functions (NAFs) will be used for density-fitting correlated calculations and also specifies the threshold for the occupation numbers of NAFs (see Ref. 84).

Options:
off

NAFs will not be constructed.

$<$any real number in the [0,1] interval$>$

A NAF basis will be constructed and NAFs with occupation numbers smaller than this number will be dropped.

on

Equivalent to naf_cor=1e-2 for local correlation methods, and to naf_cor=5e-3 otherwise.

Default:

according to the value of lcorthr for local correlation methods, naf_cor=0.1 for the reduced-cost excited-state approaches of Refs. 96 and (see keyword redcost_exc), naf_cor=off otherwise.

Example:

to use NAFs and set a threshold of $10^{-2}$ type naf_cor=1e-2

naf_scf

Specifies whether NAFs will be used for density-fitting SCF calculations and also specifies the threshold for the occupation numbers of NAFs (see Ref. 84). The syntax is analogous with that for keyword naf_cor.

naftyp

Specifies how NAFs will be constructed in the case of local correlation calculations. NAFs are constructed by diagonalizing the $\mathbf{W}=\mathbf{J}^{\mathrm{T}}\mathbf{J}$ matrix where $\mathbf{J}$ is a particular block of the three-center Coulomb integral matrix (see Refs. 84 and 106 for details).

Options:
Jai

NAFs are constructed from $J^{P}_{ai}$

Jpi

NAFs are constructed from $J^{P}_{pi}$

Jpq

NAFs are constructed from $J^{P}_{pq}$

Jmi

NAFs are constructed from $J^{P}_{\mu i}$

Default:

naftyp=Jpq for local CC and localcc=2016 or 2018 as well as for the reduced-cost excited-state approaches of Refs. 96 and (see keyword redcost_exc), naftyp=Jai for local MP2 and localcc=2016 or 2018

Example:

to construct NAFs using $J^{P}_{pq}$ set naftyp=Jpq

nchol

Number of Cholesky vectors/quadrature points for the Laplace integral in the case methods based on the decomposition of energy denominators. See also the description of keyword dendec.

Options:
auto

The number of Cholesky vectors/quadrature points will be automatically determined to achieve the required precision.

$<$any positive integer$>$

The number of Cholesky vectors/quadrature points will also be automatically determined but the maximum number of the vectors cannot exceed this number.

Default:

nchol=auto

Example:

to use ten Cholesky vectors/quadrature points give nchol=10

ndeps

Step size for the numerical differentiation in atomic units.

Options:

$<$any positive real number$>$

Default:

ndeps=1e-3

Example:

for a step size of $5\cdot 10^{-4}$ a.u. for numerical Hessian evaluation set ndeps=5e-4

nstate

Number of electronic states including the ground state and excited states. In non-relativistic calculations, for closed-shell reference determinants nstate is supposed to be the number of singlet states. See also keywords nsing and ntrip.

Options:

$<$any positive integer$>$

Default:

nstate=max(1, nsing+ntrip)

Example:

for three states give nstate=3

nsing

Number of singlet electronic states (strictly speaking the number of of states with $M_{S}=0$ and $S$ is even) including the ground state and excited states. Use this option only for non-relativistic calculations and closed-shell reference determinants, it should be zero otherwise. In the case of closed-shell reference determinants a partial spin-adaptation is possible, see Ref. 82. This enables us to search for singlet and triplet roots separately. See also keywords nstate and ntrip.

Options:

$<$any positive integer$>$

Default:

nsing=1 for closed-shell reference determinants, nsing=0 otherwise

Example:

for two singlet states give nsing=2

ntrip

Number of triplet electronic states (strictly speaking the number of of states with $M_{S}=0$ and $S$ is odd) including the ground state and excited states. Use this option only for non-relativistic calculations and closed-shell reference determinants, it should be zero otherwise. See the description of keywords nstate and ntrip.

Options:

$<$any positive integer$>$

Default:

ntrip=0

Example:

for two triplet states give ntrip=2

occ

Specifies the occupation of the Hartree–Fock determinant.

Options:
1.

If this keyword is not given, the occupation is automatically determined in the SCF calculations.

2.

For RHF calculations the occupation should be given in the following format:
occ=$,,\dots,$
where $$ is the number of occupied orbitals in irrep $i$, and $N_{ir}$ is the number of irreps.

3.

For ROHF and UHF calculations the occupation should be given as
occ=$,\dots,$/$,\dots,$
where $$ is the number of occupied $\sigma$ spinorbitals in irrep $i$.

Default:

occ is not specified, that is, the occupation is set by the SCF program.

Examples:
1.

Water, RHF calculation:
occ=3,1,1,0

2.

Water, UHF calculation:
occ=3,1,1,0/3,1,1,0

3.

Carbon atom, ROHF or UHF calculation:
occ=2,0,0,0,0,1,0,1/2,0,0,0,0,0,0,0

optalg

Specifies the optimization algorithm used for geometry and basis set optimizations. For basis set optimization, at the moment, the downhill simplex method of Nelder and Mead [108] is the only available option. A geometry optimization can be carried out by either the Broyden–Fletcher–Goldfarb–Shanno (BFGS) or the simplex algorithm.

Options:
simplex

the simplex method of Nelder and Mead.

BFGS

the BFGS algorithm.

Default:

optalg=simplex for basis set optimization, optalg=BFGS otherwise.

Example:

To run a geometry optimization with the simplex algorithm set optalg=simplex

optmaxit

Maximum number of iteration steps allowed in a geometry or basis set optimization. If the simplex algorithm is used, i.e., optalg=simplex, the maximum number of function evaluations is also controlled by the parameter optmaxit: it is set to 15$\times$optmaxit. If the optimization is terminated with a message “the maximum number of function evaluation is exceeded”, then you can increase the value of optmaxit appropriately.

Options:

$<$any positive integer$>$

Default:

optmaxit=50

Example:

to allow 60 iteration steps set optmaxit=60

optetol

Convergence threshold for geometry or basis set optimization. If the simplex algorithm is used, i.e., optalg=simplex, the optimization is terminated when the energy difference (in E${}_{h}$) becomes less then this value and the optstol criterion is also fulfilled. In addition to the criterion for the gradient (optgtol) and the step-size (optstol) the energy change between the cycles is also monitored. For a successful geometry optimization it is required that the optgtol criterion is satisfied and either the energy difference between the last two steps becomes less than this value (in E${}_{h}$) or the optstol criterion is met.

Options:

$<$any positive real number$>$

Default:

optetol=1e-6

Example:

for a convergence threshold of $5\cdot 10^{-6}$ E${}_{h}$ set optetol=5e-6

optgtol

Convergence threshold for geometry optimization, upper limit (in E${}_{h}$/bohr) for the maximum gradient component. For a successful geometry optimization this criterion must be fulfilled.

Options:

$<$any positive real number$>$

Default:

optgtol=1e-4

Example:

for a convergence threshold of $3\cdot 10^{-4}$ E${}_{h}$/bohr set optgtol=3e-4

optstol

Convergence threshold for geometry or basis set optimization. For the latter the optimization is terminated when the maximum change in the parameters becomes less then this value and the optetol criterion is also fulfilled, for the former this criterion is met when the maximum step-size from the previous step (in bohr) is lower than this value. The geometry optimization is terminated successfully if, in addition to the optgtol criterion obeyed, either this criterion is met or the energy difference between the last two steps becomes less than optetol.

Options:

$<$any positive real number$>$

Default:

optstol=1e-3 for a basis set optimization, optstol=1e-4 otherwise.

Example:

to set a threshold of $10^{-5}$ bohr for a geometry optimization type optstol=1e-5

orblocc

Specifies what type of orbital localization is performed for the core molecular orbitals.

Options:

All the options introduced for keyword orbloco also work for orblocc, see the description of keyword orbloco for details.

Default:

orblocc=orbloco if localcc=2013, or core=corr and
localcc$\neq$off; orblocc=off otherwise

Example:

to localize of core orbitals with the Pipek–Mezey algorithm specify orblocc=on

orblocguess

Initial guess for the orbital localization.

Options:
cholesky

Guess orbitals are calculated by the Cholesky decomposition of the one-particle density matrix [5].

restart

Guess orbitals for the localization are read from the
MOCOEF.LOC file produced in a previous run where orbital localization was performed.

Orbitals are read from the MOCOEF.LOC file and directly employed at later steps of the calculation without any change. The locality of the orbitals is not checked and the orbital localization is skipped entirely.

Note:

The combination of orblocguess=restart (or orblocguess=read) with scfiguess=off can be particularly useful if the result files of the converged SCF iteration and orbital localization steps are available and only the local correlation step should be repeated with different settings. See also the description of scfiguess=off for this option.

Default:

orblocguess=cholesky

Example:

to speed up the orbital localization by using the localized orbitals of a previous calculation as guess specify orblocguess=restart.

orbloco

Specifies what type of orbital localization is performed for occupied molecular orbitals.

Options:
off

No orbital localization.

boys

Boys localization is performed [36].

pm

Pipek–Mezey localization is performed [129].

IBO

intrinsic bond orbitals of Knizia are constructed [68].

cholesky

localized orbitals are calculated by the Cholesky decomposition of the one-particle density matrix [5].

Default:

orbloco=off in the general case, orbloco=boys for local correlation calculations

Example:

to carry out Pipek–Mezey localization for the occupied orbitals type orbloco=pm

orblocv

Specifies what type of orbital localization is performed for virtual molecular orbitals.

Options:

All the options introduced for keyword orbloco excepting IBO also work for orblocv, see the description of keyword orbloco for details. In addition, for local correlation calculations there is one more option:

pao

Projected atomic orbitals.

Default:

orblocv=off in the general case, orblocv=pao for local correlation calculations

Example:

to carry out Boys localization for the virtual orbitals type orblocv=boys

osveps

Threshold for the occupation numbers of orbital specific virtual orbitals (OSVs) used at the evaluation of pair correlation energies in local MP2 and dRPA calculations. See the description of keyword wpairtol for more details.

Options:
off

OSVs are not constructed and not dropped

$<$any real number in the [0,1] interval$>$

Orbitals with occupation numbers smaller than this number will be dropped.

Default:

osveps=1e-3 for localcc=2015, osveps=off for localcc=2016 or 2018

Example:

to set a threshold of $10^{-4}$ type osveps=1e-4

ovirt

This keyword controls the cost reduction approaches based on natural orbital (NO) or optimized virtual orbitals (OVOs) techniques. For a ground-state correlation calculation, if this keyword is set, the virtual MOs will be transformed to MP2 NOs or OVOs [109]. Subsequently the virtual space will be truncated on the basis of the populations of the orbitals, which can be controlled by keywords eps and ovosnorb. See Ref. 136 for more details.

Options:
off

The virtual MOs are not changed.

MP2

MP2 NOs will be used.

OVOS

Optimized virtual orbitals will be used.

Default:

ovirt=off

Example:

to use MP2 NOs for a ground-state CC calculation give ovirt=MP2

ovltol

Tolerance for the eigenvalues of the AO overlap matrix. Eigenvectors corresponding to the eigenvalues lower that ovltol will be removed to cure the linear dependence of the AO basis set.

Options:
$<$any positive real number or zero$>$

This number will be used as the threshold for the eigenvalues of the overlap matrix.

Default:

ovltol=1e-7

Example:

to keep all the basis functions set ovltol=0.0

ovosnorb

Specifies the retained percentage of virtual orbitals in an optimized virtual orbitals (OVOs) calculations. ovosnorb % of virtual orbitals will be retained.

Options:

$<$any number between 0 and 100$>$

Default:

ovosnorb=80.0

Example:

to retain only 70 % of the virtuals give ovosnorb=70.0

popul

This keyword controls the wave function analysis.

Options:
off

No wave function analysis is performed.

Mulli

A population analysis is also performed, and Mulliken and Löwdin atomic charges as well as Mayer bond orders are computed [101, 93].

IAO

In addition to the above parameters, intrinsic atomic orbitals (IAOs) are constructed and IAO partial charges are calculated [68].

Default:

popul=Mulli if dens $\neq$ 0, popul=off otherwise

Example:

to calculate IAO charges set popul=IAO

pressure

The pressure in Pa at which the thermodynamic properties are evaluated (see also keyword freq).

Options:

$<$any positive integer$>$

Default:

pressure=100000

Example:

for 1 atm set pressure=101325

ptfreq

Frequency of the perturbation for frequency-dependent properties in atomic units (available only with Cfour). See Refs. 78 and 113 for more details.

Options:

$<$any real number$>$

Default:

ptfreq=0.0

Example:

to set a frequency of 0.1 a.u. give ptfreq=0.1

qmmm

This keyword tells Mrcc that a QM/MM calculation is performed and the point charges included in the input file must be processed. This keyword is automatically added to the MINP file by the MM program conducting the QM/MM calculations, and you do need to bother with it. Use this keyword in the only case if you want to add point charges manually.

Options:
off

QM/MM calculation is not performed and no point charges are added.

Amber

Currently this is the only option, the Amber MD code will be used for the QM/MM calculation or point charges will be added.

Note:

Point charges can be manually added to the end of the input file in the following format:
pointcharges
$<$number of point charges$>$
$$ $$ $$ $$
$$ $$ $$ $$
⋮
where $x_{i}$, $y_{i}$, and $z_{i}$ are the Cartesian coordinates and $q_{i}$ is the charge for point charge $i$. The charge must be given in a.u., while for the coordinates the same unit must be used as for the specification of the molecular geometry.

Default:

qmmm=off

Example:

to add two point charges with coordinates (0, 0, 1) and (0, 1, 0) a.u. (provided that the geometry is also given in bohr) and charges of 0.5 and -0.5 a.u. qmmm=Amber should be set anywhere in the MINP file, and the following lines should be added to the end the file:
pointcharges
2
0.0 0.0 1.0 0.5
0.0 1.0 0.0 -0.5

qscf

Use this option to carry out quadratic SCF calculations. One can use either Newton or trust region iterations with optional line search.

Options:
off

No quadratic SCF, conventional SCF algorithm will be executed.

Newton

Simple Newton iteration.

NewtonL

Simple Newton iteration with line search.

AugHess

Trust region method with augmented Hessian algorithm without line search. The trust radius is updated according to the scheme described in Ref. 55.

AugHessM

Trust region method with augmented Hessian algorithm using line search. The trust radius update scheme is similar to the above one but uses different coefficients.

AugHessL

Trust region method with augmented Hessian algorithm using line search. The trust radius is updated as described in Ref. 55, but this method takes into account the step length found by the line search.

AugHessG

Trust region method with augmented Hessian algorithm using line search. The trust radius update scheme takes into account the change of the gradient in the consecutive iterations.

Note:

The simple Newton iteration schemes are only efficient in the vicinity of a minimum and not recommended in the general case. We recommend the use of the AugHessG or AugHessL options in difficult cases.

Default:

qscf=off

Example:

UKS calculation with the B3LYP functional using the AugHessG algorithm:
calc=B3LYP
scftype=UHF
qscf=AugHessG

redcost_exc

This keyword controls the cost reduction approaches based on natural orbital (NO) and natural auxiliary function (NAF) techniques for excited states. For a CIS(D${}_{\infty}$), ADC(2), and CC2 calculation, if redcost_exc$\neq$off, the reduced-cost approach of Refs. 96 and is invoked, and truncated state-averaged MP2/CIS(D) NOs as well as NAFs will be used. The lnoepso, lnoepsv, naf_cor, and naftyp keywords will be automatically set depending on the selected option for redcost_exc. For CIS, TD-HF, TDA, TD-DFT calculations the NAF approximation will be invoked if redcost_exc=8.

Options:
off

The cost reductions techniques will not be used.

on

Reduced-cost calculation will be executed with default settings (see also Table 4).

cust

Reduced-cost calculation will be executed with customized truncation thresholds. The threshold for the complete MO space NAFs (CS-NAFs), CIS coefficients, orbital energies, and linear dependency must be specified in the subsequent line, respectively, as
$<$threshold 1$>$ $<$threshold 2$>$ $<$threshold 3$>$ $<$threshold 4$>$
The default values of these thresholds are 0.1 a.u., 0.35, 0.15 a.u., $10^{-7}$, respectively. These are used with if any other option is chosen. See Ref. for the detailed description of these parameters.

$<$positive integers from 1 to 10$>$

See Table 4.

Default:

redcost_exc=off

Examples:
1.

Reduced-cost ADC(2) calculation for the lowest singlet excited state with the default settings proposed in Ref. :
nsing=2
redcost_exc=on

2.

Reduced-cost TD-DFT calculation with the PBE0 functional for the lowest 3 singlet excited states of a molecule:
calc=PBE0
redcost_exc=8
nstate=4

Table 4: Options for keyword redcost_exc. CS-NAF – complete MO space NAFs are used, NO – frozen natural orbitals are used, Can. – the NO space is augmented with canonical virtual orbitals, RS-NAF – restricted NO space NAFs are used. See Ref. for more details.
 Option CS-NAF NO Can. RS-NAF Note off no no no no Default on yes/no yes yes yes Equivalent to 1 or 6 depending on keyword localcc cust yes yes yes yes Customized thresholds 1 yes yes yes yes The approach presented in Ref. , default if localcc=off and redcost_exc=on 2 no yes no yes The approach presented in Ref. 96 3 no no no yes The approach presented in Ref. 96 without NOs 4 no yes no no The approach presented in Ref. 96 without NAFs 5 no yes yes no 6 no yes yes yes Default if localcc=2016 and redcost_exc=on 7 yes yes no yes 8 yes no no no For reduced-cost CIS, TD-HF, TDA, TD-DFT, … 9 yes yes yes no 10 yes yes no no
refdet

The reference determinant (Fermi-vacuum) for CI/CC calculations can be specified using this keyword. By default the reference determinant is identical to the HF determinant, but sometimes it is necessary to change this.

Options:
none

The reference determinant is identical to the HF determinant.

serialno

Using this option one can define the occupation of the correlated orbitals in the reference determinant specifying their serial numbers. This option requires three more lines. In the first line the serial numbers of the doubly-occupied orbitals must be given, while in the second and third lines those orbitals should be specified which are singly-occupied by an alpha or a beta electron, respectively. For the format of these lines see the description of the serialno option of the active keyword. For relativistic calculations the occupation of the spinors (i.e., not that of the Kramers-pairs) should be given. For technical reasons all electrons are treated as alpha electrons and the serial numbers of the occupied spinors must be given in the second line, the first and third lines must be left blank.

vector

Using this option one can set the occupation numbers for each correlated orbital. In the subsequent line an integer vector should be supplied with as many elements as the number of correlated orbitals (correlated spinors for relativistic calculations, not Kramers-pairs!). The integers must be separated by spaces. In the case of non-relativistic calculations type 2 for doubly-occupied orbitals, 1 for open-shell orbitals with alpha electron, -1 for open-shell orbitals with beta electron, 0 otherwise. In the case of relativistic calculations type 1 for each occupied spinor, 0 otherwise.

Notes:
1.

Frozen orbitals must not be considered here in any case.

2.

If the MO integrals are taken over from another program, the numbering of orbitals may be different from that of the parent program. Here the order of MOs: doubly occupied, open shell, virtual; and in each of this blocks the MOs are reordered according to the orbital energies (natural orbital occupations in the case of MCSCF orbitals).

3.

If the MO integrals are taken over from another program, and this line is omitted, the program will fill the orbitals with electrons from the bottom automatically. In this manner we do not need this line for closed shells or a doublet ref. det., but e.g. for high spin states the Fermi vacuum must be defined here.

4.

For relativistic calculations (Dirac interface) this line is always required. The spinors are symmetry-blocked according to the Fermion irreps of the corresponding double group. Complex conjugate irreps follow each other. Within each irrep the spinors are numbered according to orbital energies. Please note that this line is automatically printed by the dirac_mointegral_export program, and you do not have to do it by hand. However, for technical reasons, always a closed-shell occupation is generated, and you may need to remove or add some electrons.

Default:

refdet=none, that is, the reference determinant is identical to the HF determinant.

Examples:
1.

We have 20 correlated orbitals, 10 electrons, and we are interested in a high-spin triplet state. Suppose that orbitals 1 to 4 are doubly-occupied while orbitals 5 and 6 are singly occupied by alpha electrons. Using the serialno option the input should include the following four lines (note the blank line at the end):
refdet=serialno
1-4
5,6

2.

The same using the vector option:
refdet=vector
2 2 2 2 1 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0

3.

We perform a relativistic calculation for the Be atom with 20 correlated spinors. We have 6, 6, 4, and 4 spinors in the four Fermion irreps, E${}_{1/2g}$, E${}_{-1/2g}$, E${}_{1/2u}$, and E${}_{-1/2u}$ of the $C_{2h}^{*}$ double group, respectively, and two occupied spinors in both of the gerade irreps. Thus using the vector option the occupation vector should be given as:
refdet=vector
1 1 0 0 0 0 1 1 0 0 0 0 0 0 0 0 0 0 0 0

4.

The same using the serialno option (note the blank lines):
refdet=serialno

1,2,7,8

rest

Use this keyword to restart canonical (i.e., not local) CI and CC calculations from previously calculated wave function parameters (cluster amplitudes, CI coefficients, $\lambda$ amplitudes, etc.) if ccprog=mrcc. For restarting local correlation calculations see keyword lccrest.

Options:
1

The program restarts from the previously calculated parameters.

2

The program executes automatically the lower-level calculations of the same type consecutively (e.g., CCSD, CCSDT, and CCSDTQ if CCSDTQ is requested) and restarts each calculation from the previous one (this is only available for energy calculations).

3

Same as rest=1, however, only selected roots from the previous calculation will be used as initial guess. The serial number of the roots must be specified in the subsequent line as
$\dots$
where $$ is the serial number of the states. The number of states given here must be equal to nstate or nsing + ntrip. Please note that the ground state solution is not automatically selected, it should also be given here if needed. It is recommended to use root following (diag=follow) together with this option.

4

Same as rest=2 but the initial vectors are selected as in the case of rest=3.

Notes:

use the restart option, e.g.,

1.

after system crash.

2.

if the iteration procedure did not converge in the given number of steps.

3.

for geometry optimization.

4.

for potential curve calculations.

5.

if you are interested in high-order CC/CI energies. Then it is worth restarting the calculations with higher excitations using the converged vectors of the same method including lower excitations, e.g., CCSDT using the converged CCSD amplitudes, CCSDTQ using the CCSDT amplitudes, and so on. With this trick 1-3 iteration steps can be saved usually, but more ones in the case of strong static correlation (i.e., large cluster amplitudes). Use exclusively rest=2 for this purpose (that is, not rest=1)!

6.

if you are interested in calculating the energies for all methods in a hierarchy (e.g., executing all CC methods up to CCSDTQP). Use exclusively rest=2 for this purpose (that is, not rest=1)!

7.

to generate brute-force initial guess for excited-state calculations (rest=3 or 4). That is, if you do not want to bother with the initial guess for excited states, but you know approximately the energy of the excited states, then execute a low-level method for many roots. Use LR-CCS (calc=CCS) and CIS (calc=CIS), respectively, for higher-order LR-CC and CI calculations. Select the desired roots on the basis of their energies, and use them as initial guess in high-level calculations. (For other options for the initial guess for excited-state calculations see keyword ciguess.)

8.

Please note that the program always needs the file fort.16 from the previous calculation for the restart and also fort.17, if more than one root is sought or for geometry optimization.

Default:

rest=0

Examples:
1.

to restart a CC calculation after power failure set rest=1

2.

to restart a LR-CCSD calculation using the first, third, and fifth roots of a previous LR-CCS calculation the input should include the following two lines:
rest=3
1 3 5

rgrid

Specifies the radial integration grid for DFT calculations. For the details of the grid construction see the description of keyword agrid. See also the description of keyword grtol.

Options:
Log3

the Log3 quadrature of Mura and Knowles [102]

GS

Gauss–Chebyshev quadrature. A modified version of the mapping function of Ref. 72 is employed: the function is scaled by the atomic scaling parameters of Becke [7].

EM

Note:

the number of grid points is calculated by the max{20,5*(3*grtol/2+$i$-8)} formula with $i$ as the number of the row in the periodic table where the atom is located [72]. To change the number of radial integration points set the value of grtol.

Default:

rgrid=Log3

Examples:
1.

to use the Euler–Maclaurin scheme set rgrid=EM

2.

to use the default Log3 grid with more points set grtol=12

rohftype

Specifies the type of the ROHF orbitals. See also the description of keyword scftype.

Options:
standard

Standard ROHF orbitals obtained by diagonalizing the ROHF Fock-matrix.

semicanonical

Semicanonical ROHF orbitals obtained by separately diagonalizing the alpha and beta UHF Fock-matrices constructed using the converged ROHF orbitals.

Notes:
1.

rohftype=semicanonical is required for perturbative CC methods if ROHF orbitals are used, otherwise the expressions for the perturbative corrections are not correct. Iterative CC and CI methods are invariant to the choice of ROHF orbitals (if all electrons are correlated).

2.

It is very important to give this keyword if Mrcc is used together with another code and ROHF orbitals are used since this keyword tells Mrcc what type of ROHF orbitals are taken over from the other code.

Default:

rohftype=standard for iterative CC and CI methods,
rohftype=semicanonical for perturbative methods.

Example:

to use semicanonical ROHF orbitals for iterative CC methods give rohftype=semicanonical

scfalg

Specifies what type of SCF algorithm is to be used.

Options:
disk

Conventional SCF algorithm, two-electron integrals are stored on disk.

direct

Direct SCF algorithm, two-electron integrals are recalculated in each iteration step.

auto

Based on the size and geometry of the molecule the program will automatically select the more efficient one from the above options.

Default:

scfalg=auto

Example:

to run direct SCF add scfalg=direct

scfdamp

Specifies whether damping of the SCF density matrices is performed.

Options:
off

No damping.

$<$any real number in the [0,1] interval$>$

In each SCF iteration cycle the new and old SCF density matrices are mixed by factors (1-scfdamp) and scfdamp, respectively.

on

Equivalent to scfdamp=0.7

Default:

scfdamp=off

Example:

to use a damping factor of $0.8$ type scfdamp=0.8

scfdiis

Specifies if DIIS convergence acceleration is used in the SCF calculations.

Options:

on or off

Default:

scfdiis=on

Example:

to turn off DIIS convergence accelerator add scfdiis=off

scfdiis_end

Specifies the last iteration step in which the DIIS convergence acceleration is applied.

Options:

$<$any positive integer$>$

Default:

scfdiis_end=scfmaxit, that is, the DIIS procedure is not turned off.

Example:

to turn off the DIIS convergence accelerator after iteration step 20 give scfdiis_end=20

scfdiis_start

Specifies the first iteration step in which the DIIS convergence acceleration is applied.

Options:

$<$any positive integer$>$

Default:

scfdiis_start=1, that is, the DIIS procedure is active from the first iteration.

Example:

to turn on the DIIS convergence accelerator in iteration step 5 give scfdiis_start=5

scfdiis_step

Specifies the frequency of DIIS extrapolations. The extrapolation will be carried out in every scfdiis_step’th iteration cycle.

Options:

$<$any positive integer$>$

Default:

scfdiis_step=1, that is, the DIIS extrapolation is performed in each iteration step.

Example:

to carry out DIIS extrapolation only in every second iteration step give scfdiis_step=2

scfdtol

Convergence threshold for the density matrix in SCF calculations. The RMS change in the density matrix will be smaller than $10^{-{\tt scfdtol}}$.

Options:

$<$any integer$>$

Default:

scfdtol=scftol+2 for frequency calculations, otherwise
scfdtol=scftol+1 for correlation calculations, scfdtol=scftol for SCF calculations

Example:

for an accuracy of $10^{-8}$ one must give scfdtol=8

scfext

Specifies the number of Fock-matrices used for the DIIS extrapolation in SCF calculations.

Options:

$<$any positive integer$>$

Default:

scfext=10

Example:

to increase the number of DIIS vectors to 15 give scfext=15

scfiguess

Initial guess for the SCF calculation.

Options:

Superpositions of atomic densities. For each atom a density-fitting UHF calculation is performed, and the initial one-particle density matrix is constructed from the averaged alpha and beta atomic densities.

ao

Atomic density initial guess. The initial one-particle density matrix is constructed from diagonal atomic densities derived from the occupation of the atoms. It is efficient for Dunning’s basis sets.

core

Core Hamiltonian initial guess. The initial MOs are obtained by diagonalizing the one-electron integral matrix.

mo

The SCF calculation will use the MO coefficients obtained in a previous calculation and stored in the MOCOEF file. The calculation can only be restarted from the MOs computed with the same basis set.

restart

The SCF calculation will use the density matrices obtained in a previous calculation and stored in the SCFDENSITIES file. If the calculation is restarted from the densities obtained with another basis set, the VARS file is also required.

off

No SCF calculation will be performed, but the Fock-matrix and the MO coefficients obtained in a previous calculation will be used in the correlation calculations. This requires the FOCK, MOCOEF, and VARS files from the previous calculation.

min

A density fitting SCF calculation will be performed using the cc-pVTZ-min minimal basis set (see the description of keyword basis), and the resulting density will be used as initial guess. In the minimal-basis SCF calculation the AO basis set is used as the auxiliary basis, and loose convergence thresholds are employed, consequently, the energy is unreliable and should not be used for any purpose.

small

A density fitting SCF calculation will be performed using a smaller basis set which must be specified by keyword basis_sm (see the description of keyword basis_sm), and the resulting density will be used as initial guess.

Notes:
1.

Restarting from densities obtained with a bigger basis set is not allowed.

2.

To restart SCF runs from the results of DFT embedding calculations with the Huzinaga-equation- or projector-based approaches use MO coefficients, i.e., scfiguess=mo, since only the subsystem densities are stored at the end of the embedding calculation but the MOs are available for the entire system.

Default:

Examples:
1.

For a core Hamiltonian initial guess set scfiguess=core

2.

For restarting the SCF calculation from the results of a calculation performed with the same basis set type scfiguess=restart. Note that you need the SCFDENSITIES file from the previous run.

3.

You would like to generate a good initial guess for an aug-cc-pVTZ SCF calculation. First, run a calculation with the cc-pVTZ basis set (cc-pVTZ-min is also a good option), that is, your input file should contain the
basis=cc-pVTZ
line. Then, restart your aug-cc-pVTZ calculation from the cc-pVTZ density matrix. To that end the MINP file should include the following lines:
basis=aug-cc-pVTZ
scfiguess=restart

Note that the SCFDENSITIES and the VARS files from the cc-pVTZ run must be copied to the directory where the aug-cc-pVTZ calculation is executed.

4.

The calculations in the previous example can be run more simply, in one step using the small option and the basis_sm keyword as
basis=aug-cc-pVTZ
basis_sm=cc-pVTZ
scfiguess=small

scflshift

Level shift parameter for the SCF calculation.

Options:
off

No level shifting.

$<$any real positive number $>$

The value of the level shift parameter in a.u.

on

Equivalent to scflshift=0.2

Default:

scflshift=off

Example:

To use a level shift value of 0.5 a.u. give scflshift=0.5

scfmaxit

Maximum number of iteration steps in SCF calculations.

Options:

$<$any positive integer$>$

Default:

scfmaxit=50

Example:

to increase the maximum number of SCF iterations to 200 give scfmaxit=200

scftol

Convergence threshold for the energy in SCF calculations. The energy will be accurate to $10^{-{\tt scftol}}$ E${}_{h}$.

Options:

$<$any integer$>$

Default:

scftol=max(8,cctol) for property calculations,
scftol=max(6,cctol) otherwise

Example:

for an accuracy of $10^{-8}$ E${}_{h}$ one must give scftol=8

scftype

Specifies the type of the Hartree–Fock/Kohn–Sham SCF procedure, or the type of the molecular orbitals if the MO integrals are computed by other programs. See also the description of keyword rohftype.

Options:

RHF, ROHF, UHF, or MCSCF

Notes:
1.

scftype=MCSCF is only available if Mrcc is used together with Columbus or Molpro. In that case the MCSCF calculation is performed by the aforementioned codes and the transformed MO integrals are passed over to Mrcc.

2.

It is very important to give this keyword if Mrcc is used together with another code and ROHF or MCSCF orbitals are used since this keyword tells Mrcc that the orbitals are not canonical HF orbitals. Please also set keyword rohftype in this case.

3.

If a HF-SCF calculation is run, the type of the SCF wave function can also be controlled by keyword calc. See the description of calc.

4.

For DFT calculations only the RHF and UHF options can be used, which, in that case, instruct the code to run RKS or UKS calculations, respectively.

Default:

scftype=RHF for closed-shell systems, scftype=UHF for open shells.

Example:

to use ROHF for open-shell systems type scftype=ROHF

scsps

Scaling factor for the antiparallel-spin component of the correlation energy in spin-component scaled MP2 (SCS-MP2) calculations [46].

Options:
$<$any real number$>$

the antiparallel-spin component of the correlation energy will be scaled by this number.

Default:

scsps=6/5

Example:

to set a scaling factor of 1.5 type scsps=1.5

scspt

Scaling factor for the parallel-spin component of the correlation energy in spin-component scaled MP2 (SCS-MP2) calculations [46].

Options:
$<$any real number$>$

the parallel-spin component of the correlation energy will be scaled by this number.

Default:

scspt=1/3

Example:

to set a scaling factor of 0.5 type scspt=0.5

spairtol

Threshold for the selection of strong pairs in local MP2, dRPA, and CC methods. For each orbital pair an estimate of the pair correlation energy is calculated (see the description of keyword wpairtol). An orbital pair will be considered as strong pair if the absolute value of the pair correlation energy estimate is greater than spairtol. In the subsequent calculations strong pairs will be treated at a higher level, while for the other pairs (weak and distant) the corresponding pair correlation energy estimates will be added to the correlation energy. See also Refs. 85 and 106 for more details.

Options:
off

the local MP2 or dRPA pair correlation energy estimate is not calculated, an orbital pair will be considered as strong pair in this case if the absolute value of the available pair correlation energy estimate is greater than wpairtol. See also the description of keyword wpairtol and Ref. 105.

$<$any positive real number$>$

Orbital pairs with pair correlation energy estimates greater than this number (in E${}_{h}$) will be considered as strong pairs.

Default:

spairtol=1e-4 for localcc=2015, spairtol=off for
localcc=2016 and localcc=2018

Example:

to set a threshold of $10^{-5}$ E${}_{h}$ type spairtol=1e-5

symm

Spatial symmetry (irreducible representation) of the state. See Sect. 13 for the implemented point groups, conventions for irreps, etc.

Options:
0

No symmetry adaptation, that is, all calculations will use the $C_{1}$ point group

off

Equivalent to symm=0

1, 2, …, 8

Serial number of the irrep (see Sect. 13).

$<$irrep label$>$

Label for the irrep (see Sect. 13).

Note:

Irreps can only be specified by their serial numbers if Mrcc is used with another program. In that case please check the manual or output of the other program system for the numbering of irreps.

Default:

by default the state symmetry is determined on the basis of the occupation of the HF determinant.

Examples:
1.

for the second irrep of the point group type symm=2

2.

for the B${}_{1u}$ irrep of the $D_{2h}$ point group type symm=B1u

talg

Specifies the algorithm for the calculation of the (T) correction in the case of the CCSD(T) method.

Options:
occ

The outmost loops run over the occupied indices of the triples amplitudes.

virt

The outmost loops run over the virtual indices of the triples amplitudes.

lapl

Laplace transformed (T) energy expression for localcc=2016 or 2018 according to Ref. 104. See also the laptol keyword to set the accuracy of the numerical Laplace transform.

topr

The T0${}^{\prime}$ semi-canonical approximation of the local (T) expression according to Ref. 104.

Default:

talg=occ for conventional CCSD(T) calculations, talg=lapl for the local CCSD(T) scheme of localcc=2016 or 2018, and talg=virt for local CCSD(T) with localcc=2015 or 2013.

Notes:
1.

For algorithmic reasons in the case of previous local CCSD(T) schemes (localcc=2015 or 2013) talg=virt is the only option. For localcc=2016 or 2018, the default is talg=lapl, and talg=virt is used if lcorthr=0 is set.

2.

For conventional CCSD(T) calculations talg=occ is recommended since the algorithm is somewhat faster than the other one. In turn, its memory requirement is higher. The program checks automatically if the available memory is sufficient for the first algorithm (i.e., talg=occ). If this is not the case, talg will be automatically set to virt.

3.

In the case of localcc=2016 or 2018 the talg=topr algorithm is approximately three times faster than the talg=lapl with its default settings, but considerably less accurate. For quick exploratory calculations the talg=lapl algorithm is recommended in combination with laptol=0.1.

Example:

to change the default for a conventional CCSD(T) calculation set talg=virt

temp

The temperature in K at which the thermodynamic properties are evaluated (see also keyword freq).

Options:

$<$any positive real number$>$

Default:

temp=298.15

Example:

for 300 K set temp=300.0

test

A keyword for testing Mrcc. If an energy value is specified using this keyword, it will be compared to the energy calculated last time [e.g., the CCSD(T) energy and not the CCSD or HF energy if calc=CCSD(T)] in the Mrcc run. An error message will be displayed and the program exits with an error code if the test energy and the calculated energy differ. This keyword is mainly used by the developers of the program to create test jobs to check the correctness of the computed energies. (See Sect. 8 for the further details.)

Options:
off

No testing.

$<$any real number$>$

The energy to be tested.

Default:

test=off

Example:

to set a test energy of -40.38235315 E${}_{h}$ type test=-40.38235315

tprint

Controls the printing of converged cluster amplitudes/CI coefficients if ccprog=mrcc.

Options:
off

No printing.

$<$any real number$>$

Cluster amplitudes/CI coefficients whose absolute value is greater than this number will be printed.

Note:

The value of the cluster amplitude/CI coefficient and the corresponding spin-orbital labels (serial number of the orbital + a or b for alpha or beta spin orbitals, respectively) will be printed. The numbering of the orbitals corresponds to increasing orbital energy order. Note that orbital energies are printed at the end of the SCF run if verbosity$\geq$3. You can also identify the orbitals using Molden (see Sect. 14.1).

Default:

tprint=off

Example:

to set a threshold of 0.01 give tprint=0.01

uncontract

Uncontract contracted basis sets.

Options:

on or off

Default:

uncontract=off

Example:

to uncontract the basis set add uncontract=on

unit

Specifies the units used for molecular geometries.

Options:
angs

Ångströms will be used

bohr

Atomic units will be used

Default:

unit=angs

Example:

to use bohrs rather than ångströms the user should set unit=bohr

verbosity

Controls the verbosity of the output.

Options:

0, 1, 2, 3. The verbosity of the output increases gradually with increasing value of the option. Error messages are not suppressed at any level.

Default:

verbosity=2

Example:

to increase the amount of information printed out give
verbosity=3

wpairtol

Threshold for the selection of weak pairs in local MP2, RPA, and CC methods. For each orbital pair the estimate of the pair correlation energy is calculated with a multipole approximation [57, 134, 105]. An orbital pair will be considered as distant pair if the absolute value of the multipole-based pair correlation energy estimate is smaller than wpairtol. For the distant pairs the corresponding multipole-based pair correlation energy estimates will be added to the correlation energy, and distant pairs will be neglected in the subsequent calculations.

In the case of localcc=2015, for the remaining pairs a more accurate pair correlation energy estimate will be calculated using orbital specific virtuals (OSVs) controlled by keyword osveps, and these pairs will be further classified as weak and strong pairs controlled by keyword spairtol, see the description of keyword spairtol. The extended domain of an occupied orbital will include those orbitals for which the latter accurate pair correlation energy estimate is greater than spairtol. See also Ref. 85 for more details.

In the case of localcc=2016 or 2018, spairtol=off is set as default, and the extended domain of an occupied orbital will include those orbitals for which the multipole-based pair correlation energy is greater than wpairtol. See also Refs. 105 and 106 for more details.

Note that for local CC methods if spairtol$\neq$off is specified as a non-default option in the case of localcc=2016 or 2018, accurate MP2 pair energies are computed in the extended domains for the remaining non-distant pairs. Then the non-distant pairs are further divided into weak and strong categories according to the value of spairtol, as discussed above. In this case the MP2 pair energies of the weak pairs are added to the correlation energy, and new, somewhat smaller extended domains are constructed to proceed with the higher-level computation as above using solely the strong pair list. See also Ref. 106 for more details.

Options:
$<$any positive real number$>$

Orbital pairs with multipole-based pair correlation energy estimates smaller than this number (in E${}_{h}$) will be considered as distant pairs.

Default:

wpairtol=1e-5 for local MP2 and CC if localcc=2018,
wpairtol=min(1e-6, 0.01*spairtol) for localcc=2015

Note:

For defaults with other than the above settings, see the description of lcorthr

Example:

to set a threshold of $5\cdot 10^{-6}$ E${}_{h}$ type wpairtol=5e-6