2.0cm 2.0cm 2.0cm 2.0cm
Manual for the MRCC
Program System
Release date: August 28, 2023
Department of Physical Chemistry and Materials Science
Budapest University of Technology and Economics
https://www.mrcc.hu/
Mrcc is a suite of ab initio and density functional quantum chemistry programs for high-accuracy electronic structure calculations developed and maintained by the quantum chemistry research group at the Department of Physical Chemistry and Materials Science, TU Budapest [66]. Its special feature, the use of automated programming tools enabled the development of tensor manipulation routines which are independent of the number of indices of the corresponding tensors, thus significantly simplifying the general implementation of quantum chemical methods. Applying the automated tools of the program several quantum chemistry models and techniques of high complexity have been implemented so far including arbitrary single-reference coupled-cluster (CC) and configuration interaction (CI) methods, multi-reference CC approaches, CC and CI energy derivatives and response functions, arbitrary perturbative CC approaches. Many features of the package are also available with relativistic Hamiltonians allowing for accurate calculations on heavy element systems. The developed cost-reduction techniques and local correlation approaches also enable high-precision calculations for medium-sized and large molecules.
In the following, words set in typewriter font denote file names, shell commands, environment variables, and input records of the input file. These must be typed as shown. Variables, that is, numbers, options, etc. which must be specified by the user will be given as $$variable$\mathrm{>}$. These must be replaced by the corresponding values of the variables. Optional items are denoted by brackets, e.g., as [$$variable$\mathrm{>}$].
The main authors of the Mrcc code and their major contributions are as follows.
general design; driver program (dmrcc); input analyzer (minp); automated, string-based many-body code (goldstone, xmrcc, mrcc); domain construction for local correlation calculations (mulli); particular integral evaluation algorithms (integ); direct and density-fitting (DF) Hartree–Fock (HF) algorithms, DFT algorithms (scf); density-fitting MP2 and RPA algorithms (drpa); explicitly correlated models (mp2f12)
local correlation methods and multi-level approaches (drpa); closed-shell coupled-cluster singles and doubles with perturbative triples [CCSD(T)] code (ccsd); parallelization of CCSD(T)
excited-state approaches (cis), range-separated double hybrid functionals (dft), bond functions (integ)
closed-shell CCSD(T) code (ccsd); parallelization of CCSD(T)
quadratic and multi-configurational self-consistent field algorithms, reduced-cost SCF algorithms (scf); analytic gradients for local density fitting SCF (prop)
open-shell CCSD(T) code (uccsd), open-shell local correlation approaches (drpa)
integral transformation and orbital optimization code (ovirt); domain construction for local correlation calculations (mulli)
density-fitting integrals (integ)
installation script (build.mrcc), geometry optimization, basis set optimization, the Mrcc homepage
QM/MM and embedding approaches (qmmod)
solvation models
Hartree–Fock self-consistent field code (scf)
closed-shell CCSD(T) code (ccsd)
four-center atomic orbital integral code (integ)
orbital relaxation contribution to gradients for density-fitting methods, density-fitting MP2 gradient (prop)
orbital relaxation contribution to gradients for coupled-cluster methods (prop)
double hybrid and van der Waals density functionals (dft)
explicitly correlated models (mp2f12)
In addition, Bence Kornis, Levente Dojcsák, Huliyar S. Nataraj, and Sanghamitra Das have also contributed to the development of the Mrcc code.
If results obtained with the Mrcc code are published, the
publication must acknowledge the following two references.
“M. Kállay,
P. R. Nagy, D. Mester, Z. Rolik, G. Samu, J. Csontos, J. Csóka,
P. B. Szabó, L. Gyevi-Nagy, B. Hégely, I. Ladjánszki, L. Szegedy,
B. Ladóczki, K. Petrov, M. Farkas, P. D. Mezei, and Á. Ganyecz:
The MRCC program system: Accurate quantum chemistry from water to
proteins, J. Chem. Phys. 152, 074107 (2020).”
“Mrcc, a quantum chemical program suite written by M. Kállay,
P. R. Nagy, D. Mester, L. Gyevi-Nagy, J. Csóka, P. B. Szabó,
Z. Rolik, G. Samu, J. Csontos, B. Hégely, Á. Ganyecz,
I. Ladjánszki, L. Szegedy, B. Ladóczki, K. Petrov, M. Farkas,
P. D. Mezei, and R. A. Horváth. See www.mrcc.hu.”
In addition, credit must be given to the corresponding papers which describe the implementation and the underlying methodological developments. The corresponding references are given in Sect. 6 of the manual.
If Mrcc is used combined with other program systems or libraries, the users are also requested to include appropriate citations to those packages as required by their authors.
Mrcc can be used as a standalone code, but interfaces have been developed to the Cfour, Columbus, Dirac, Molpro, Orca, and Psi quantum chemistry packages, the xTB and Mopac semi-empirical programs as well as the Amber molecular dynamics (MD) code. Mrcc, in standalone mode, can currently be used for single-point energy calculations, evaluation of first-order properties, geometry optimizations, and harmonic vibrational frequency calculations with the standard nonrelativistic Hamiltonian and effective core potentials. The interfaces enable the calculation of further molecular properties as well as several other features, such as the use of relativistic Hamiltonians and MCSCF orbitals or quantum mechanics/molecular mechanics (QM/MM) simulations. If Mrcc is used together with the aforementioned quantum chemistry packages, the integral, property integral, HF, MCSCF, and CPHF calculations, the integral and density-matrix transformations, etc. are performed by these program systems. Transformed MO (property) integrals are passed over to Mrcc, which carries out the correlation calculation and returns unrelaxed MO density matrices if necessary.
In the following we describe the use of the Cfour, Columbus, Dirac, Molpro, xTB, Mopac, and Amber interfaces and the features that they enable. For a complete list of available features, see Sect. 6. See also the description of keyword iface. For the Orca and Psi interfaces see the manual of these packages.
Most of the implemented features are available via the Cfour interface using RHF, ROHF, and UHF orbitals: single-point energy calculations, geometry optimizations, first-, second-, and third-order property calculations, electronic excitation energies, excited-state and transition properties, diagonal Born–Oppenheimer correction (DBOC) calculations. Most of the properties implemented in Cfour are also available with Mrcc. The interface also enables the use of several relativistic Hamiltonians.
The Cfour interface is very user-friendly. You only have to prepare the input file ZMAT for Cfour with the keyword CC_PROG=MRCC, and run Cfour. The Mrcc input file is then written automatically and Mrcc is called directly by Cfour, and you do not need to write any input file for Mrcc. Most of the features of Mrcc can be controlled by the corresponding Cfour keyword, see Cfour’s homepage at www.cfour.de. If you use the Cfour interface, you can safely ignore the rest of this manual.
You also have the option to turn off the automatic construction of the Mrcc input file by giving INPUT_MRCC=OFF in the Cfour input file ZMAT. However, it is only recommended for expert users.
Single-point energies, equilibrium geometries, ground- and excited-state first-order properties, and transition moments can be computed with RHF, ROHF, and MCSCF reference states using the Columbus interface. Evaluation of harmonic vibrational frequencies is also possible via numerically differentiated analytical gradients.
Running Mrcc with Columbus requires three additional programs, colto55, coldens, and runc_mrcc, which are available for Columbus licensees from the authors of Mrcc upon request.
To use this interface for single-point energy calculations first prepare input files for Columbus using the colinp script. It is important to set a calculation in the input file which requires a complete integral transformation (e.g., CISD and not just MCSCF). Execute Columbus. If you do not need the results of the Columbus calculations, you can stop them after completing the integral transformation. Run the colto55 program in the WORK directory created by Columbus. This will convert the Columbus integral files to the Mrcc format. Prepare input file MINP for Mrcc as described in Sect. 11. Run dmrcc as described in Sect. 9. It is recommended to execute first some inexpensive calculation (e.g., CISD) with Mrcc and compare the HF and CISD energies in order to test your input files.
For property calculations create the Columbus and the Mrcc input files. In the Columbus input set the corresponding MRCI property calculation. Copy the Mrcc input file MINP to the WORK directory of Columbus. If the directory does not exist, create it. Then execute the runc_mrcc script.
The interface to the Dirac code enables four-component relativistic calculations with the full Dirac–Coulomb Hamiltonian and several approximate variants thereof. Single-point energy calculations are possible with all CC and CI methods implemented in Mrcc using Kramers-paired Dirac–Fock orbitals. First-order property (unrelaxed) calculations are available with iterative CC and CI methods. See Refs. 109 and 67 for more details.
If you use Dirac, you should first prepare input files for
Dirac (see http://diracprogram.org/).
It is important to run a full integral transformation with Dirac
(see the description of the MOLTRA keyword in Dirac’s
manual), and to use Abelian symmetry (that is, the
${D}_{2h}$ group and its subgroups).
Execute the pam script saving the MRCONEE and
MDCINT files, e.g., running it as
pam --get="MRCONEE MDCINT" --inp=Y.inp --mol=X.mol
where X.mol and Y.inp should be replaced by your input files
as appropriate. Then run the dirac_mointegral_export interface
program, which generates the files needed by Mrcc.
It also creates a sample input file MINP for Mrcc, which
contains the input for a closed-shell CCSD calculation.
If you intend to run another type of
calculation, please edit the file as described in Sect. 11.
Please also note that you may need to modify the occupation vector
under the refdet keyword (see the description of the keyword on
page 12), and you should set hamilton=x2c if exact
2-component Hamiltonians are used.
Then run dmrcc as described in Sect. 9.
For relativistic property calculations define the corresponding
operator in the Dirac input file (see the description of the
PROPERTIES and MOLTRA keywords in Dirac’s manual).
Then execute the pam script as
pam --get="MRCONEE MDCINT MDPROP" --inp=Y.inp --mol=X.mol
and edit the MINP file, in particular, set the dens
keyword (page 12). The CC property code currently does not
work with double-group symmetry, and you need to turn off symmetry for
CC property calculations, that is, set symm=off in the MINP
file. Finally run dmrcc.
With Molpro single-point energy calculations are available using RHF, UHF, ROHF, and MCSCF orbitals. The interface also enables the use of Douglas–Kroll–Hess Hamiltonians as well as effective core potentials.
The Molpro interface is very user-friendly. You only have to prepare the input file for Molpro with a line starting with the mrcc label followed by the corresponding keywords, and run Molpro. The Mrcc input file is then written automatically and Mrcc is called directly by Molpro, and you do not need to write any input file for Mrcc. Most of the features of Mrcc can be controlled by the corresponding Molpro keywords. If you use Molpro, you also have the option to install Mrcc with the makefile of Molpro.
For a detailed description of the interface point your browser to the Molpro User’s Manual at www.molpro.net and then click “34 The MRCC program of M. Kallay (MRCC)”.
If you use the Molpro interface, you can safely ignore the rest of this manual.
The Amber/Mrcc interface enables QM/MM and other multiscale calculations. The interface of the Amber MD code and Mrcc is based on the work of Götz et al. [35], which facilitates the integration of QM codes into Amber as external modules. The detailed description of the Amber/Mrcc interface is documented in Ref. 47. Currently for the separation of the QM and MM subsystems only the link atom approach is supported.
With the Amber/Mrcc program the projection-based embedding techniques, namely the projector-augmented operator of Manby and Miller [84] and our Huzinaga-operator [48] approaches, can also be employed for the QM region. The multilevel approach based on our local correlation methods is also supported [48]. The latter approaches enables the embedding of wave function theory (WFT) or density functional theory (DFT) methods into lower-level WFT or DFT methods and also the combinations thereof (DFT-in-DFT, WFT-in-DFT, WFT-in-WFT, and WFT-in-WFT-in-DFT). Consequently, with the Amber interface you can also define three (QM/QM/MM) or four (QM/QM/QM/MM) layers for the multilevel calculations. At the present stage of development, the multilevel methods with three or four layers are only available for single point energy calculations. The traditional multilevel approach, “Our own n-layered integrated molecular orbital and molecular mechanics” (ONIOM) of Morokuma et al. [86], is also supported with two QM layers by our in-house interface for single point, geometry optimization, and (BornâOppenheimer) molecular dynamics calculations. The unofficial Mrcc module, which supports ONIOM calculations, is available on user-request.
Using the RISM modules of Amber, the embedded cluster reference interaction site model (EC-RISM) method of Kast et al. [kloss2008quantum, kast2010prediction] was implemented as a solvation model [RISM]. These calculations are controlled by Mrcc, and Amber is used for the 3D-RISM part of the calculations. Based on the integral equation theory, 3D-RISM can provide solvation thermodynamics and distribution functions like molecular dynamics (MD), but in a few minutes, while its accuracy is also satisfying. EC-RISM is a QM/MM-like combination of 3D-RISM and QM to describe systems in solutions. See keyword rism for further details.
To use the Amber/Mrcc interface you need a properly installed version of Amber (version 2017 or later), see the Amber homepage at ambermd.org. The detailed description of the usage of the interface is well documented in the Amber manual, see section “11.2.6.6. AMBER/MRCC” in the manual of Amber21.
The ONIOM method [86] can be also applied with the xTB program [xtbprog], which implements the density functional tight binding (DFTB) methods of Grimme and co-workers [GFN1, GFN2]. These theoretical models are parameterized for geometries, vibrational frequencies, and non-covalent interactions (GFN), and the necessary parameters are available for most of the periodic system. In the ONIOM implementation, single-point and gradient calculations are supported with the mechanical embedding approach, single-point calculations with various point charge models can be utilized for electronic embedding, including the Mulliken charges of the GFN-xTB models, and single-point and gradient calculations can be carried out with electronic embedding if the point charges are user-defined and independent of the geometry. Note that the use of the force-field version of the GFN family (GFN-FF) [GFN-family] and the implicit solvent models of xTB [GFN-ALPB] is also supported in the ONIOM framework. Finally, the GFN1-xTB and GFN2-xTB can also be used for initial guess calculations as it is described in Ref. guess.
The Mrcc program is also interfaced to the Mopac program of Stewart [MOPAC2016], which allows various semi-empirical molecular orbital methods, for instance, the PM6 and PM7 versions of the neglect of diatomic differential overlap (NDDO) approach, to be used with the mechanical embedding version of the ONIOM implementation. Single-point and gradient calculations are supported along with the COSMO implicit solvent model of the Mopac program.
In this section the available features of the Mrcc code are summarized. We also specify what type of reference states (orbitals) can be used, and if a particular feature requires one of the interfaces or is available with Mrcc in standalone mode. We also give the corresponding references which describe the underlying methodological developments.
Single-point energy calculations are possible with the following methods. See also Tables LABEL:GrTable and LABEL:ExTable for the available ground- and excited-state models, respectively, as well as Table LABEL:RCTable for the reduced-cost approaches. For a detailed list of the available density functional approximations, see Table LABEL:FuncTable.
conventional and density-fitting (resolution-of-identity) Hartree–Fock SCF (Ref. 138): restricted HF (RHF), unrestricted HF (UHF), and restricted open-shell HF (ROHF)
conventional and density-fitting (resolution-of-identity) multi-configurational SCF (MCSCF)
conventional and density-fitting (resolution-of-identity) Kohn–Sham (KS) density functional theory (DFT) (Ref. 71): restricted KS (RKS), unrestricted KS (UKS), restricted open-shell KS (ROKS); local density approximation (LDA), generalized gradient approximation (GGA), meta-GGA (depending on kinetic-energy density and/or the Laplacian of the density), hybrid, range-separated hybrid, double hybrid (DH), third-order double hybrid, and range-separated double hybrid (RSDH) functionals (for the available functionals see the description of keyword dft); dispersion corrections
time-dependent HF (TD-HF), time-dependent DFT (TD-DFT), TD-DFT in the Tamm–Dancoff approximation (TDA); currently only for closed-shell molecules with RHF/RKS reference using density-fitting, for LDA and GGA functionals and their hybrids and double hybrids even with range separation (Refs. 91, 90, RSDH, SCS_RSDH, and ADC_RSDH)
density-fitting (resolution-of-identity) MP2, spin-component scaled MP2 (SCS-MP2), and scaled opposite-spin MP2 (SOS-MP2); currently only for RHF and UHF references (Ref. 71)
density-fitting (resolution-of-identity) random-phase approximation (RPA, also known as ring-CCD, rCCD), direct RPA (dRPA, also known as direct ring-CCD, drCCD), second-order screened exchange (SOSEX), renormalized second-order perturbation theory (rPT2), and approximate RPA with exchange (RPAX2); currently the dRPA, SOSEX, rPT2, and RPAX2 methods are available for RHF/RKS and UHF/UKS references, while the RPA method is only implemented for RHF/RKS (Refs. 71 and 72)
density-fitting (resolution-of-identity) second-order coupled-cluster singles and doubles (CC2), configuration interaction singles with perturbative correction for double excitations [CIS(D)], iterative doubles correction to configuration interaction singles [CIS(D${}_{\mathrm{\infty}}$)], and second-order algebraic diagrammatic construction [ADC(2)] approaches; spin-component scaled CC2, CIS(D), CIS(D${}_{\mathrm{\infty}}$), and ADC(2) [SCS-CC2, SCS-CIS(D), SCS-CIS(D${}_{\mathrm{\infty}}$), and SCS-ADC(2)]; scaled opposite-spin CC2, CIS(D), CIS(D${}_{\mathrm{\infty}}$), and ADC(2) [SOS-CC2, SOS-CIS(D), SOS-CIS(D${}_{\mathrm{\infty}}$), and SOS-ADC(2)] (Refs. 93, 94, 90, and 92)
arbitrary single-reference coupled-cluster methods (Ref. 69): CCSD, CCSDT, CCSDTQ, CCSDTQP, …, CC($n$)
arbitrary single-reference configuration-interaction methods (Ref. 69): CIS, CISD, CISDT, CISDTQ, CISDTQP, …, CI($n$), …, full CI
arbitrary perturbative coupled-cluster models (Refs. 63, 8, and 65):
CCSD[T], CCSDT[Q], CCSDTQ[P], …, CC($n-1$)[$n$]
CCSDT[Q]/A, CCSDTQ[P]/A, …, CC($n-1$)[$n$]/A
CCSDT[Q]/B, CCSDTQ[P]/B, …, CC($n-1$)[$n$]/B
CCSD(T), CCSDT(Q), CCSDTQ(P), …, CC($n-1$)($n$)
CCSDT(Q)/A, CCSDTQ(P)/A, …, CC($n-1$)($n$)/A
CCSDT(Q)/B, CCSDTQ(P)/B, …, CC($n-1$)($n$)/B
CCSD(T)${}_{\mathrm{\Lambda}}$, CCSDT(Q)${}_{\mathrm{\Lambda}}$, CCSDTQ(P)${}_{\mathrm{\Lambda}}$, …, CC($n-1$)($n$)${}_{\mathrm{\Lambda}}$
CCSDT-1a, CCSDTQ-1a, CCSDTQP-1a, …, CC($n$)-1a
CCSDT-1b, CCSDTQ-1b, CCSDTQP-1b, …, CC($n$)-1b
CC2, CC3, CC4, CC5, …, CC$n$
CCSDT-3, CCSDTQ-3, CCSDTQP-3, …, CC($n$)-3
explicitly correlated MP2 (MP2-F12) using the 2B ansatz, the F + K commutator approximation, fixed amplitudes, spin-flipped geminals in open-shell calculations (Ref. F12pT)
explicitly correlated CCSD and CCSD(T) utilizing the CCSD(F12*)(T+) model and the above MP2-F12 approximation (Ref. F12pT)
multi-reference CI approaches, Davidson correction for MRCI (Ref. 70)
multi-reference CC approaches using a state-selective ansatz (Ref. 70)
arbitrary single-reference linear-response (equation-of-motion, EOM) CC methods (Ref. 62): LR-CCSD (EOM-CCSD), LR-CCSDT (EOM-CCSDT), LR-CCSDTQ (EOM-CCSDTQ), LR-CCSDTQP (EOM-CCSDTQP), …, LR-CC($n$) [EOM-CC($n$)]
linear-response (equation-of-motion) MRCC schemes (Ref. 62)
DFT/WFT embedding (Ref. 48): DFT-in-DFT, WFT-in-DFT, WFT-in-WFT, and WFT-in-WFT-in-DFT embedding (where WFT stands for wave function theory); currently LDA, GGA, meta-GGA, or hybrid functionals can be used as the DFT method; any WFT method implemented in Mrcc can by used in WFT-in-DFT type embedding calculations; WFT-in-WFT multi-level methods (via keyword corembed) is only available for local correlation methods.
ONIOM approach with arbitrary number of layers [86]. Currently, only the integrated MO+MO (IMOMO) [57] part of ONIOM is available in the standalone version of Mrcc with mechanical embedding, electrostatic embedding, and automated link atom handling. Note that all types of the aforementioned DFT/WFT embedding techniques are also supported as top-level ONIOM calculations.
The CI and CC approaches listed above are available with RHF, UHF, standard/semi-canonical ROHF, MCSCF orbitals. Density-fitting with CI or high-order CC is currently enabled only for RHF references. For the perturbative CC approaches with ROHF reference determinant, for theoretical reasons, semi-canonical orbitals are used (see Ref. 65).
The CI and CC approaches listed above are also available with the following interfaces and references.
Cfour, Columbus, and Molpro
Cfour, Columbus, and Molpro
Cfour
Cfour, and Molpro
Columbus and Molpro
QM/MM calculations can be performed by the Amber interface using any method implemented in Mrcc for the QM region.
Semi-empirical quantum mechanical calculations are available through the Mopac and xTB interfaces in the ONIOM framework with mechanical embedding using automated link atom handling. The electrostatic embedding version of the ONIOM scheme is only available with the xTB program.
Single-point calculations are also possible with several types of relativistic Hamiltonians and reference functions, see Sect. 6.7 for more details.
CC$n$ calculations with ROHF orbitals are not possible for theoretical reasons, see Ref. 65 for explanation.
Single-point CI and CC calculations are, in principle, possible with RKS, UKS, ROKS orbitals.
Solvation effects can be modeled by the polarizable continuum model (PCM) [147] via an interface to the PCMSolver library [16, 17, 1]. The PCM treatment is self-consistent for HF and KS SCF calculations, while, in the post-SCF steps, the potential of the solvent optimized at the SCF level is frozen. See the descriptions of keywords pcm*. An improved treatment of solvent effects is possible using the embedded cluster reference interaction site model (EC-RISM) [RISM], which is available via the Amber interface, see Sect. 5.5 and keyword rism for further details. Note that the implicit solvation models of the Mopac and xTB programs are also available in the ONIOM framework.
Geometry optimizations and first-order property calculations can be performed using analytic gradient techniques with the following methods. See also Tables LABEL:GrTable and LABEL:ExTable for the availability of analytic gradients for ground- and excited-state models, respectively.
conventional and DF (RI) HF-SCF (Ref. 138): RHF and UHF
conventional and DF (RI) DFT (Ref. 71): RKS and UKS with LDA, GGA, meta-GGA (depending only on kinetic-energy density), and hybrid functionals as well as dispersion corrections
double hybrid density functional methods (Ref. 71), such as B2PLYP, B2PLYP-D3, B2GPPLYP, etc. (current limitations: only MP2 correlation, closed shell RKS, no spin-component scaling, no meta-GGA functionals, only DH functionals for which the DFT contribution to the energy is stationary with respect to the variation of the MO coefficients)
DF-MP2 (RI-MP2), currently only for RHF references (Ref. 71)
CIS and TD-HF methods using RHF references
The ONIOM approach with arbitrary number of layers [86]. Currently, only the IMOMO [57] part of ONIOM is available in the standalone version of Mrcc with mechanical embedding and automated link atom handling. ONIOM gradients with electronic embedding are also available with perturbation independent point charges.
Currently only unrestricted geometry optimizations are possible, and electric dipole, quadrupole, and octapole moments as well as the electric field at the atomic centers can be evaluated. In addition, Mulliken, Löwdin, intrinsic atomic orbital (IAO), CHELPG, and Merz–Singh–Kollman atomic charges, and Mayer bond orders can be computed using the SCF wave functions. Analytic gradients for the CI and CC methods listed above are available with RHF, UHF, and standard ROHF orbitals without density fitting.
The following keywords are available to control the optimization process
– to select an algorithm for the optimization
– maximum number of iterations allowed
– convergence criterion for energy change
– convergence criterion for the gradient change
– convergence criterion for the step-size
The optimization will be terminated and regarded as successful when the maximum gradient component becomes less than optgtol and either an energy change from the previous step is less than optetol or the maximum displacement from the previous step is less than optstol. For their detailed description see Sect. 12.
The implemented analytic gradients for the CI and CC approaches listed above can also be utilized via the Cfour and Columbus interfaces with the following references.
Cfour and Columbus
Cfour and Columbus
Cfour
Columbus
In addition to geometries, most of the first-order properties (dipole moments, quadrupole moments, electric field gradients, relativistic contributions, etc.) implemented in Cfour and Columbus can be calculated with Mrcc.
QM/MM geometry optimizations and MD calculations can be performed by the Amber interface using any method implemented in Mrcc for which analytic gradients are available.
Analytic gradients are also supported for the semi-empirical quantum mechanics methods of the Mopac and xTB programs with the mechanical embedding of the ONIOM framework. The electronic embedding with perturbation independent point charges is also supported for the DFTB methods of the xTB program.
Geometry optimizations and first-order property calculations can also be performed via numerical differentiation for all methods available in Mrcc using the Cfour interface.
Analytic gradients are are also available with several types of relativistic Hamiltonians and reference functions, see Sect. 6.7 for more details.
The implicit solvation models of the xTB and Mopac programs and their analytic gradient are available in the ONIOM implementation. See the descriptions of keywords pcm* and oniom_pcm.
Harmonic vibrational frequencies, infrared (IR) intensities, and ideal gas thermodynamic properties can be evaluated using numerically differentiated analytic gradients for all the methods listed in Sect. 6.2.
CC and CI harmonic frequency and second-order property calculations for RHF and UHF references can also be performed using analytic second derivatives (linear response functions) with the aid of the Cfour interface. Analytic Hessians (LR functions) are available for the following approaches.
In addition to harmonic vibrational frequencies [61], the analytic Hessian code has been tested for NMR chemical shifts [61], static and frequency-dependent electric dipole polarizabilities [64], magnetizabilities and rotational $g$-tensors [32], electronic $g$-tensors [31], spin-spin coupling constants, and spin rotation constants. These properties are available via the Cfour interface.
Using the Cfour interface, harmonic frequency calculations are also possible via numerical differentiation of energies for all implemented methods with RHF, ROHF, and UHF orbitals.
Using the Cfour or the Columbus interface, harmonic frequency calculations are also possible via numerical differentiation of analytic gradients for all implemented methods for which analytic gradients are available (see Sect. 6.2 for a list of these methods). With Cfour, the calculation of static polarizabilities is also possible using numerical differentiation.
NMR chemical shifts can be computed for closed-shell molecules using gauge-including atomic orbitals and RHF reference function.
Third-order property calculations can be performed using analytic third derivative techniques (quadratic response functions) invoking the Cfour interface for the following methods with RHF and UHF orbitals.
The analytic third derivative code has been tested for static and frequency-dependent electric-dipole first (general, second-harmonic-generation, optical-rectification) hyperpolarizabilities [114] and Raman intensities [113]. Please note that the orbital relaxation effects are not considered for the electric-field. These properties are available via the Cfour interface.
Using the Cfour interface, anharmonic force fields and the corresponding spectroscopic properties can be computed using numerical differentiation techniques together with analytic first and/or analytic second derivatives at all computational levels for which these derivatives are available (see Sect. 6.2 and 6.3 for a list of these methods).
Diagonal Born–Oppenheimer correction (DBOC) calculations can be performed using analytic second derivative techniques via the Cfour interface for the following methods with RHF and UHF references.
Excitation energies, first-order excited-state properties, and ground to excited-state transition moments can be computed as well as excited-state geometry optimizations can be performed using linear response theory and analytic gradients with the following methods.
Excitation energy and property calculations for the aforementioned methods are available with RHF, UHF, and standard ROHF orbitals. Density-fitting is only possible for RHF-based single point calculations. So far electric and magnetic dipole transition moments, both in the length and the velocity gauge, as well as the corresponding oscillator and rotator strengths have been implemented. For the list of implemented first-order properties see Sect. 6.2.
Excitation energies can also be computed for closed-shell systems using the density-fitting approximation and RHF (RKS) orbitals with the following methods (Refs. 93, 94, 91, 90, RSDH, SCS_RSDH, and ADC_RSDH).
CIS, time-dependent HF (TD-HF), Tamm–Dancoff approximation (TDA), time-dependent DFT (TD-DFT)
configuration interaction singles with perturbative correction for double excitations [CIS(D)]
second-order coupled-cluster singles and doubles (CC2) method
iterative doubles correction to configuration interaction singles [CIS(D${}_{\mathrm{\infty}}$)] method
second-order algebraic diagrammatic construction [ADC(2)] approach
spin-scaled versions of the latter approaches: SCS-CIS(D), SCS-CC2, SCS-CIS(D${}_{\mathrm{\infty}}$), SCS-ADC(2), SOS-CIS(D), SOS-CC2, SOS-CIS(D${}_{\mathrm{\infty}}$), and SOS-ADC(2)
CIS(D)- and ADC(2)-based double hybrid TD-DFT methods even with range separation
Analytic gradients are implemented for CIS and TD-HF. Ground to excited-state transition moments are available for TD-HF, TDA, TD-DFT, hybrid and double-hybrid TD-DFT even with range separation, CIS(D), SCS-CIS(D), SOS-CIS(D), ADC(2), SCS-ADC(2), and SOS-ADC(2). For the CIS, TD-HF, hybrid TDA and TD-DFT, CC2, CIS(D${}_{\mathrm{\infty}}$), ADC(2), SCS-CC2, SCS-CIS(D${}_{\mathrm{\infty}}$), SCS-ADC(2), SOS-CC2, SOS-CIS(D${}_{\mathrm{\infty}}$), and SOS-ADC(2) methods efficient reduced-cost approaches are also implemented, see Sect. 6.8 for details.
Excitation energies, first-order excited-state properties, and ground to excited-state transition moments can also be calculated as well as excited-state geometry optimizations can also be carried out using the following interfaces and reference states.
Cfour, Columbus, and Molpro (only excitation energy)
Cfour, Columbus, and Molpro (only excitation energy)
Cfour and Molpro (only excitation energy)
Columbus and Molpro (only excitation energy)
Please note that for excitation energies and geometries LR-CC methods are equivalent to the corresponding EOM-CC models. It is not true for first-order properties and transition moments.
With CI methods, excited to excited-state transition moments can also be evaluated.
Excited-state harmonic frequencies can be evaluated for the above methods with the help of numerical differentiation of analytical gradients, see Sect. 6.3.
Excited-state harmonic frequencies can also be calculated for the above methods via numerical differentiation using the Cfour or Columbus interface.
Excited-state harmonic frequencies and second-order properties can be evaluated for CI methods using analytic second derivatives and the Cfour interface.
Natural transition orbitals (NTOs) can be computed for the CIS, TDHF, TDA, TD-DFT, ADC(2), CC2, and CIS(D${}_{\mathrm{\infty}}$) methods, while state-averaged natural orbitals (NOs) can be computed for the second-order methods. See the description of keyword nto for details.
Treatment of special relativity in single-point energy calculations is possible for all the CC and CI methods listed in Sect. 6.1 using various relativistic Hamiltonians with the following interfaces (Refs. 109 and 67).
With Molpro, relativistic calculations can be performed with Douglas–Kroll–Hess Hamiltonians using RHF, UHF, ROHF, and MCSCF orbitals. The interface also enables the use of effective core potentials (see Molpro’s manual for the specification of the Hamiltonian and effective core potentials).
With Cfour, exact two-component (X2C) and spin-free Dirac–Coulomb (SF-DC) calculations can be performed. The evaluation of mass-velocity and Darwin corrections is also possible using analytic gradients for all the methods and reference functions listed in Sect. 6.2. (See the description of the RELATIVISTIC keyword in the Cfour manual for the specification of the Hamiltonian.)
Treatment of special relativity in analytic gradient calculations is possible for all the CC and CI methods listed in Sect. 6.2 using various relativistic Hamiltonians with the following interfaces.
With Cfour, analytic gradient calculations can be performed with the exact two-component (X2C) treatment.
The methods for which reduced-cost algorithms are available are collected in Table LABEL:RCTable.
The different types of SCF calculations can be sped up utilizing various tricks. In the dual basis set methods [81, 49], the SCF equations are solved with a smaller basis set, the resulting density matrix is projected into a large basis set, and then, the Fock matrix is constructed and diagonalized in the large basis (see keyword dual). In a similar manner, dual auxiliary bases, dual grids, or dual fitting metrics can also be used [DualDF] (see keywords dual_df, grid_sm, and fitting).
For the reduction of the computational expenses of the exchange contribution to the Fock matrix, several techniques have been implemented. The evaluation of three-center Coulomb integrals and their contribution can be accelerated with the multipole approximation [Multipole] (see keywords fmm, fmmord). The exchange contribution can also be computed with semi-numerical techniques exploiting pseudo-spectral approximations [PSSPFMM, 110, DualDF], for which also a dual grid approach is available [DualDF] (see keywords pssp, agrid_pssp, agrid_pssp_sm). Probably the most efficient method to speed up SCF calculations is the local density fitting approach, where products of orbitals are expanded in local domains of fitting functions [127, 138, 107] (see keywords scfalg, excrad, excrad_fin). The efficiency of the latter can further be improved with the OCC-RI-K algorithm [occri, DualDF] algorithm, which can also be used together with non-local direct DF-SCF calculations (see keyword occri).
The computational expenses of the CC and CI methods listed in Sect. 6.1 can be reduced via orbital transformation techniques (Ref. 136). In this framework, to reduce the computation time, the dimension of the properly transformed occupied and/or virtual one-particle space is truncated. Currently optimized virtual orbitals (OVOs) or MP2 natural orbitals (NOs) can be chosen. This technique is recommended for small to medium-size molecules. This scaling reduction approach is available using RHF or UHF orbitals for virtual orbitals. Occupied NOs are only implemented for RHF orbitals and with density fitting. Spatial symmetry for the frozen NO (FNO) methods can also be utilized with the general-order CC and CI methods of the mrcc program. For the correction of the dropped orbitals, besides the conventional $\mathrm{\Delta}$MP2 approach, novel PPL+ and (T+) corrections [FNOExt] can also be used. See the description of keywords ovirt, eps, lnoepsv, and ovosnorb for more details. Reduced-cost ground-state correlation methods can also be activated by adding the FNO- prefix to the acronym of the method [e.g., calc=FNO-CCSD(T)], see notes for keyword calc.
The cost of density-fitting methods can be reduced using natural auxiliary functions (NAFs) introduced in Ref. 71. The approach is very efficient for dRPA, but considerable speedups can also be achieved for MP2, CC2, and ADC(2). See the description of keywords naf_cor and naf_scf for more details.
The computational expenses of CIS, TD-HF, hybrid TDA and TD-DFT, CC2, CIS(D${}_{\mathrm{\infty}}$), ADC(2), SCS-CC2, SCS-CIS(D${}_{\mathrm{\infty}}$), SCS-ADC(2), SOS-CC2, SOS-CIS(D${}_{\mathrm{\infty}}$), and SOS-ADC(2) excited-state calculations can be efficiently reduced using local fitting domains as well as state-averaged NOs and NAFs (Refs. 93, 94, and 91). The scaling of CC2 and ADC(2) calculations can also be decreased by local correlation approaches (Ref. 90). See the description of keywords redcost_exc and redcost_tddft for more details.
The cost of MP2, dRPA, SOSEX, corresponding double- and dual-hybrid DFT, as well as single-reference iterative and perturbative coupled-cluster calculations can be reduced for large molecules by the local natural orbital CC (LNO-CC) approach (Refs. 135, 138, 72, 107, ROLocalMP2, 105, 108, and 106). This method combines ideas from the cluster-in-molecule approach of Li and co-workers [80], the incremental approach of Stoll et al. [144], domain- and pair approximations introduced first by Pulay et al. (see, e.g., Ref. 129) with frozen natural orbital, natural auxiliary function, and Laplace transform techniques. All of the above methods are available for closed-shell molecules using RHF (RKS) orbitals, while currently only LMP2 [ROLocalMP2] is implemented using ROHF (ROKS) orbitals or QROs formed from unrestricted orbitals. The LMP2 and LNO-CCSD(T) implementations are highly competitive, applicable up to a thousand atoms or about 50,000 atomic orbitals [108, 106] with accurate settings and basis sets. Highly-converged LNO settings and quadruple- or quintuple-$\zeta $ basis sets are feasible up to a few hundred atoms [L7CCDMC, 26, MeLysine] depending on the system type.
Further unique features of our local correlation methods include the exceptionally small memory and disk use, frequent checkpointing and restarting capability (see lccrest), support for general-order LNO-CC approaches, and treatments of point group symmetry (see localcorrsymm), and near-linear dependent AO basis sets [108, 106]. See the description of keywords localcc, lnoepso, lnoepsv, domrad, lmp2dens, dendec, nchol, osveps, spairtol, wpairtol, laptol, lccrest, and lcorthr for further details. Extensive benchmarks regarding the accuracy and efficiency of the local correlation methods are also provided in Refs. 72, 107, ROLocalMP2, 108, and 106 and Section II.G of Ref. 66. The definition, performance analysis, and the illustration of the systematic convergence of the lcorthr settings, as well as of the corresponding LNO extrapolation scheme is documented in Ref. 106. Extensive independent LNO-CCSD(T) benchmarks and comparisons to alternative local CCSD(T) approaches are also available [MobiusCC, LCCcompRuMartin, 117, LCC_MOBH35re, LCCcomp_anionbind].
Utilizing the above local correlation techniques a multi-level scheme is defined in which the LMOs are classified as active or environment (Ref. 48, 49). The contributions of these LMOs to the total correlation energy are evaluated using different models for the two subsystems, for instance, one can choose a LNO-CC model for the active subsystem and LMP2 for the environment. See the description of keyword corembed for further details.
The optimization of basis set’s exponents and contraction coefficients can be performed with any method for which single-point energy calculations are available (see Sect. 6.1). The implementation is presented in Ref. 89. The related keywords are
– to turn on/off basis set optimization
– to select an algorithm for the optimization
– maximum number of iterations allowed
– convergence criterion for energy change
– convergence criterion for parameter (exponent, contraction coefficient) change
For their detailed description see Sect. 12.
For the optimization of basis sets it is important to know the format for the storage of the basis set parameters. In Mrcc the format used by the Cfour package is adapted. The format is communicated by the following example.
actual lines | description | |||
---|---|---|---|---|
C:6-31G | $\hookleftarrow $ Carbon atom:basis name | |||
Pople’s Gaussian basis set | $\hookleftarrow $ comment line | |||
$\hookleftarrow $ blank line | ||||
2 | $\hookleftarrow $ number of angular momentum types | |||
0 | 1 | $\hookleftarrow $ 0$\to $s , 1$\to $p | ||
3 | 2 | $\hookleftarrow $ number of contracted functions | ||
10 | 4 | $\hookleftarrow $ number of primitives | ||
$\hookleftarrow $ blank line | ||||
3047.5249 | 457.36952 | … | $\hookleftarrow $ exponents for s functions | |
$\hookleftarrow $ blank line | ||||
0.0018347 | 0.0000000 | 0.0000000 | $\hookleftarrow $ contraction coefficients | |
0.0140373 | 0.0000000 | 0.0000000 | for s functions | |
0.0688426 | 0.0000000 | 0.0000000 | ||
0.2321844 | 0.0000000 | 0.0000000 | ||
0.4679413 | 0.0000000 | 0.0000000 | ||
0.3623120 | 0.0000000 | 0.0000000 | ||
0.0000000 | 0.1193324 | 0.0000000 | ||
0.0000000 | 0.1608542 | 0.0000000 | ||
0.0000000 | 1.1434564 | 0.0000000 | ||
0.0000000 | 0.0000000 | 1.0000000 | ||
$\hookleftarrow $ blank line | ||||
7.8682724 | 1.8812885 | … | $\hookleftarrow $ exponents for p functions | |
$\hookleftarrow $ blank line | ||||
0.0689991 | 0.0000000 | $\hookleftarrow $ contraction coefficients | ||
0.3164240 | 0.0000000 | for p functions | ||
0.7443083 | 0.0000000 | |||
0.0000000 | 1.0000000 |
In a basis set optimization process you need two files in the working directory: the appropriate MINP file with the basopt keyword set and a user supplied GENBAS file that contains the basis set information in the above format. You do not need to write the GENBAS file from scratch, you can use the files in the BASIS directory of Mrcc to generate one or you can use the Basis Set Exchange [5, 22, 142, 128] to download a basis in the appropriate form (CFOUR format). Note that you can optimize several basis sets at a time: all the basis sets which are added to the GENBAS file will be optimized.
You can perform unconstrained optimization when all the exponents and contraction coefficients are optimized except the ones which are exactly 0.0 or 1.0. Alternatively, you can run constrained optimizations when particular exponents/coefficients or all exponents and coefficients for a given angular momentum quantum number are kept fixed during the optimization. The parameters to be optimized can be specified in the GENBAS file as follows.
Unconstrained optimization: no modifications are needed—by default all exponents and contraction coefficients will be optimized except the ones which are exactly 0.0 or 1.0.
Constrained optimization: by default all the exponents and coefficients will be optimized just as for the unconstrained optimization. To optimize/freeze particular exponents or coefficients special marks should be used:
use the “--” mark (without quotes) if you want to keep an exponent or coefficient fixed during the optimization. You should put this mark right after the fixed parameter (no blank space is allowed). If this mark is attached to an angular momentum quantum number, none of the exponents/coefficients of the functions in the given shell will be optimized except the ones which are marked by “++”.
use the “++” mark (without quotes) if you want a parameter to be optimized. Then you should put this mark right after it (no blanks are allowed). You might wonder why this is needed if the default behavior is optimization. Well, this makes life easier. If you want to optimize just a few parameters, it is easier to constrain all parameters first then mark those, which are needed to be optimized (see the example below).
Examples:
To reoptimize all parameters in the above basis set but the exponents
and coefficients of s-type functions you should copy the basis set to
the GENBAS file and put mark “--” after the angular
momentum quantum number of 0. The first lines of the GENBAS file:
C:6-31G
Pople’s Gaussian basis set
2
0--
1
3
2
10
4
3047.5249
457.36952
…
Both s- and p-type functions are fixed but the first s-exponent:
C:6-31G
Pople’s Gaussian basis set
2
0--
1--
3
2
10
4
3047.5249++
457.36952
…
During the optimization the GENBAS file is continuously updated, and if the optimization terminated successfully, it will contain the optimized values (in this case it is equivalent to the GENBAS.opt file, see below, the only difference is that the file GENBAS.opt may contain the special marks, i.e., “++”, “--”). Further files generated in the optimization are:
GENBAS.init – the initial GENBAS file saved
GENBAS.tmp – temporary file, updated after each iteration, can be used to restart conveniently a failed optimization process
GENBAS.opt – this file contains the optimized parameters after a successful optimization.
After registration at the Mrcc homepage, pre-built binaries are
available in the download area. These binaries were compiled with the Intel
compiler (and Intel MKL) version 2021.10.0 and should utilize optimal instruction sets
(e.g., AVX-512) of modern CPUs. If your CPU is still not correctly identified,
the MKL_ENABLE_INSTRUCTIONS environment variable enables you to
use an architecture-specific code path of your choice in Intel MKL routines. To
install these executables, Linux operating system and the 2.27 or later version of the GNU C Library (glibc) is required.
To use the MPI-parallel executables, Intel MPI 2021 also has to be installed. The precompiled binaries are linked against the Intel MPI 2021.10 library, and
it is highly recommended to use this particular Intel MPI version with the
precompiled binaries. It is also strongly suggested to install the newest stable
version of
Libfabric (1.9.0 or later) from the https://github.com/ofiwg/libfabric
repository as some of the previous versions provided via the Intel MPI package could
cause irregular runtime behavior. The binaries are provided in a gzipped tar file,
mrcc.YYYY-MM-DD.binary.tar.gz, where YYYY-MM-DD is the
release date of the program. Note that you will
find several program versions on the homepage. Unless there are
overriding reasons not to do so, please always download the last
version. To unpack the file type
tar xvzf mrcc*.binary.tar.gz
Please do not forget to add the name of the directory where the
executables are placed to your PATH environment variable.
To install Mrcc from source code some version of the Unix operating system, Fortran 90 and C compilers as well as BLAS (basic linear algebra subprograms) and LAPACK (linear algebra package) libraries are required. The tested BLAS and LAPACK implementations are Intel MKL (oneMKL), AOCL-BLIS, and AOCL-libFLAME. Optionally, Mrcc can also be linked with the Libxc library of density functionals [85, 56] and the PCMSolver library for continuum solvation [17, 1]. For an MPI-parallel build, a working MPI installation is also required. Please be sure that the directories where the compilers are located are included in your PATH environment variable. Please also check your LD_LIBRARY_PATH environment variable, which must include the directories containing the BLAS and LAPACK and, if linked against, the Libxc, PCMSolver, and MPI libraries.
After registration at the Mrcc homepage the program can be
downloaded as a gzipped tar file, mrcc.YYYY-MM-DD.tar.gz, where
YYYY-MM-DD is the release date of the program. Note that you will
find several program versions on the homepage. Unless there are
overriding reasons not to do so, please always download the last
version. To unpack the file type
tar xvzf mrcc*.tar.gz
To install Mrcc add the current working directory to PATH
(e.g., export PATH=.:$PATH on bash) and run the build.mrcc script as
build.mrcc [$$compiler$\mathrm{>}$] [-i$$option1$\mathrm{>}$] [-p$$option2$\mathrm{>}$] [-g] [-d] [-s]
[-f$$folder$\mathrm{>}$]
[-l$$library$\mathrm{>}$]
[-L$$lspath$\mathrm{>}$]
[-nomkl]
$$compiler$\mathrm{>}$ specifies the compiler to be used. The options are
as follows. Note that currently only the Intel compiler system is
supported. The freely available Intel oneAPI compilers (with the base and hpc
toolkits) are recommended if an existing Intel compiler installation
is not available.
Previously, the other listed compiler systems were also
used to generate working binaries.
Intel
Intel compiler
GNU
GNU compiler (g77 or gfortran)
PGF
Portland Group Fortran compiler
G95
G95 Fortran 95 compiler
PATH
Pathscale compiler
HP
HP Fortran Compiler
DEC
Compaq Fortran Compiler (DEC machines)
XLF
XL Fortran Compiler (IBM machines)
Solaris10
Sun Solaris10 and Studio10 Fortran Compiler (AMD64)
If the build.mrcc script is invoked without specifying the
$$compiler$\mathrm{>}$ variable, a help message is displayed.
Optional arguments:
specifies if 32- or 64-bit integer variables are used.
Accordingly, $$option1$\mathrm{>}$ can take the value of 32 or
64. The 32-bit integer version is not supported any more.
Default: 64 for 64-bit machines, 32 otherwise.
generates parallel code using message passing interface
(MPI) or open multi-processing (OpenMP) technologies. Accordingly, $$option2$\mathrm{>}$ can take
the MPI[=$$MPI implementation$\mathrm{>}$] or OMP values.
OpenMP and MPI parallelizations have been tested with the Intel compiler
and the Intel MPI and Open MPI implementations.
If -pMPI is specified, the library given by $$MPI implementation$\mathrm{>}$
will be linked to Mrcc.
The default value of $$MPI implementation$\mathrm{>}$ is IntelMPI
with the Intel compiler and OpenMPI for other compilers.
Please note that currently the
two parallelization schemes can only be combined for the scf,
mrcc, and ccsd programs,
other executables will be compiled with only OpenMP parallelization
even if -pOMP and -pMPI are both set.
Default: no parallelization.
source codes are compiled with debugging option
(use this for development purposes)
Default: no debugging option.
source codes are compiled for development, no
optimization is performed (use this for development purposes)
Default: codes are compiled with highest level optimization.
specifies the installation folder. Executables, basis set libraries, and test jobs will be copied to directory $$folder$\mathrm{>}$. If this flag is not used, you will find the executables, etc. in the directory where you perform the installation.
Mrcc is linked with an external library. The libxc and pcm options for $$library$\mathrm{>}$ require the installation of the Libxc and the PCMSolver libraries, respectively. See notes below. Other values for -l$$library$\mathrm{>}$ are passed to the linker unchanged.
Adds $$lspath$\mathrm{>}$ to the library search path while linking.
Turns off linking with Intel MKL library when using the Intel compiler. If this option is used, the BLAS and LAPACK libraries must be specified manually using -l, see example below. It is recommended to use Intel MKL library for best performance.
Mrcc is linked statically. Non-MPI-parallel executables are linked entirely statically, while MPI-parallel executables link Intel-provided libraries statically.
Notes:
After the installation please do not forget to add the directory where the Mrcc executables are located to your PATH environment variable. This is the $$folder$\mathrm{>}$ directory if you used the -f flag at the installation, otherwise the directory where you executed the build.mrcc script.
The build.mrcc script has been tested on several platforms with several versions of the supported compilers and libraries. Nevertheless you may need to customize the compiler flags, names of libraries, etc. These data can be found in the build.mrcc.config file, please edit this file if necessary. Please do not change build.mrcc.
To ensure the best performance of the software, the use of Intel compiler is recommended. The current release of Mrcc has been tested with the 2021.10 version of the Intel compiler (freely available from Intel).
If you use Mrcc together with Molpro, you can also use the Molpro installer to install Mrcc, please follow the instructions in the Molpro manual (www.molpro.net).
If Mrcc is linked with the Libxc library, Libxc must be installed before starting the installation of Mrcc. Note that you must compile Libxc with the same Fortran compiler as used for the installation of Mrcc. At the installation of Libxc, it is recommended to set the installation path of Libxc (--prefix=<install dir> option of configure) to the directory where the installation of Mrcc is carried out, otherwise please set the LIBS_LIBXC environment variable to the installation path of Libxc (i.e. export LIBS_LIBXC=<install dir> in bash, where <install dir> contains lib/libxcf03.a and lib/libxc.a) before running build.mrcc. See the manual for the Libxc project for details [56], as well as the examples below. The current release of Mrcc has been tested with the 6.2.2 version of Libxc.
If Mrcc is linked with the PCMSolver library, PCMSolver must be installed before starting the installation of Mrcc. At the installation of PCMSolver, it is recommended to set the installation path of PCMSolver (--prefix=<install dir> option of setup.py) to the directory where the installation of Mrcc is carried out, otherwise please set the LIBS_PCM environment variable to the installation path of PCMSolver (i.e., export LIBS_PCM=<install dir> in bash, where <install dir> contains lib/libpcm.so, or lib/libpcm.a in the case of static linking) before running build.mrcc. See the manual for the PCMSolver project for details [1], as well as the examples below. The current release of Mrcc has been tested with the 1.3.0 version of PCMSolver.
If the program is compiled for multi-node parallel execution, MPI-parallel executables (*_mpi) are generated. The compilation is performed by mpiifort when $$MPI implementation$\mathrm{>}$ = IntelMPI is set, otherwise the mpifort compiler wrapper is used.
For the compilation of MPI-parallel executables, a working MPI installation is necessary. Currently Open MPI version 4 and Intel MPI (2017 or later) implementations are supported. Open MPI has to be patched with commits 07830d0 and 51acbf7 from https://github.com/open-mpi/ompi. Please, consult Sect. 9.3 for additional MPI settings required at runtime.
Examples:
Compile Mrcc for OpenMP parallel execution with the Intel compiler
(recommended):
build.mrcc Intel -pOMP
Compile Mrcc for OpenMP parallel execution with the Intel compiler
and install it to the /prog/mrcc directory (recommended):
build.mrcc Intel -pOMP -f/prog/mrcc
Compile Mrcc for OpenMP and combined OpenMP-MPI parallel execution with the Intel compiler and Intel MPI, and link with the Libxc and the PCMSolver libraries supposing that Mrcc is compiled in the /prog/mrcc directory. This will enable all features of Mrcc and is highly recommended.
Installation of the Libxc library:
Download the Libxc library (libxc*.tar.gz) from the homepage
of the Libxc project [56].
tar xvzf libxc*.tar.gz
cd libxc*
./configure --prefix=/prog/mrcc/ FC=ifort --enable-kxc
make
make check
make install
cd /prog/mrcc
Installation of the PCMSolver library:
curl -L https://github.com/PCMSolver/pcmsolver/archive/v1.3.0.tar.gz | tar -xz
cd pcmsolver-1.3.0
Replace -openmp by -qopenmp in file cmake/downloaded/autocmake_omp.cmake
cp cmake/custom/compilers/Intel/C.make
cmake/custom/compilers/Intel/C.cmake
./setup.py --cxx=icpc --cc=icc --fc=ifort --int64 --omp --prefix=/prog/mrcc/
cd build
make
make install
export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:/prog/mrcc/lib
cd /prog/mrcc
cp pcmsolver-1.3.0/api/pcmsolver.f90 .
Compiling and linking Mrcc:
build.mrcc Intel -pOMP -pMPI=IntelMPI -llibxc -lpcm
Compile Mrcc for serial execution with the Intel compiler:
build.mrcc Intel
Compile Mrcc for parallel execution using MPI environment
with the Intel compiler for 32-bit machines:
build.mrcc Intel -i32 -pMPI
Compile Mrcc with the Intel compiler for parallel execution
using OpenMP and MPI parallelization through the Open MPI library:
build.mrcc Intel -pOMP -pMPI=OpenMPI
Compile Mrcc with the Intel compiler using BLAS and LAPACK libraries
other than Intel MKL installed in a standard directory:
build.mrcc Intel -nomkl -lblas -llapack
Compile Mrcc with the Intel compiler using BLAS and LAPACK libraries
other than Intel MKL, the shared object files libblas.so and liblapack.so can be found in
$$path-to-blas$\mathrm{>}$ and $$path-to-lapack$\mathrm{>}$, respectively:
build.mrcc Intel -nomkl -L$$path-to-lapack$\mathrm{>}$ -llapack -L$$path-to-blas$\mathrm{>}$ -lblas
Under the Windows operating system the pre-built binaries cannot be directly executed, and the direct compilation of the source code has not been attempted so far. For Windows users we recommend the use of virtualization software packages, such as VirtualBox, which allow Linux as a guest operating system. In that environment Mrcc can be installed in the normal way as described in the previous subsections.
Once you have successfully installed Mrcc, you may wish to test the correctness of the installation. For that purpose, numerous test jobs are at your disposal. The corresponding input files can be found in the MTEST directory created at the installation, where a test script, mtest, is also available. Your only task is to change to the MTEST directory and execute the mtest script. (Please do not forget to add the directory where the Mrcc executables are located to your PATH environment variable.) The test jobs will be automatically executed and you will receive feedback about the results of the tests. The corresponding output files will be left in the MTEST directory, and you can also check them. If all the tests complete successfully, your installation is correct with high probability.
The execution of the test jobs will take for a couple of hours. If you want to run the test on another machine, e.g., on a node of a cluster, you should copy the entire MTEST directory to that machine and start the mtest script there.
Please note that there are some test jobs that allocate a small amount of memory to test the out-of-core algorithms of the program (MINP_*smallmem). These test jobs run with four threads by default when testing OpenMP-parallel executables (i.e., the build.mrcc script was run with the -pOMP switch) with the mtest script. If you run these test jobs with more than four threads, some of them will fail since the memory requirement for OpenMP-parallel runs grows with the number of threads. In this case, the failure of these tests does not indicate a problem with your installation.
Testing of the MPI-parallel executables is also supported. The mtest script automatically detects MPI-parallel executables in PATH if the build.mrcc script was run with the -pMPI switch and runs the appropriate MPI-parallel test jobs (MINP_*_MPI).
To test other external libraries, two switches can be used. With the -l switch, Libxc test jobs (containing _Libxc in the name) are also run, while with the -p switch, the PCMSolver library is also tested (_PCM in the file name). The two options can also be used at the same time, i.e., if mtest is invoked as mtest -l -p, all the test jobs will be executed.
Testing of the features available with external semi-empirical and MM programs, such as ONIOM models including semi-empirical layers or RISM solvation calculations, is also automated. The mtest script automatically searches for those programs, and the corresponding test jobs are executed if they are available.
Please note that you can also create your own test jobs, e.g., if you modify the code, compile the program with new compiler versions, or use unusual combination of keywords. To that end you should calculate a reliable energy for your test job (e.g., using a stable compiler version), include the test keyword and the calculated energy to the MINP file (see the description of the keyword for more details), and copy the MINP file to the MTEST directory renaming it as MINP_$$job_name$\mathrm{>}$. Then the new job will be automatically executed when the mtest script is invoked next time.
Please be sure that the directory where the Mrcc executables are located are included in your PATH environment variable. Note that the package includes several executables, and all of them must be copied to the aforementioned directory, not only the driver program dmrcc. Please also check your LD_LIBRARY_PATH environment variable, which must include the directories containing the libraries linked with the program. This variable is usually set before the installation, but you should not change by removing the names of the corresponding directories. Please do not forget to copy the input file MINP (see Sect. 11) to the directory where the program is invoked.
To run Mrcc in serial the user must invoke the driver of the
package by simply typing
dmrcc
on a Unix console. To redirect the input one should execute dmrcc
as
dmrcc $\mathrm{>}$ out
where out is the output file.
Several executables of the package can be run in OpenMP parallel mode, hence it is recommended to use this option on multiprocessor machines.
The pre-built binaries available at the Mrcc homepage support OpenMP-parallel execution. If you prefer source-code installation, to compile the program for OpenMP parallel execution you need to invoke the build.mrcc script with the -pOMP option at compilation (see Sect. 7). The OpenMP parallelization has been tested with the Intel compiler. Please be careful with other compilers, run, e.g., our test suite (see Sect. 8) with the OpenMP-complied executables before production calculations.
To run the code with OpenMP
you need to set the environment variable OMP_NUM_THREADS
to the number of cores you want to use. E.g., under Bourne shell (bash):
export OMP_NUM_THREADS=4
Then the program should be executed as described above.
The provided binaries are linked with threaded Intel MKL routines, thus,
when those are executed, the environment variable MKL_NUM_THREADS
should also be set, e.g.:
export MKL_NUM_THREADS=4
If source-code installation is preferred, it is recommended to link the
Mrcc objects with threaded BLAS and LAPACK libraries and
employ the required runtime settings of the employed libraries (e.g., define
MKL_NUM_THREADS for Intel MKL).
The binding of threads to hardware elements might affect the performance on certain systems. The thread affinity can be specified with the OpenMP environment variables OMP_PLACES and OMP_PROC_BIND or the Intel MKL specific variable KMP_AFFINITY when the precompiled binaries are executed or Intel MKL is used. For the computations involving the ccsd program, such as DF-, FNO-, and LNO-CCSD computations, nested parallelism is utilized. It is suggested to try setting the OMP_PLACES=cores and OMP_PROC_BIND=spread,close for improved performance. However, in some cases, such thread binding may negatively effect the performance (e.g., when the number of processes exceeds the number of physical CPUs).
Currently executables scf, mrcc, and ccsd can be
run in parallel using the MPI technology.
To compile the program for MPI-parallel execution, you need to invoke the
build.mrcc script with the -pMPI option at compilations
(see Sect. 7). It has been tested with the Intel
compiler and the Open MPI (version 4) and Intel MPI (2017 or later) environments.
If the precompiled binaries are used or Intel MPI 2021 or newer is linked to Mrcc,
it is strongly recommended to install and use the newest stable version of Libfabric
(1.9.0 or later) as some of the previous versions provided via the Intel MPI package
could cause irregular runtime behavior.
The Libfabric library can be downloaded from https://github.com/ofiwg/libfabric.
If not the Intel provided Libfabric library is used, the environment for Intel MPI
should be set using the -ofi_internal=0 option of mpivars.sh (e.g., if
Intel MPI 2019 is installed, source
$$Intel MPI install dir$\mathrm{>}$/parallel_studio_xe_2019/
compilers_and_libraries_2019/linux/mpi/intel64/bin/mpivars.sh
release_mt -ofi_internal=0.
For the MPI-parallel execution, the mpitasks keyword has to be set. Then, it is sufficient to execute the dmrcc binary as usual. The program will spawn the number of scf, mrcc, or ccsd processes specified with the mpitasks keyword and copy the necessary input files to the compute nodes, therefore the input files need to be present only in the directory where dmrcc is executed. Note that the working directory can be the same for all MPI processes, e.g., a directory in the network file system of a computer cluster. Alternatively, process-specific working directories are also supported to exploit local hard drives within compute nodes. In both cases the spawned MPI process will start the execution in the directory with the same path, which might be on a separate file system.
If you wish to run Mrcc with other mpirun options, the MPI-parallel dmrcc_mpi executable should be launched as mpirun -np 1 <options> dmrcc_mpi. You should not run dmrcc using mpirun since it will result in launching mpirun recursively, and your job might fail to start. Please note that the total number of processes will be higher than mpitasks, so you might need to oversubscribe nodes using the appropriate mpirun or scheduler option (e.g., sbatch --overcommit … with SLURM or mpirun --oversubscribe … with Open MPI). For optimal performance, please set mpitasks at the total number of available CPUs, non-uniform memory access (NUMA) nodes, nodes, cores, etc., as the additional number of processes on top of mpitasks are driver processes running mostly in the background and do not require dedicated resources.
On systems consisting of more than one NUMA node (e.g., containing more than one CPU), the performance may be increased by running one process on each NUMA node of the compute nodes. This strategy is beneficial, for instance, when the number of OpenMP processes would otherwise surpass a few tens. Instead, the number of MPI tasks can be increased for better parallel efficiency. Note, however, that in this case the total memory requirement is increased because each process allocates the amount of memory specified in the input file as all MPI algorithms currently available in Mrcc rely on replicated memory strategies.
Pinning processes to CPU cores in MPI parallel runs might affect the performance. When Open MPI is linked to Mrcc, binding can be set by the -bind-to option of mpirun, via modular component architecture (MCA) parameters (e.g., --mca hwloc_base_binding_policy core), or setting the environment variable OMPI_MCA_hwloc_base_binding_policy. It is also suggested to set the Open MPI MCA parameter rmaps_base_inherit to $1$. In the case Mrcc is linked with Intel MPI or the precompiled binaries are used, pinning can be controlled by the I_MPI_PIN and I_MPI_PIN_PROCESSOR_LIST environment variables. If the internode connection is established via an InfiniBand network, other MCA parameters might need to be set as well (e.g., btl_openib_allow_ib to true).
To run Mrcc in parallel using MPI, the environment needs to be set up correctly. You need a working MPI installation with the Mrcc executables available on all nodes used for the job. MPI can be set up by setting the appropriate environment variables (e.g., PATH, LD_LIBRARY_PATH), sourcing the setup script that comes with MPI (e.g., if Intel MPI is used) or maybe in other ways (e.g., modules).
MPI might need to be compiled with your chosen resource manager’s integration if you plan to use its own tools to launch MPI processes. If you do not use a resource manager, you should be able to ssh/rsh between the compute nodes. In this case, you might need to set the proper startup method with the appropriate option (such as ssh, slurm, pbs, ...) of your MPI implementation (e.g., the I_MPI_HYDRA_BOOTSTRAP environment variable or the -bootstrap option of mpirun in the case of Intel MPI). In the case, you use Intel MPI with the SLURM resource manager, you should set the I_MPI_PMI_LIBRARY environment variable pointing to SLURM’s own PMI library. Alternatively, it is suggested that you do not apply a third party PMI library (unset I_MPI_PMI_LIBRARY) and use the mpirun command with ssh bootstrap.
Depending on the hardware and software setup of the network on your cluster, you may need to adjust the fabric interface provider. For instance, for Intel MPI, use the FI_PROVIDER variable, e.g., with verbs, tcp, etc.
For further details, please refer to the documentation of your chosen MPI implementation. You may also find useful tips to solve MPI-related issues on the Mrcc user forum (thread no. 1031, 1032, 1092, 1094). In the case of runtime issues (e.g., the program hangs in an MPI parallel job step or runs with error), it may be helpful to increase the verbosity level, e.g., with Intel MPI via setting the environment variables as I_MPI_DEBUG=5, I_MPI_HYDRA_DEBUG=1, and I_MPI_OFI_PROVIDER_DUMP=1.
In order to test if the problem is related to the MPI settings, not optimal but usually
accessible options can be tried. For instance, for Intel MPI you may try to set one or more of the
following:
unset I_MPI_PMI_LIBRARY
export I_MPI_HYDRA_BOOTSTRAP=ssh
export I_MPI_OFI_PROVIDER=tcp
and depending on the Intel MPI version:
source
$$Intel MPI install dir$\mathrm{>}$/parallel_studio_xe_2019/
compilers_and_libraries_2019/linux/mpi/intel64/bin/mpivars.sh
release_mt -ofi_internal=1
or
source $$Intel MPI install dir$\mathrm{>}$/oneapi/mpi/latest/env/vars.sh
-i_mpi_library_kind=release_mt -i_mpi_ofi_internal=1.
In this section, we discuss the major characteristics of the programs of the Mrcc package, and also provide some information about their use and the corresponding outputs.
Driver for the program system. It calls the programs of the suite (except build.mrcc). It is recommended to run always dmrcc, but advanced users may run the programs one-by-one (e.g., for the purpose of debugging). See also Sect. 9 for further details.
Input reader and analyzer. This program reads the input file MINP, checks keywords, options, and dependencies; sets default values for keywords.
An open-ended atomic orbitals integral code. This code reads and analyzes the molecular geometry, reads the basis sets, and calculates one- and two-electron integrals as well as property integrals over Gaussian-type atomic orbitals. Both the Obara–Saika and the Rys quadrature schemes are implemented for the evaluation of two-electron integrals. In principle integrals over basis functions of arbitrary high angular momentum can be evaluated using the Obara–Saika algorithm.
Hartree–Fock and Kohn–Sham SCF code. It solves the RHF, UHF, ROHF, RKS, UKS, or ROKS equations using either conventional or direct SCF techniques. It also performs the semi-canonicalization of orbitals (if requested) for ROHF wave functions.
Orbital localization program. It performs the localization of MOs using the Cholesky, Boys, or generalized Boys procedures. It also constructs the domains for local correlation calculations.
An efficient three-index integral transformation, density-fitting MP2, RPA, dRPA, SOSEX, and RPAX2 code. The dRPA method is implemented using the modified algorithm of Ref. 54, which scales as the fourth power of the system size, see Ref. 71. drpa is also the driver for local correlation calculations.
Domain construction for local correlation calculations. It assigns the localized MOs (LMOs) to atoms using the Boughton–Pulay method, and for each occupied LMO it constructs a domain of occupied and virtual LMOs on the basis of their spatial distance. Projected atomic orbitals (PAOs) are also constructed if requested.
Integral transformation and orbital optimization code. This program performs the four-index integral transformations of AO integral for correlation calculations. It also carries out the construction of optimized virtual orbitals (OVOs) or MP2 natural orbitals in the case of reduced-cost CC calculations.
A very fast, hand-coded, MO-integral-based (DF) CCSD and CCSD(T) code. The code has been designed for local CC calculations but is also very efficient for conventional CC calculations. Currently it only functions for closed-shell systems, and the spatial symmetry is not utilized.
A fast, hand-coded, MO-integral-based open-shell (DF) CCSD and CCSD(T) code. The code has been designed for local CC calculations but is also efficient for conventional CC calculations. Currently, the spatial symmetry is not utilized.
An efficient explicitly correlated MP2 (MP2-F12) program. It also computes the intermediates for explicitly correlated CCSD and CCSD(T) calculations.
A very fast, hand-coded, integral direct DF CIS, TDA, TD-DFT, ADC(2), CIS(D), CIS(D${}_{\mathrm{\infty}}$), and CC2 code. Currently it only works for closed-shell systems, and the spatial symmetry is not utilized.
This program solves the CPHF/KS equations, constructs relaxed density matrices, calculates first-order molecular properties and Cartesian gradients.
Interface program for QM/MM and embedding calculations.
This program generates the formulas for mrcc. The program also estimates the memory requirement of the calculation. This is a very crude (the symmetry and spin is not treated exactly) but quick estimate. The real memory requirement, which is usually much smaller, is calculated by xmrcc after the termination of goldstone.
It calculates the exact memory requirement for mrcc. Note that it may take a couple of minutes for complicated wave functions (e.g., MRCC derivatives). It prints out five numbers at the end (in MBytes):
Minimal and optimal memory for double-precision (real*8) arrays.
Memory allocated by mrcc for integer arrays.
(= Real*8 + Integer) The minimal and the optimal amount of total required memory. It is not worth starting the calculation if the real physical memory of the machine is smaller than the Minimal value. The performance of the program is optimal if it can use at least as much memory as the Optimal value. If the memory is between the Minimal and Optimal values, out-of-core algorithms will be executed for particular tasks, and it may result in slow down of the code. Please note that the memory available to the program can be specified by keyword mem (see page 12).
Automated, string-based many-body code. It performs the single-point energy as well as derivative calculations for general CC, LR-CC, and CI methods. Abelian spatial symmetry is utilized and a partial spin adaptation is also available for closed-shell systems.
Installation script of the suite. See Sect. 7 for a detailed description.
The input file of the Mrcc package is the MINP file. This file must be placed in the directory where the program is invoked. In addition, if you use your own basis sets (see keyword basis), angular integration grids for DFT calculations (see keyword agrid), or Laplace-quadrature for Laplace transform calculations (see keyword dendec), you may also need the GENBAS file, and then it must be also copied to the above directory.
In general, the execution of Mrcc is controlled by keywords.
The list of the keywords is presented in Sect. 12.
The keywords and the corresponding options must be given
in the MINP file as
…
$$keyword$\mathrm{>}$[=$$option$\mathrm{>}$]
…
You can add only one keyword per line, but there are keywords which
require multiple-line input, and the corresponding variables must be
specified in the subsequent lines as
…
$$keyword$\mathrm{>}$[=$$option$\mathrm{>}$]
$$input record 1$\mathrm{>}$
$$input record 2$\mathrm{>}$
…
$$input record n$\mathrm{>}$
…
The input is not case-sensitive.
Any number of lines can be left blank between two items, however, if
a keyword requires multiple-line input, the lines including the keyword
and its input records cannot be separated. Under similar conditions
any line can be used for comments, but the beginning of a comment line
must not be identical to a keyword because that line may be identified
as a keyword by the input reader and misinterpreted. Thus it is
recommended to start comment lines with some special character, e.g.,
hash mark.
Please note that you can find input files for numerous test jobs in the MTEST directory created at the installation of Mrcc (see Sect. 8). The input files have self-explanatory names and also include a short description at the beginning. You should look at these files for examples for the structure of the input file and the use of various keywords. You can use these files as templates, but please note that these files have been created to thoroughly check the correctness of the code and the installation, and thus some of them contain very tight convergence thresholds as well as unusual combination of (auxiliary) basis sets. In production calculations, you should use the default convergence thresholds (i.e., delete the lines including keywords itol, scftol, cctol, etc.), select the basis set carefully (i.e., set the appropriate option for keyword basis), and use the default auxiliary basis sets (i.e., delete the lines including keywords dfbasis_scf or dfbasis_cor). Please also do not forget to remove keyword test and to specify the amount of memory available to the program by setting the mem keyword.
In this section, the keywords of the Mrcc input file are listed in alphabetical order. To help the reader, the list of the keywords can also be found at the end of the manual.
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. Note also that this keyword only sets the active orbitals for the post-SCF calculation, the MCSCF active orbitals can be specified by keyword mact.
All orbitals are inactive (i.e., single-reference calculation).
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.
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.
active=none
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
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
Specifies the angular integration grid for DFT calculations. 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 [6] 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. See also the description of keywords grid, rgrid, and grtol.
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.
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 are MMMM and NNNN, respectively. The angular grids are selected for each radial point so that the error in the angular integrals will not be larger than ${10}^{1-\mathrm{\U0001d690\U0001d69b\U0001d69d\U0001d698\U0001d695}}$. If grid=TA*, the above procedure is not carried out, but the angular grids are taken from Ref. AhlrichsGrid.
agrid=LD0006-LD0590
for a 974-point Lebedev grid set agrid=LD0974
to use an adaptive grid with at least 110 and at most 974 angular points set agrid=LD0110-LD0974
for a very fine grid use
agrid=LD0110-LD0974
grtol=12
Specifies the angular integration grid for COSX calculations. For more information about the usage and available options, see keyword agrid. To learn more about COSX calculations, see keywords pssp and dual_df.
Same as for agrid
agrid_pssp=LD0006-LD0194
for a fine grid use agrid=LD0006-LD0590
Specifies the small grid for dual grid COSX calculations. For more information about COSX calculations, see keywords pssp and dual_df.
Same as for agrid and agrid_pssp
agrid_pssp_sm=LD0006-LD0074
for a fine grid use agrid=LD0006-LD0194
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.
If the same basis set is used for all atoms, the label of the basis set must be given.
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$\mathrm{>}$:$$basis set$\mathrm{>}$ .
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$\mathrm{>}$
$$basis set label 1$\mathrm{>}$
$$,$$,…,$$-$$,…
$$basis set label 2$\mathrm{>}$
$$,$$,…,$$-$$,…
…
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.
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).
It is the same as embed, but the partitioning defined by keyword corembed will be used.
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.
By default the following basis sets are available for elements H to Kr in Mrcc:
the def2 Gaussian basis sets of Weigend and Ahlrichs [152]: def2-SV(P), def2-SVP, def2-TZVP, def2-TZVPP, def2-QZVP, def2-QZVPP
the augmented def2 Gaussian basis sets of Rappoport and Furche [130]: def2-SVPD, def2-TZVPD, def2-TZVPPD, def2-QZVPD, def2-QZVPPD
F12 basis sets for explicitly correlated wave functions developed by Peterson et al. [120]: cc-pV$X$Z-F12 ($X$ = D, T, Q)
the Gaussian basis sets of Dunning and Hay (LANL2DZ) [19]
the auxiliary basis sets of Hellweg and Rappoport for the augmented def2 Gaussian basis sets [53]: 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 [155]: cc-pV$X$Z-RI-JK, aug-cc-pV$X$Z-RI-JK ($X$ = D, T, Q, 5), def2-QZVPP-RI-JK
Weigend’s universal Coulomb fitting auxiliary basis sets for density fitting/resolution of the identity DFT calculations [Weigend_JFIT]: def2-QZVPP-RI-J
the complementary auxiliary basis sets (CABSs) for explicitly correlated calculations developed by Yousaf and Peterson [OPTRI1, OPTRI2]: cc-pV$X$Z-F12-OPTRI ($X$ = D, T, Q), aug-cc-pV$X$Z-OPTRI ($X$ = D, T, Q, 5)
From Na to Og the following basis sets are available, which must be used together with the corresponding ECP (see also the description of keyword ECP):
the def2 Gaussian basis sets of Weigend and Ahlrichs [152]: def2-SV(P), def2-SVP, def2-TZVP, def2-TZVPP, def2-QZVP, def2-QZVPP
the augmented def2 Gaussian basis sets of Rappoport and Furche [130]: def2-SVPD, def2-TZVPD, def2-TZVPPD, def2-QZVPD, def2-QZVPPD
the auxiliary basis sets of Hellweg and Rappoport for the augmented def2 Gaussian basis sets [53]: def2-SVPD-RI, def2-TZVPD-RI, def2-TZVPPD-RI, def2-QZVPD-RI, def2-QZVPPD-RI
the auxiliary basis sets of Hättig and Hill for correlation calculations with the core-valence and PP basis sets [cc-pwCVXZ-RI, aug-cc-pwCVXZ-PP-RI]: aug-cc-pwCV$X$Z-RI, cc-pV$X$Z-PP-RI and aug-cc-pV$X$Z-PP-RI, aug-cc-pwCV$X$Z-PP-RI ($X$ = D, T, Q, 5)
Please note that some of the above basis sets are not available for all elements.
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.
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.
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.
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$\mathrm{X}$Z. Then the diffuse functions will
be automatically removed from the hydrogen atoms.
Only the conventional AO basis set can be specified with this keyword. For the auxiliary basis sets used in density-fitting or F12 approximations see the description of keywords dfbasis_*.
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.
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.
If the (aug-)cc-pV$X$Z-PP basis set does not exist for an element with $Z\le 28$, the program will automatically attempt to use the corresponding (aug-)cc-pV$X$Z basis instead.
If the (aug-)cc-pV($X$+d)Z basis set does not exist for an element (i.e., $Z\le 12$ or $Z\ge 19$), the program will automatically attempt to use the corresponding (aug-)cc-pV$X$Z basis instead.
none, that is, the basis set must be specified (excepting the case when Mrcc is used together with another code, that is, iface $\ne $ none).
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
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
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
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
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
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.
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.
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
Specifies the small basis set used in dual basis-set calculations as well as for generating SCF initial guess (scfiguess=small).
the options are the same as for keyword basis, but there is an additional one, none, which means that no small basis is defined.
basis_sm=none
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
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
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.
on, off, or $$any real number$\mathrm{>}$
basopt=off
To optimize a basis set variationally set basopt=on
To optimize a basis set minimizing the difference of the calculated energy and -76.287041 E${}_{h}$ set basopt=-76.287041
Specifies the bond function (BF) basis (see Ref. 89 for details).
No BFs are used.
name of the BF basis to be used.
The format of the name of the BF basis, $$BF basis name$\mathrm{>}$, is $$AO basis name$\mathrm{>}$-$$BF type$\mathrm{>}$. 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.
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$\mathrm{>}$:$$BF basis name$\mathrm{>}$. E.g., BCH:6-31G-1s1p is 6-31G-1s1p BF basis optimized for the C–H bond.
If BF bases are used, the geometry must be given in mol format (see the description of keyword geom)
bfbasis=none
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
Specifies the number of gradient vectors used for the BFGS update in quadratic SCF calculations using the BFGS algorithm (qscf=BFGS).
$$any positive integer$\mathrm{>}$
bfgsmem=10
to increase the number of vectors to 15 set bfgsmem=15
Threshold (in E${}_{h}$) for starting the BFGS algorithm in quadratic SCF calculations using the BFGS algorithm (qscf=BFGS). The calculation starts with the conventional DIIS-based convergence acceleration, and the BFGS algorithm will be switched on if the maximum norm of the gradient is smaller than bfgstol.
$$any positive real number$\mathrm{>}$
bfgstol=1e-3
for a convergence threshold of ${10}^{-4}$ E${}_{h}$ set bfgstol=1e-4
Specifies the optimization algorithm used for the Boys–Foster localization.
Uses standard Jacobi rotations.
Uses trust-region optimization technique.
boysalg=Jacobi
for the trust-region optimization set boysalg=Newton
Boughton–Pulay completeness criterion [9] 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.
This number will be used as the completeness criterion.
bpcompo=0.95 for local excited-state calculations, bpcompo=0.985 otherwise
Atom domains determined by bpcompo are also utilized to construct local fitting domains in the case of localcc$\mathrm{\ge}$2016 according to Ref. 107.
to set a threshold of 0.99 type bpcompo=0.99
Boughton–Pulay completeness criterion [9] for virtual orbitals (projected atomic orbitals). See also keyword bpcompo.
This number will be used as the completeness criterion.
bpcompv=0.98
to set a threshold of 0.95 type bpcompv=0.95
Boughton–Pulay completeness criterion [9] for the occupied orbitals of an extended domain. See also keyword bpcompo.
This number will be used as the completeness criterion.
bpedo=bpcompo is set if bpedo is not specified and not employed in the local correlation calculation
to set a threshold of 0.9998 type bpedo=0.9998
Boughton–Pulay completeness criterion [9] for the virtual orbitals (projected atomic orbitals) of an extended domain. See also keyword bpcompo.
This number will be used as the completeness criterion.
bpedv=0.995 is set as default in the case of localcc$\mathrm{\ge}$2016 according to Ref. 107.
bpedv=bpcompv is set if bpedv is not specified and not employed in the local correlation calculation
to set a threshold of 0.99 type bpedv=0.99
Boughton–Pulay completeness criterion [9] for the occupied orbital of a primary domain. See also keyword bpcompo.
This number will be used as the completeness criterion.
bppdo=0.999 is set as default in the case of localcc$\mathrm{\ge}$2016 according to Ref. 107.
bppdo=bpcompo is set if bppdo is not specified and not employed in the local correlation calculation
to set a threshold of 0.99 type bppdo=0.99
Boughton–Pulay completeness criterion [9] for virtual orbitals (projected atomic orbitals) of a primary domain. See also keywords bppdo and bpcompo.
This number will be used as the completeness criterion.
bppdv=bpcompv
to set a threshold of 0.99 type bppdv=0.99
Specifies the type of the calculation.
Hartree–Fock or Kohn–Sham SCF calculation, the type of the calculation can be controlled by keyword scftype (see also keyword scftype).
Restricted, unrestricted, restricted open-shell Hartree–Fock SCF; restricted, unrestricted, restricted open-shell Kohn–Sham SCF; or multi-configurational SCF calculation, respectively. The type of the HF/KS/MCSCF procedure 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.
Kohn–Sham SCF calculation with the specified density functional. The type of the Kohn–Sham procedure (i.e., RKS, UKS, or ROKS) can be controlled by keyword scftype (see also keyword scftype). The options are identical to those of keyword dft (except for off, user, userd, and rsdh), 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.
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.
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. For excited-state calculations with double hybrid functionals see also keyword dhexc.
TD-DFT in the Tamm–Dancoff approximation (TDA). For HF reference it is equivalent to CIS. It is only available with density fitting. For excited-state calculations with double hybrid functionals see also keyword dhexc.
Second-order Møller–Plesset (MP2) calculation. Note that efficient MP2 calculations are only possible with the density-fitting (resolution-of-identity) approximation, which can be performed with the DF-MP2 or the RI-MP2 keywords or by setting dfbasis_cor$\mathrm{\ne}$none. In the case of a DF-MP2 calculation, the spin-component scaled MP2 (SCS-MP2) [39] and the scaled opposite-spin MP2 (SOS-MP2) [59] energy will also be computed (see also keywords scsps and scspt).
Third-order Møller–Plesset (MP3) calculation [41]. Note that MP3 calculations without the density-fitting (resolution-of-identity) approximation are not possible with RHF reference wave functions.
Scaled opposite-spin second-order Møller–Plesset (SOS-MP2) calculation [59] 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).
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.
Direct random-phase approximation (dRPA) calculation (see Eqs. 7 and 8 in Ref. 148). 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).
Random-phase approximation (RPA) calculation (see Eqs. 10 and 13 in Ref. 148, 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).
Second-order screened exchange (SOSEX) [40] calculation (see Eqs. 7 and 9 in Ref. 148), the dRPA and the rPT2 [132] energies are 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).
The scaled-equation and down-scaled dRPA and SOSEX methods of Ref. 100. See also the notes for options dRPA and SOSEX as well as the descriptions of keywords scspe and scsph.
RPAX2 calculation (see Eqs. 17 to 19 in Ref. 54). 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).
Configuration interaction singles (CIS) calculation [91]. 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.
Configuration interaction singles with perturbative correction for double excitations [CIS(D)] calculation [46, 90]. Note that CIS(D) calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-CIS(D) [$\equiv $ RI-CIS(D)] calculation is performed [that is, options CIS(D), DF-CIS(D), and RI-CIS(D) are synonyms].
Iterative doubles correction to configuration interaction singles [CIS(D${}_{\mathrm{\infty}}$)] calculation [46, 93, 94]. Note that CIS(D${}_{\mathrm{\infty}}$) calculations are only possible with the density-fitting (resolution-of-identity) approximation, and, by default, a DF-CIS(D${}_{\mathrm{\infty}}$) [$\equiv $ RI-CIS(D${}_{\mathrm{\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, 43, 93, 94]. 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].
Second-order coupled-cluster singles and doubles (CC2) calculation [10, 43, 93, 94]. 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.
Scaled opposite-spin CC2, CIS(D), CIS(D${}_{\mathrm{\infty}}$), and ADC(2) [SOS-CC2, SOS-CIS(D), SOS-CIS(D${}_{\mathrm{\infty}}$), SOS-ADC(2)] calculation [10, 156, 93, 94, 90]. An ${N}^{4}$-scaling algorithm based on the Cholesky decomposition/Laplace transform of energy denominators is executed. It is only available with density fitting. Note that two different parameterizations are implemented for SOS-CIS(D) (see Refs. 52 and 133). See also keywords scsps, scspt, scsps_t, and scspt_t.
Spin-component scaled CC2, CIS(D), CIS(D${}_{\mathrm{\infty}}$), and ADC(2) [SCS-CC2, SCS-CIS(D), SCS-CIS(D${}_{\mathrm{\infty}}$), SCS-ADC(2)] calculation [10, 156, 93, 94, 90]. It is only available with density fitting. Note that two different parameterizations are implemented for SCS-CIS(D) (see Refs. 52 and 133). See also keywords scsps, scspt, scsps_t, and scspt_t.
The corresponding single-reference CC calculation with perturbative corrections (see Ref. 63).
The corresponding single-reference CC calculation with perturbative corrections (see Ref. 63).
The corresponding CCSD(T)${}_{\mathrm{\Lambda}}$, CCSDT(Q)${}_{\mathrm{\Lambda}}$, etc. calculation (see Ref. 63).
The corresponding iterative approximate single-reference CC calculation (see Ref. 63).
The corresponding iterative approximate single-reference CC calculation (see Ref. 63).
The corresponding iterative approximate single-reference CC calculation (see Ref. 63).
The corresponding iterative approximate single-reference CC calculation (see Ref. 63).
The corresponding single-reference CC calculation with perturbative corrections using ansatz A (see Ref. 65).
The corresponding single-reference CC calculation with perturbative corrections using ansatz B (see Ref. 65).
The corresponding single-reference CC calculation with perturbative corrections using ansatz A (see Ref. 65).
The corresponding single-reference CC calculation with perturbative corrections using ansatz B (see Ref. 65).
The explicitly correlated MP2 (MP2-F12) method using the 2B ansatz, the F + K commutator approximation, fixed amplitudes, spin-flipped geminals in open-shell calculations [MP2-R12_Noga, MP2-F12_TM, F12pT]. It is only available with DF.
The explicitly correlated CCSD and CCSD(T) utilizing the CCSD(F12*) [CCSD(F12)] and (T+) [F12pT] models and the above MP2-F12 approximation. It is only available with DF.
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$\mathrm{>}$) 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).
For excited-state calculations with the TD-HF, TDA, TD-DFT (including double hybrid approaches), CIS, CIS(D), CIS(D${}_{\mathrm{\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$\mathrm{\ne}$off) calculations will be carried out for the excited states.
The active orbitals can be selected and the MRCI/CC calculations can be controlled by keywords nacto, nactv, active, maxex, and maxact. Note that, by default, MRCI/CC calculations are executed using HF reference. To use MCSCF orbitals scftype must be set to MCSCF, and the MCSCF wave function must be defined by keywords docc and mact. In this case nacto and nactv are taken over from the MCSCF calculation, and a MRCI/CC calculation will be run automatically. You should only set the above keywords if you want to run the post-MCSCF calculation with an active space different from the MCSCF one.
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.
Local correlation methods, for both ground and excited states, can be run if the prefix “L” is added to the corresponding option of the keyword, e.g., as LMP2, LdRPA, LCCSD(T), LADC(2), LCC2, 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., and for ADC(2) [LNO-ADC(2)]. Both options are equivalent to setting localcc=on.
For the dRPA, RPA, and SOSEX methods the use of PBE orbitals is recommended (this is set by default).
For the RPAX2 method the use of PBEx orbitals is recommended (this is set by default).
If FNO-CCSD or FNO-CCSD(T) is set, a CCSD or CCSD(T) calculation is performed utilizing the FNO and NAF approximations [FNONAF] and the PPL+ and (T+) corrections [FNOExt]. In these cases, ovirt=ppl, naf_cor=on, and dfbasis_cor=auto are set by default. CCSD-PPL+ and CCSD-PPL+(T+) are equivalent to FNO-CCSD and FNO-CCSD(T).
If the FNO- prefix is added to the acronym of an explicitly correlated CC method, that is, if FNO-CCSD-F12, FNO-CCSD(T)-F12, FNO-CCSD(F12*)(T+), etc. is set, the corresponding calculation is performed utilizing the FNO, NAF, and NAB approximations [RedCostF12]. In these cases, ovirt=mp2, naf_cor=on, naf_f12=on, nab=on, and dfbasis_cor=auto are set by default.
calc=SCF
To run a CCSD(T) calculation, the user should set calc=CCSD(T)
For DF-HF (RI-HF) calculations type:
calc=DF-HF
which is equivalent to the following input:
calc=SCF
dfbasis_scf=auto
For a local CCSD(T) calculation using the local natural orbital approximation set calc=LCCSD(T) or calc=LNO-CCSD(T)
For a RI-MP2 calculation set calc=MP2
For a DFT calculation with the B3LYP functional set calc=B3LYP
Direct RPA calculation with Kohn–Sham orbitals calculated with the
PBE functional:
calc=dRPA
dft=PBE
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
Maximum number of iteration steps in correlated calculations (CC, CI, RPA, …).
$$any positive integer$\mathrm{>}$
ccmaxit=50
to increase the maximum number of CC iterations to 100 give ccmaxit=100
Specifies the CC program to be used.
The automated, string-based CC program mrcc will be called.
The very fast, hand-coded CCSD(T) codes, ccsd or uccsd, will be executed (currently the spatial symmetry cannot be utilized).
The very fast, hand-coded, integral direct DF-CIS code cis will be executed (currently the spatial symmetry cannot be utilized).
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.
ccprog=ccsd for conventional and explicitly correlated CCSD and CCSD(T) calculations, ccprog=cis for CIS, CIS(D), CIS(D${}_{\mathrm{\infty}}$), ADC(2), and CC2 calculations, ccprog=mrcc otherwise.
to use the mrcc code for CCSD or CCSD(T) calculations give ccprog=mrcc
Specifies the CCSD algorithm employed in the ccsd or uccsd program if ccprog=ccsd.
An MO-based, hand-coded, OpenMP-parallel CCSD algorithm is invoked [138]. All the necessary four-center integrals (including the ones with four virtual orbitals) are stored on disk after the integral transformation (and integral assembly) steps. The four-center integrals are read in each CCSD iteration from these files. Careful optimization was performed significantly improving both the CPU usage and the OpenMP parallelization efficiency of the implementation (see Ref. 41 for ccsd).
For ccsd: A highly-optimized, hand-coded, OpenMP-parallel, extremely low-memory, in-core, integral-direct ${t}_{1}$-transformed DF-CCSD algorithm is invoked [41]. The necessary four-center integrals and intermediates are constructed completely in memory in an integral direct manner in each iteration, and disk I/O is avoided. This implementation, including the cost of the repeated integral assembly, is still more efficient than the optimized ccsdalg=disk one, especially if the OpenMP parallelization is employed. For uccsd: currently the most storage-intensive four-external four-center integrals are assembled on-the-fly from the corresponding two- and three-center DF integrals in each CCSD iteration. The remaining four-center integrals are stored on disk.
The dfdirect algorithm is currently only available with the ccsd or uccsd programs, i.e., via ccprog=ccsd.
ccsdalg=dfdirect if density fitting is applied for the correlation energy calculation, e.g., for calc=DF-CCSD, calc=DF-CCSD(T), calc=LNO-CCSD(T), etc. If density fitting is not employed, i.e., dfbasis_cor=none, then ccprog=disk is selected automatically.
ccsdalg=disk switches to the CCSD algorithm that stores the four-center integrals on disk.
Controls the parallel execution of Intel MKL subroutines in program ccsd.
Threaded/parallel version of MKL is used
Sequential MKL procedures are called in parallel
ccsdmkl=seq
to run BLAS subroutines on more than $1$ thread in program ccsd set ccsdmkl=thr
Use this keyword to restart canonical or reduced-cost (i.e., not local) CC calculations from previously calculated amplitudes if ccprog=ccsd and talg=occ. For restarting local correlation calculations, see keyword lccrest. For restarting canonical CI or CC calculations with ccprog=mrcc, see keyword rest.
Turns off checkpointing and performs the CC calculation normally
Turns on checkpointing creating the ccsd.rst and pt.rst files and performs the CC calculation normally
Turns on checkpointing and restarts the CC calculation from the last completed iteration for CCSD or the next unfinished $ijk$ index triplet for the (T) correction
In the case of restarting a CC calculation, program ccsd is run directly skipping the SCF procedure and the integral transformation. Therefore, the VARS, 55, and DFINT_* files, as well as the checkpoint files ccsd.rst and pt.rst from the previous run need to be present in the working directory.
ccsdrest=trf
to restart a CCSD/CCSD(T) calculation set ccsdrest=ccsd
Sets the number of outer OpenMP threads in program ccsd while performing the CCSD iteration.
$$any positive integer$\mathrm{>}$
ccsdthreads=2
to reduce the memory requirement of a CCSD calculation by turning off nested OpenMP set ccsdthreads=1
Convergence threshold for the energy in correlated calculations (CC, CI, dRPA, RPA, etc.). The energy will be accurate to ${10}^{-\mathrm{\U0001d68c\U0001d68c\U0001d69d\U0001d698\U0001d695}}$ E${}_{h}$.
$$any integer$\mathrm{>}$
cctol=8 for property calculations, cctol=[-${\mathrm{log}}_{10}$(optetol)]+2 for geometry optimizations, cctol=5 for localcc$\mathrm{\ge}$2016, cctol=6 otherwise
for an accuracy of ${10}^{-8}$ E${}_{h}$ one must give cctol=8
Charge of the system.
$$any integer$\mathrm{>}$
charge=0
for the Cl${}^{-}$ ion one should give charge=-1
Specifies what type of algorithm is to be used in CIS, TDA, TD-HF, and TD-DFT calculations.
Conventional algorithm, two-electron integrals are stored on disk
Completely I/O-free, integral-direct algorithm, two-electron integrals are recalculated in each iteration step.
Partially I/O-free, integral-direct algorithm; recommended if the I/O is fast and/or few states are required.
Variant of direct2, but usually slower.
Based on the size of the molecule the program will automatically select the most efficient one from the above options.
cialg=auto
to use disk-based algorithm set cialg=disk
The initial guess vectors for CI and LR-CC calculations can be specified using this keyword.
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:
$$
where $$ is the level of excitation, and the electrons are
promoted from occupied orbitals $$ to virtual
orbitals $$ with spins
$$ ($$ is 1 for alpha and 0
for beta), respectively. $$ is the corresponding coefficient.
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.
ciguess=off
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.
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.
The guess vector is not required to be normalized, it is done automatically.
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).
Specifies the computational point group. All calculations will use the specified Abelian group. See Sect. 13 for more details.
The molecular symmetry is automatically recognized.
Schönflies symbol of the Abelian point group such as C1, Ci, Cs, C2, C2v, C2h, D2, D2h
cmpgrp=C1 is equivalent to symm=off
cmpgrp=auto
to use ${C}_{2v}$ point group for benzene set cmpgrp=C2v
Controls the restart procedure for composite calculations. Currently, only ONIOM calculations can be restarted with the keyword.
The ONIOM calculation is restarted and halted after the layer boundaries are defined. This option is designed for the setup of layer boundaries after a successful calculation on the full system.
The composite calculation is restarted after the $n$th subcalculation is completed. This option is designed for the modification of composite calculations after successful subcalculations of a composite method. See the keyword subminp to modify the initial settings.
The composite calculation is restarted from the last successful subcalculation of a composite method. This option is designed, e.g., for power failures.
No restart is requested
The comprest=0 option requires the MOCOEF file for Mrcc, the corresponding auxiliary file (aux) for Mopac, and the wbo file for xTB calculations.
In the case of single point ONIOM calculations with option comprest $\mathrm{>}$ 0, the ONIOMSPEC and iface* files of the corresponding subcalculations are checked if they are available. Additionally, the GRAD* and VARS* files are also searched for if gradient or electronic embedding calculations are requested. Calculations using the Mopac or the xTB programs are automatically checked by the external codes.
Except of comprest=0, the boundary settings are read from the ONIOMSPEC file, while the boundary specification in the MINP file is ignored.
In the case of comprest=auto, the keyword scfiguess is set to restart, and the SCFDENSITIES file is reused if the upcoming calculation uses Mrcc.
comprest=off
If a series of ONIOM calculations are stopped due to power failure at unknown points of the composite procedure, type comprest=auto to continue from the last successful subcalculation.
To reconfigure the boundary settings after the first ONIOM subcalculation, type comprest=0. After the setup is made, type comprest=1 to continue with the second subcalculation.
Specifies whether the core electrons are correlated.
Frozen core approximation
All core electrons are correlated
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.
core=frozen
to correlate all core electrons set core=corr or core=0
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.
Conventional case, a single model defined by calc is used for the entire system.
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.
corembed=off
Local correlation methods available with localcc$\mathrm{\ge}$2015 (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. 48). Note that models with KS-DFT reference, such as dRPA, SOSEX, etc., are not available for multi-level local correlation calculations.
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 2021 (or for previous versions according to lcorthr=loose and localcc=2015 or localcc=2016) are employed for the low-level model of the environment.
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
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
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
Switches on the approximate long-range spin polarization approach of Ref. ROLocalMP2 for the restricted open-shell local correlation methods with localcc$\mathrm{\ge}$2021. This option treats the orbital domains without singly occupied orbitals via a closed-shell formalism with negligible approximations and considerable gain in efficiency.
on or off
csapprox=on
Switches on the core-valence separation (CVS) for excited-state calculations [CVS, CVS2]. Only excitations including at least one core orbital specified by this keyword will be included in the excited-state wave function.
The CVS approximation will not be used.
The serial number of the core orbitals (i.e., absolute serial number according to orbital energy order) must be given separated by commas. Serial numbers separated by dash means that all the orbitals in the specified interval will be included. For unrestricted orbitals, only the serial numbers of the alpha orbitals must be given.
A projected eigenvalue equation [CVS] is used for ccprog=ccsd. It is recommended for higher-order linear-response CC and CI calculations.
A decoupled Hamiltonian [CVS2] is used for ccprog=cis. It is recommended for CIS, TDHF, CIS(D), ADC(2), and various (DH-)TDDFT calculations [CVS_DH]. In these cases, the reduced-cost and reduced-scaling approximations are also available (see the description of keyword redcost_exc and Ref. RedCostCVS).
If the CVS approximation is employed, it is recommended to give an appropriate initial guess using the ciguess keyword.
cvs=off
to consider orbitals 1, 3, 4, and 5 set:
cvs=1,3-5
Diagonal Born–Oppenheimer correction (DBOC) (available only with Cfour).
on or off
dboc=off
for a DBOC calculation set dboc=on
Selects the algorithm for the decomposition of energy denominators, Cholesky-decomposition or Laplace transform, for canonical SOS-MP2, dRPA (also required for SOSEX), SOS-CC2, SOS-CIS(D), SOS-CIS(D${}_{\mathrm{\infty}}$), and SOS-ADC(2) as well as for local MP2 and dRPA calculations. The dRPA calculation is performed using the modified algorithm of Heßelmann [54] 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. 72). 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.
for SOS-CC2, SOS-CIS(D), SOS-CIS(D${}_{\mathrm{\infty}}$), and SOS-ADC(2) calculations the decomposition will not be used, but a fifth-power scaling algorithm will be executed
Cholesky decomposition will be used
Laplace transform will be used
dendec=Laplace for SOS-MP2, SOS-CC2, SOS-CIS(D) SOS-CIS(D${}_{\mathrm{\infty}}$), and SOS-ADC(2); dendec=Cholesky otherwise
The algorithms based on the Laplace-transformed technique use minimax quadratures obtained from Ref. 146.
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. 146). The subsequent NN lines must contain, respectively, the weights and quadrature points.
to use Laplace transform give dendec=Laplace
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. 60, 61, 62, 64, 113, 114 for more details.
Density-matrix calculation (for geometry optimizations, first-order properties, etc.)
Density-matrix first derivatives (for second-order property calculations, available only with Cfour)
Transition density matrices (for transition moment calculations)
Second and third derivatives of the density-matrix (for third-order property calculations, available only with Cfour)
dens=2 for geometry optimizations and QM/MM calculations, dens=0 otherwise
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).
If dens $\ne $ 0, a population analysis is also performed, and Mulliken and Löwdin atomic charges as well as Mayer bond orders are computed.
for the calculation of both one- and two-particle density matrices set dens=2
Specifies whether the density-based basis-set correction proposed by Giner, Toulouse, and co-workers [83] is calculated.
The density-based basis-set correction is not calculated.
Both the CCSD(T) energy and the correction are calculated.
Only the density-based basis-set correction is calculated.
To calculate the correction, even if the CCSD(T) energy is not required, the following
format should be used (see keyword dft):
denscorr=2
calc=DF-CCSD(T)
dft=userd
1
1.00 HFx
1
1.0 ECMD_PBEc
In the case of ROHF calculations, the semicanonical orbitals are used by default (see keywords rohfcore and rohftype). When comparing the results from other programs, it might be important to verify the ROHF orbitals.
denscorr=0
for only the density-based basis-set correction, the user should set denscorr=2
Specifies how the inverse of the two-center Coulomb integral matrix is decomposed in density fitting direct SCF calculations.
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.
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 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.
dfalg=InvSqrt for property calculations, dfalg=LinEq otherwise
to use Cholesky decomposition set dfalg=Cholesky
Specifies the complimentary auxiliary basis set (CABS) for explicitly correlated methods.
No CABS.
The specified basis set is used as CABS. For the specification of the basis the same rules apply as for keyword basis, see the description of keyword basis.
This option can only be used if Dunning’s aug-cc-pV$X$Z or Peterson’s cc-pV$X$Z-F12 basis sets are used as the normal basis set. In this case, if dfbasis_cab= auto, the corresponding “OPTRI” basis set will be used as CABS.
For the available CABSs, see the notes for keyword basis.
dfbasis_cab=auto for all the explicitly correlation methods, dfbasis_cab=none otherwise.
To use the cc-pVTZ-F12-OPTRI CABS with the aug-cc-pVDZ AO basis set:
basis=aug-cc-pVDZ
dfbasis_cab=cc-pVTZ-F12-OPTRI
Specifies whether the density fitting approximation will be used in the correlated calculations and also specifies the fitting basis set.
The density fitting approximation is not used for the correlated calculation.
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.
This option can only be used if Dunning’s (aug-)cc-pV$X$Z or (aug-)cc-pV($X$+d)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) and (aug-)cc-pV($X$+d)Z 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.
For the available fitting basis sets see the notes for keyword basis on page 12.
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.
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 $\ne $ off), dfbasis_cor=none otherwise.
To use the cc-pVTZ-RI fitting basis in the correlated calculation for all atoms the input must include dfbasis_cor= cc-pVTZ-RI
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
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
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
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. Exceptions are the KS-DFT calculations with non-hybrid functionals, where the def2-QZVPPD-RI-J basis set will be used in the following cases: for H with the cc-pVDZ, aug-cc-pVDZ, cc-pVDZ-F12, def2-SV(P), def2-SVP, def2-TZVP, def2-SVPD, def2-TZVPD, and LANL2DZ AO basis sets; for other atoms with the cc-pVDZ, cc-pVTZ, cc-pV(D+d)Z, cc-pV(T+d)Z, aug-cc-pVDZ, aug-cc-pVTZ, aug-cc-pV(D+d)Z, aug-cc-pV(T+d)Z, cc-pVDZ-F12, def2-SV(P), def2-SVP, def2-TZVP, def2-TZVPP, def2-SVPD, def2-TZVPD, def2-TZVPPD, and LANL2DZ AO basis sets.
dfbasis_scf=auto if dfbasis_cor$\mathrm{\ne}$none and for DFT calculations, dfbasis_scf= none otherwise.
Specifies the small fitting basis set for dual auxiliary basis set calculations. Any basis set can be used as a small auxiliary basis set which is available in Mrcc or you can use your own sets using the GENBAS file. For a small set of atoms (H, C, N, O, F, P, S, and Cl), preoptimized cc-pV(T,Q)Z-RI-JK-dual basis sets are also available [DualDF]. These are the recommended ones for calculations with the (aug-)cc-pV(T,Q)Z AO basis sets. See also the description of keyword dual_df for further details and examples.
Same as dfbasis_scf.
dfbasis_scf_sm=dfbasis_scf
For more information on dual auxiliary basis calculations see the keyword dual_df.
Specifies the integral transformation program to be used for the transformation of three-center Coulomb integrals.
the drpa program will be called
the ovirt program will be called
dfintran=ovirt if ovirt$\mathrm{\ne}$off, dfintran=drpa otherwise.
to use the ovirt code set dfintran=ovirt
Use this keyword to perform DFT calculations and to specify the functional.
No DFT calculation is carried out.
The name of the functional, see Table LABEL:FuncTable for the available functionals.
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 [56]).
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 [56]).
the HF exchange, denoted by HFx
the long-range HF exchange, denoted by lrHFx. The range-separation parameter ($\omega $) should be specified in a.u. after the lrHFx flag separated by space, see the example below. If the parameter is not set, $\omega $ $=$ 0.15 a.u. will be used.
the short-range HF exchange, denoted by srHFx. The syntax is same as for lrHFx, see the example below.
the MP2, MP3, dRPA, and SOSEX correlation, denoted, respectively, by MP2, MP3, 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.
the long- and short-range dRPA correlation denoted, respectively, by lrdRPA and srdRPA. The corresponding spin components, lrdRPAs, srdRPAs, lrdRPAt, and srdRPAt, are also available; the syntax is same as for lrHFx.
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$\mathrm{>}$
$$coefficient 1$\mathrm{>}$ $$functional name 1$\mathrm{>}$
$$coefficient 2$\mathrm{>}$ $$functional name 2$\mathrm{>}$
$$coefficient 3$\mathrm{>}$ $$functional name 3$\mathrm{>}$
…
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$\mathrm{>}$
$$coefficient 1$\mathrm{>}$ $$functional name 1$\mathrm{>}$
$$coefficient 2$\mathrm{>}$ $$functional name 2$\mathrm{>}$
$$coefficient 3$\mathrm{>}$ $$functional name 3$\mathrm{>}$
…
$$number of entries for energy$\mathrm{>}$
$$coefficient 1’$\mathrm{>}$ $$functional name 1’$\mathrm{>}$
$$coefficient 2’$\mathrm{>}$ $$functional name 2’$\mathrm{>}$
$$coefficient 3’$\mathrm{>}$ $$functional name 3’$\mathrm{>}$
…
See option user for the possible values of
$$functional name n$\mathrm{>}$ and $$functional name n’$\mathrm{>}$.
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.
Range-separated double hybrid functionals
of Mester and Kállay [RSDH] based on the ansatz of Kalai and
Toulouse [Toulouse_RSDH].
The entries should be specified in the
subsequent lines as follows (see also the examples below):
ssr$$functional name 1$\mathrm{>}$
ssr$$functional name 2$\mathrm{>}$
…
$$coefficient 1$\mathrm{>}$
$$coefficient 2$\mathrm{>}$
Here, $$coefficient 1$\mathrm{>}$ and $$coefficient 2$\mathrm{>}$
stand for the parameters $\lambda $ and $\omega $, respectively.
ssr$$functional name i$\mathrm{>}$ is
the range-separated variant of the corresponding functional evaluated
utilizing the local-scaling approximation
proposed by Scuseria and co-workers [30].
This framework works with arbitrary combination of
LDA/GGA/meta-GGA exchange and correlation functionals using either the built-in ones
or those from the Libxc library.
dft=PBE for dRPA, RPA, SOSEX, and for scaled-equation and down-scaled dRPA and SOSEX; dft=PBEx for RPAX2; otherwise dft=off
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 [36, 37] 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,
CAM-B3LYP-D3, wB97X-D3, M06-L-D3, MN15-D3,
M06-2X-D3, LC-wPBE-D3, PW91-D3, SCAN-D3,
dRPA75-D3, HF-D3.
See also the description of keyword edisp.
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).
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.
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.
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.
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.
For the VV10 van der Waals functional you can modify parameters $b$ and $C$ (see Ref. 150) if it is used with the user or userd options. For that purpose the two parameters should be specified after the VV10 or VV10NL flag separated by spaces, see the example below. If the parameters are not set, those of Ref. 150 will be used.
The range-separated variant of a functional is available with the ssr prefix utilizing the local-scaling approximation proposed by Scuseria and co-workers [30]. This framework works for arbitrary LDA/GGA/meta-GGA exchange and correlation functionals using the either built-in ones or those from the Libxc library.
ECMD functionals proposed by Giner and co-workers [83] are available with the ECMD_ prefix. This framework works for the built-in PBEc functional or arbitrary LDA/GGA/meta-GGA correlation functionals from the Libxc library.
To perform a DFT calculation with the B3LYP functional give dft=B3LYP or calc=B3LYP
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
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
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
SOSEX calculation with Kohn–Sham orbitals calculated with the
LDA exchange functional:
calc=SOSEX
dft=LDA
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
B2PLYP calculation, the MP2 contribution is evaluated using local
MP2 approximation:
calc=LB2PLYP
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
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
Local dRPA calculation with Kohn–Sham orbitals calculated with the
PBE functional:
calc=LdRPA
dft=PBE
To perform a DFT calculation with the B3LYP functional using its Libxc implementation set calc=HYB_GGA_XC_B3LYP5
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
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 VV10 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 VV10 8.0 0.0093
DFT calculation with a user-defined CAM-B3LYP functional. A
range-separation parameter ($\omega $) of 0.33 a.u. is set. If you want
to use the default $\omega $ $=$ 0.15 a.u., simply drop the last number
for the lrHFx entry.
dft=user
3
1.00 HYB_GGA_XC_CAM_B3LYP
0.19 HFx
0.46 lrHFx 0.33
An alternative implementation of CAM-B3LYP:
dft=user
3
1.00 HYB_GGA_XC_CAM_B3LYP
0.65 HFx
-0.46 srHFx 0.33
Alternative definition for the RS-PBE-P86 functional utilizing built-in
functionals:
dft=rsdh
ssrPBEx
ssrP86
0.5
0.7
Alternative definition for the RS-PBE-P86 functional utilizing
functionals from the Libxc library:
dft=rsdh
ssrGGA_X_PBE
ssrGGA_C_P86
0.5
0.7
To define third-order double hybrid functionals containing a
perturbative third-order correction [XYG3MP3, DSDdRPA_DSD3], the
MP3 option can be used. For example, the XYG9[f2]@B44LYP
functional of Ref. XYG3MP3 can be specified as
dft=userd
4
0.56 B88
0.44 HFx
0.19 VWN5
0.81 LYP
7
0.0881 LDA
-0.0644 B88
0.9692 HFx
0.2730 VWN5
0.6447 MP2s
0.6467 MP2t
0.2761 MP3
Density-based basis-set correction, see keyword denscorr,
using the built-in PBEc-based functional:
denscorr=2
calc=DF-CCSD(T)
dft=userd
1
1.00 HFx
1
1.0 ECMD_PBEc
Density-based basis-set correction, see keyword denscorr,
utilizing correlation functionals from the Libxc library:
denscorr=2
calc=DF-CCSD(T)
dft=userd
1
1.00 HFx
1
1.0 ECMD_GGA_C_PW91
This keyword can be used to select the wave function method for the evaluation of the second-order contributions in the case of excited-state double hybrid TD-DFT approaches. See Ref. 90 for further details.
The corresponding methods will be used for the evaluation of second-order terms.
dhexc=ADC(2)
If a spin-scaled DH is employed, the corresponding spin-scaled ADC(2) or CIS(D) method will be used even if ADC(2) or CIS(D) is set. For a non-spin-scaled DH, if a spin-scaled approach is selected, the default spin-scaling parameters described in Ref. 90 will be set. To change the spin-scaling parameters use keywords scsps and scspt. For the CIS(D)-based ansatz, the keywords scsps_t and scspt_t can also be varied.
Type of diagonalization algorithm used for the CI and LR-CC calculations.
Standard Davidson diagonalization
Davidson diagonalization with root-following, recommended for excited-state calculations if the initial guess is given manually or the calculation is restarted
diag=david
for root-following type diag=follow
Specifies the number of doubly occupied orbitals in an MCSCF calculation. See also the description of keyword mact.
The number of doubly occupied orbitals per irrep should be given in the
following format:
docc=$$
where $$ is the number of doubly occupied orbitals in irrep $i$,
and ${N}_{ir}$ is the number of irreps.
There is no default, the occupation must be set in the case of an MCSCF calculation.
Water, $2\times 2$ CAS calculation, the active space includes
an A${}_{1}$ and a B${}_{1}$ orbital:
docc=3,1,0,0
mact=1,0,1,0
Oxygen, $4\times 4$ CAS calculation, the active space includes
four orbitals of B${}_{2g}$, B${}_{3g}$, B${}_{2u}$, and B${}_{3u}$ symmetry.
Note that the multiplicity and the symmetry of the
state must also be specified, and, if the MCSCF calculation is
restarted from HF MOs, the HF occupation should also be given.
occ=3,0,1,1,0,2,1,1/3,0,0,0,0,2,1,1
docc=3,0,0,0,0,2,0,0
mact=0,0,1,1,0,0,1,1
mult=3
symm=2
Radius of atom domains for the local correlation method of Ref. 138 (localcc=2013). For each localized MO (LMO), using the Boughton–Pulay procedure [9], 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.
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.
Infinite radius will be applied, i.e., there is only one atom domain including all atoms.
domrad=10.0
To use the local CC methods as defined in Ref. 135 set domrad=inf, that is, use only one atom domain including all atoms.
to set a threshold of 12.0 bohr type domrad=12.0
Specifies the type of the algorithm for the solution of the dRPA equations or the calculation of SOS-MP2 energies. See Ref. 72 for more details.
The algorithm of Ref. 54 will be used, the fitting of integral lists will be performed before the dRPA iterations (SOS-MP2 calculation).
The algorithm of Ref. 72 will be executed, the fitting of the integrals is not performed. This algorithm is efficient for large molecules.
The dRPA correlation energy is calculated using the plasmon formula.
The algorithm is automatically selected on the basis of the size of the molecule (canonical dRPA) or the HOMO-LUMO gap (local dRPA).
For SOSEX calculations drpaalg=fit is the only option, which is forced by the program.
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.
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}$.
drpaalg=fit and drpaalg=auto for canonical and local dRPA, respectively.
to set the second option give drpaalg=nofit
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. 49 for more details. This feature can be combined with dual auxiliary basis, dual metric, dual grid, etc. approximations, see keywords grid_sm and dual_df.
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. 49. 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. 49. Similar to e1, but there is also an iterative step with the large basis.
dual=off
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
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
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
Set this keyword to perform dual auxiliary basis, dual metric DF-SCF, or dual grid COSX calculations [DualDF]. For dual grid DFT calculations, see keyword grid_sm.
Performs dual auxiliary basis SCF calculations. To specify the small auxiliary basis, set the dfbasis_scf_sm keyword.
Performs a dual metric DF-SCF calculation.
Performs a dual grid COSX calculation. To specify the small angular integration grid, set the agrid_pssp_sm keyword.
Turns off this feature.
dual_df=off for all calculations except COSX, where dual_df=grid is the default. For more information about COSX calculations, see the keyword pssp.
The dm option can be used in conjunction with the non-local direct SCF method only. The dual auxiliary basis algorithm can be applied to both local and non-local SCF calculations.
Dual auxiliary basis calculation with local fitting with the
cc-pVTZ AO basis set, the default large fitting basis set, and the
cc-pVTZ-RI-JK-dual small fitting basis set:
basis=cc-pVTZ
dfbasis_scf_sm=cc-pVTZ-RI-JK-dual
calc=DF-HF
scfalg=locfit2
dual_df=dab
Dual metric calculation with the cc-pVTZ AO basis set and
the default fitting basis set:
basis=cc-pVTZ
calc=DF-HF
scfalg=direct
dual_df=dm
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.
No ECPs will be used.
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.
If the same ECP is used for all atoms, the label of the ECP can be given here.
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$\mathrm{>}$:$$ECP label$\mathrm{>}$ .
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.
By default, the following ECP are available for elements Na to Rn in Mrcc:
Please note that some of the above ECPs are not available for all elements.
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.
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.
ecp=auto
To use the MCDHF-ECP-10 pseudopotential for all atoms the input must include ecp=MCDHF-ECP-10
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
This keyword controls the calculation of empirical dispersion corrections for DFT and HF calculations using the DFT-D3 approach of Grimme and co-workers [36, 37]. 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 environment variable.
No dispersion correction will be computed.
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).
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.
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 [37].
edisp=off
to calculate the D3 dispersion correction including BJ damping to the B3LYP energy give calc=B3LYP-D3
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
This keyword controls the calculation of empirical dispersion corrections for low-level methods (environment) in embedding calculations using the DFT-D3 approach of Grimme and co-workers [36, 37]. See the keywords edisp and embed for more details. Note that the empirical corrections are only calculated with the embedded atom list in the case of the subsystem calculations (see keyword embed). Note also that edisp_embed is not usable in the case of energy-based subsystem separation.
No dispersion correction will be computed.
The dispersion energy correction will be automatically evaluated for the low-level theory of the embedding approach. Other considerations are same as for edisp=auto.
You can directly give any options of the dftd3 code. Other considerations are same as for edisp=auto. Note that the coordinates will be taken from the COORD.xyz and the COORD_SUBSYSA.xyz files generated by Mrcc for the full system and the subsystem calculations, respectively.
If edisp_auto=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 [37].
edisp_embed=off
To calculate the D3 dispersion correction (including BJ damping) for
the interaction of the embedded subsystem and the environment
in the case of a LNO-CCSD(T)-in-BLYP type embedding give
calc=LNO-CCSD(T)
embed=huzinaga
1-2
BLYP-D3
0
To calculate the D3 dispersion correction to the
BLYP energy without the BJ damping, the input should include:
calc=LCCSD(T)
edisp_embed=-func b-lyp -zero
embed=huzinaga
1-2
BLYP
0
The electronic temperature in K at which the Fermi-smearing of the occupied orbitals is made in the case of scfiguess=gfn1 and scfiguess=gfn2.
$$any positive real number$\mathrm{>}$
temp=300.0
for 10000 K set etemp=10000.0
This keyword controls DFT embedding calculations. Currently, it is only available using density-fitting. Both restricted closed-shell determinants and unrestricted open-shell determinants are applicable for embedding calculations, however, only those open-shell cases are supported where the environmental subsystem has singlet multiplicity. See also keywords corembed and oniom for further embedding approaches.
No embedding.
The Huzinaga-equation-based embedding approach [48, 49] 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). Note that the user option of keyword dft is also supported. 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, $i$, should be given which is the number of the embedded orbitals or zero. If zero is given, then the number of embedded orbitals are 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, nmul, spad, nspad, ecore, bopu, and serialno. The amul and nmul mean the Mulliken population-based orbital selection scheme of Refs. 48 and 49. For amul, you can also specify the Mulliken population threshold (see the appendix of Ref. 49), the default threshold is 0.3. The nmul scheme selects $i$ pieces of active orbitals with the largest Mulliken populations. The spad and nspad techniques employ the orbital selection scheme of Ref. 12, where the active orbitals are selected based on the difference of the MO populations. Note that spad does not require a threshold for the selection, algorithm nspad selects $i$ orbitals which have the largest eigenvalues of the SPADE projector, and both schemes require orbloce=spade. The energy-based ecore algorithm is also supported, which does not require orbital localization and selects the high-energy orbitals as active orbitals: the ecore algorithm automatically selects the valence orbitals as active, while if $i$ is nonzero, then $i$ pieces of orbitals with the largest orbital energies are selected. Note that in the case of the ecore algorithm, the embedded atom list is not used in the calculation. bopu means the Boughton–Pulay algorithm-based scheme, see Ref. 48 for more details. For the default MO completeness criteria of the bopu algorithm, see keywords bpcompo and bpcompv. Lastly, if serialno is specified, the active orbitals are selected manually by giving the appropriate serial numbers of the orbitals after the keyword serialno. Note that the ordering of the orbitals are based on their orbital energies and this scheme does not require a list of active atoms. The default automatic selection scheme is spad if the orbloce=spade, and amul otherwise. If $i$ is not zero, nspad is the default scheme when orbloce=spade, and nmul is the default algorithm otherwise.
Finally, the embedding calculations can also be executed with selected active virtual orbitals if orblocv is defined by using orbloce=special. The orbital selection scheme of the virtual orbitals is specified in the fifth line, and all the previously described schemes are supported with a single modification: the ecore algorithm does not have an automatic selection mode and it will select the $i$ number of virtual orbitals with the lowest eigenvalues. In this case, a modified Huzinaga operator is constructed which includes the projector of the inactive virtual orbitals as well. Note that the subsequent correlation calculations will only use the active virtual subspace. See also examples below.
The projector-based embedding approach of Manby and co-workers [84] will be used. The embedded atoms and the low-level DFT approach can be specified as described above.
embed=off
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
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
The same as example 1, but the bopu algorithm is
used with 0.95 MO completeness criteria:
calc=DF-CCSD(T)
embed=huzinaga
1-2
PBE
0 bopu 0.95
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
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
B3LYP-in-BLYP embedding with the Huzinaga-equation-based
approach, the number
of the embedded orbitals will be the same as the number of
valance orbitals. Note that the active atoms are not specified
because of the ecore algorithm, and a blank line
is used in the second line after the keyword embed:
calc=B3LYP
embed=huzinaga
BLYP
0 ecore
B3LYP-in-BLYP embedding with the Huzinaga-equation-based
approach, the number
of the embedded orbitals will be the 4 MOs with the highest
orbital energies:
calc=B3LYP
embed=huzinaga
BLYP
4 ecore
B3LYP-in-BLYP embedding with the Huzinaga-equation-based
approach, orbitals with the 1, 2, 3, 4, and 13 serial numbers
are selected manually as active orbitals. Again, note that
the list of active atoms are not required in this case:
calc=B3LYP
embed=huzinaga
BLYP
0 serialno 1,2-4,13
CCSD(T)-in-PBE embedding with the Huzinaga-equation-based
approach, atoms 1 and 2 are included in the embedded region.
The spade
algorithm is used for the automatic selection of the active occupied
orbitals, while 4 active virtual orbitals are selected
by the nmul scheme. Note that the embedded occupied MOs
of the second SCF run are only allowed to mix with the
4 selected virtual MOs, and subsequently, the CCSD(T) calculation
is carried out with the active occupied and active virtual MOs:
calc=DF-CCSD(T)
orbloce=special
spade
spade
spade
embed=huzinaga
1-2
PBE
0
4 nmul
Same as example 9, but after the
embedded HF solution is obtained, projected atomic orbitals (PAO) are constructed
for the correlation calculation.
Note that the atomic orbitals of
atoms 1 and 2 will be used to calculate the
norm of the PAOs, and those PAOs that have
a norm lower than pao_subsys_tol will be omitted
in the subsequent correlation calculation.
Note the blank line below orblocv.
calc=DF-CCSD(T)
orblocv=pao-subsys
orbloce=special
spade
spade
spade
embed=huzinaga
1-2
PBE
0
4 nmul
Same as example 10, but this
time, the norm of the PAOs are calculated using an
extended number of atoms thus generating more
virtual orbitals.
calc=DF-CCSD(T)
orblocv=pao-subsys
1-2,8-10
orbloce=special
spade
spade
spade
embed=huzinaga
1-2
PBE
0
4 nmul
Same as example 11, except that
all virtual orbitals are used in the high-level SCF run,
and only the CCSD(T) calculation utilizes a reduced number
of virtual orbitals. Again, note the blank line
below orblocv, which can be used to
override the atom list which is used for the norm
calculation of PAOs. Here, the atom list
is not overridden.
calc=DF-CCSD(T)
orblocv=pao-subsys
orbloce=special
spade
spade
off
embed=huzinaga
1-2
PBE
0
Use this option to add an external perturbation to the Hamiltonian, e.g., an external electric dipole field.
No perturbations are added.
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$\mathrm{>}$ $$coefficient 1$\mathrm{>}$
$$operator 2$\mathrm{>}$ $$coefficient 2$\mathrm{>}$
$$operator 3$\mathrm{>}$ $$coefficient 3$\mathrm{>}$
…
where the operator can be x, y, z, xx, yy,
zz, xy, xz, yz, xxx, xxy,
xxz, xyy, xyz, xzz, yyy, yyz,
yzz, zzz.
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).
epert=none
to add the $\widehat{y}$ and $\widehat{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
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. If unset, keyword lnoepsv is used as the threshold for MP2 NOs. The use of keyword lnoepsv is preferred instead of keyword eps. See Ref. 136 for more details.
Virtual orbitals with cumulative populations of higher than this number will be dropped.
None. By default keyword lnoepsv is used.
This keyword can not be used with explicitly correlated calculations.
If eps is specified, the spatial symmetry is not used in correlation calculations.
to set a threshold of 0.95 type eps=0.95
Select the desired method for calculating atomic charges based on the electrostatic potential around the molecule.
No atomic charge calculation based on electrostatic potential.
Calculate CHELPG charges [breneman1990determining].
Calculate Merz–Singh–Kollman charges [singh1984approach].
espcharge=off
If espcharge$\mathrm{\ne}$off and dens=0, then it will be changed to dens=1.
Atomic radii used for the determination of the charges are listed in Table LABEL:radii. The values are based on the original articles and the Multiwfn implementation. For atoms not shown, the UFF radii divided by 1.2 are utilized.
to calculate CHELPG charges use espcharge=chelpg
Sets the radius of local fitting domains for the exchange contribution in direct density-fitting SCF calculations [127, 107]. 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.
Threshold for the integrals (in E${}_{h}$).
A threshold of zero will be applied, i.e., a conventional direct DF-SCF calculation will be executed.
Local fitting domains are available for RHF, UHF, ROHF, RKS, UKS, and ROKS calculations.
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.
excrad=1.0 if scfalg=locfit1 or scfalg=locfit2, excrad=0 otherwise
to set a threshold of 0.1 E${}_{h}$ type excrad=0.1
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.
See the description of keyword excrad.
excrad_fin=excrad/1000
to avoid the use of local fitting domains in the
extra iteration step give
excrad_fin=0.0
Specifies the fitting metric for direct DF-SCF calculations.
Coulomb-metric DF-SCF calculation.
Overlap-metric DF-SCF calculation.
fitting=Coulomb
to use the overlap metric, set fitting=Overlap
Controls the multipole method for three-center Coulomb integrals for density fitting SCF calculations. See also keyword fmmord.
No multipole approximation is used for the integrals.
The exchange contribution will be evaluated using multipole approximations.
Both the Coulomb and the exchange contributions will be evaluated using multipole approximations.
The multipole method can be employed either with or without local fitting approximations (see keyword scfalg). If no local fitting is used, the speedup is quite good, the computation time can roughly be halved. With local fitting approximations, the speedup is moderate.
fmm=off
To use the multipole method to speed up a DF-HF calculation, set fmm=on
The keyword controls the order of the multipole approximation if the multipole method is used for three-center Coulomb integrals (see also keyword fmm). If the multipole approximation is invoked to calculate an integral of the product of two Gaussians, the contribution of multipole moments up to order ${l}_{1}$ $+$ ${l}_{2}$ $+$ fmmord will be considered, where ${l}_{1}$ and ${l}_{2}$ are the angular momentum quantum numbers of the Gaussians.
$$any positive integer$\mathrm{>}$
fmmord=8
for lower accuracy but increased speed, give fmmord=6
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 $\ne $ off, the frequencies are calculated at the optimized geometry.
on or off
You should also set this keyword if you are interested in thermodynamic properties of atoms.
freq=off
for a frequency calculation set freq=on
Exponent for the Slater-type correlation factor ($\gamma $) used in explicitly correlated calculations [CorrFactTew, MP2-F12_TM].
Value of the correlation factor in bohr${}^{-1}$.
In Mrcc, six Gaussians are used to fit the Slater-type correlation factor. The exponents and coefficients were taken from Ref. CorrFactTew.
gamma=0.9 for cc-pVDZ-F12, gamma=1.1 for cc-pVQZ-F12 or aug-cc-pVTZ, gamma=1.4 for aug-cc-pVQZ or aug-cc-pV5Z, gamma=1.0 otherwise.
to set an exponent of 1.3 bohr${}^{-1}$ give gamma=1.3
Specifies whether spherical harmonic or Cartesian Gaussian basis functions will be used.
Spherical harmonic Gaussians will be used
Cartesian Gaussians will be used
For calculations using the density fitting (DF) approximation, if intalg=os or intalg=auto, the Coulomb integrals are evaluated by algorithms [2, 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).
The derivative integrals are evaluated by the solid-harmonic Hermite scheme [131] (see the description of keyword intalg), consequently, differentiated integrals, and thus energy derivatives cannot be evaluated with Cartesian Gaussian basis sets.
gauss=spher
for Cartesian Gaussians the user should set gauss=cart
Specifies the format of molecular geometry. The geometry must be given in the corresponding format in the subsequent lines.
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”.
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).
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.
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.
For the use of ghost atoms see the description of keyword ghost.
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.
the following five geometry inputs for H${}_{2}$O${}_{2}$ are equivalent
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
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
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
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
.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 atoms can be specified using this keyword, e.g., for the purpose of basis set superposition error (BSSE) calculations with the counterpoise method.
There are no ghost atoms.
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.
ghost=none
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
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
Requests geometry optimization. Currently only the full geometry optimization is supported, geometrical parameters cannot be frozen.
no geometry optimization.
full geometry optimization.
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.
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.
gopt=off
to carry out a full geometry optimization set gopt=full
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.
Densities are not evaluated.
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})$
grdens=off
to save the densities give grdens=on
Specifies the integration grid for DFT calculations. The grid construction follows the design principles of Becke [6] in any case.
An automatically constructed adaptive integration grid will be used. The adaptive scheme is motivated by Ref. 77. The important difference is that the grids are optimized for each atom separately to avoid discontinuous potential energy surfaces. Both the angular and the radial grids are selected automatically to ensure the accuracy of the angular and the radial integrations as described at keywords agrid and rgrid. See also the description of keyword grtol. The smoothing function for the Voronoi polyhedra are adopted from Ref. 104 with ${m}_{\mu}$ = 10.
Grid no. $n$ of Treutler–Ahlrichs [AhlrichsGrid] using the Log3 radial quadrature of Mura and Knowles [103].
grid=TA3 for KS-DFT calculations with non-hybrid LDA and GGA functionals, grid=auto otherwise
to use grid 2 of Treutler–Ahlrichs, set grid=TA2
If this keyword is set and different from keyword grid, a dual grid DFT calculation will be carried out. The SCF iterations will be performed using the grid specified by grid_sm, and an extra iteration will be carried out with a larger grid set by grid. A first-order correction is also evaluated [DualDF]. This can also be combined with dual (auxiliary) basis, dual metric, etc. approximations, see keywords dual and dual_df.
see keyword grid.
grid_sm=grid
To use grid 1 and 3 of Treutler–Ahlrichs in the SCF cycles and in the
final iteration, respectively, (recommended for KS-DFT calculations with
non-hybrid functionals) set:
grid_sm=TA1
grid=TA3
To use grid 1 of Treutler–Ahlrichs and a Coulomb fitting basis in the
SCF cycles and the default adaptive grid and a Coulomb-exchange fitting
basis in the final iteration (recommended for DFT calculations with a
hybrid functional; note that for the latter, you do not need to set
dfbasis_scf with def2-type AO basis sets, because this is the
default) set:
grid_sm=TA1
dual_df=DAB
dfbasis_scf_sm=def2-QZVPP-RI-J
dfbasis_scf=def2-QZVPP-RI-JK
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}^{-\mathrm{\U0001d690\U0001d69b\U0001d69d\U0001d698\U0001d695}}$, while the number of radial grid points increases linearly with grtol. See also the description of keywords agrid and rgrid.
$$any positive integer$\mathrm{>}$
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.
grtol=10
for a fine integration grid give grtol=12
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}^{-\mathrm{\U0001d690\U0001d69d\U0001d698\U0001d695}}$ bohr.
$$any integer$\mathrm{>}$
gtol=7
for a tolerance of ${10}^{-4}$ bohr give gtol=4
Use this keyword to calculate the Guido-type $\mathrm{\Delta}\mathbf{r}$ charge-transfer descriptor (see Ref. [Guido_CT]).
Calculate the corresponding descriptor.
The calculation of the descriptor is turned off.
guido_ct=off
To calculate the descriptor use guido_ct=on
Specifies what type of Hamiltonian is used in relativistic calculations. This keyword has only effect if iface=Dirac.
exact 2-component molecular-mean-field Hamiltonian [143]
other types of relativistic Hamiltonians such as the full Dirac–Coulomb Hamiltonian or the exact 2-component Hamiltonian
hamilton=DC
if you use the exact 2-component
molecular-mean-field Hamiltonian, set
hamilton=X2Cmmf
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.
Transformed MO integrals are calculated by Mrcc.
Mrcc is interfaced to Cfour.
Mrcc is interfaced to Columbus.
Mrcc is interfaced to Dirac.
Mrcc is interfaced to Molpro.
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.
If you use Mrcc together with Columbus or Dirac, this keyword must be always given.
iface=none, that is, all calculations will be performed by Mrcc.
to carry out four-component relativistic calculations using the Dirac interface give iface=Dirac
Specifies the algorithm used for the evaluation of two-electron integrals over primitive Gaussian-type orbitals.
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.
The integrals over contracted Gaussians are evaluated by the solid-harmonic Hermite scheme of Reine et al. [131].
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.
For DF methods option herm is not available.
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).
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.
intalg=auto
to use the Obara–Saika scheme for all angular momenta add intalg=os
Specifies whether standard excitation energies, ionization potentials, or electron affinities are calculated. This option is available at the ADC(2) and CIS(D) levels, as well as for the corresponding TDDFT functionals [IP_EA].
Standard excited-state calculations are carried out.
Ionization potentials are calculated.
Electron affinities are calculated.
ip_ea=off
for ionization potentials, the user should set ip_ea=ip
Threshold for integral calculation. Integrals less than ${10}^{-\mathrm{\U0001d692\U0001d69d\U0001d698\U0001d695}}$ E${}_{h}$ will be neglected.
$$any integer$\mathrm{>}$
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)
for an accuracy of ${10}^{-15}$ E${}_{h}$ one must give itol=15
Specifies the accuracy of the numerical Laplace transform in the (T) correlation energy term of local CCSD(T) calculations with localcc$\mathrm{\ge}$2016. See also the description of keyword talg and Ref. 105.
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}(\text{\U0001d695\U0001d68a\U0001d699\U0001d69d\U0001d698\U0001d695})|$.
laptol=${\mathrm{10}}^{\mathrm{-}\mathrm{2}}$, but laptol=${\mathrm{10}}^{\mathrm{-}\mathrm{3}}$ is set if lcorthr=tight. See also the description of lcorthr for further details.
to use a threshold of ${10}^{-3}$ type laptol=1e-3
Specifies the order of operations in LNO-CC calculations.
The integral transformation steps yielding the LNO integral lists of a LIS are immediately followed by the CC calculation in the same LIS. The three-center LNO integrals are stored in memory and passed directly to the ccsd program (unless ccsdalg=disk was requested). Hence, in combination with the default ccsdalg=dfdirect, the LNO integrals can be kept in memory throughout the entire LNO-CC calculation, and disk I/O is completely avoided in these steps. This option is only available with localcc$\mathrm{\ge}$2018 and ccprog=ccsd.
The integral transformation steps to the LNO bases of the LISs are completed first for all LISs, which is followed by the CC calculations performed for each LIS in their respective LNO basis. This option requires storing the integral lists of all LISs on disk. It is usually not limiting to store the three-center LNO integrals required for ccprog=ccsd or the four-center LNO integrals required for ccprog=mrcc.
lccoporder=lccfirst for localcc$\mathrm{\ge}$2018,
lccoporder=trffirst for localcc$$2018 or ccprog=mrcc.
to force storing the integrals on disk for localcc$\mathrm{\ge}$2018 set lccoporder=trffirst
Use this keyword to restart local CC calculations in the case of localcc$\mathrm{\ge}$2015, e.g., after power failure. If localcc$\mathrm{\ge}$2018 (and lccoporder=lccfirst) the calculation can be restarted at any stage of the local correlation calculation via lccrest=restart. In this case all the intermediate files located in the folder of execution have to be kept intact, and only the updated MINP file (containing lccrest=restart) should be overwritten before the calculation is restarted. The HF calculation and the orbital localization is skipped via the scfiguess=off and the orblocguess=read options, therefore it is important not to modify the VARS, SCFDENSITIES, FOCK, and MOCOEF* files either. 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.
For the restart with lccrest=on, compatible with options localcc=2015, localcc=2016, or localcc$\mathrm{\ge}$2018 and lccoporder=trffirst, the LMP2 calculation and the integral transformations have to be completed, and for the remaining domains the localcc.restart, DFINT_AI.*, DFINT_AB.*, DFINT_IJ.*, ajb.*, 55.*, 56, and VARS files are required, and only the updated MINP file (containing lccrest=on) should be overwritten before the calculation is restarted. In the case of localcc$\mathrm{\ge}$2016, talg=lapl or topr, and lccoporder=trffirst the laplbas.* files are also needed.
In the case of an LMP2 calculation or the LMP2 part of an LNO-CC calculation with lccoporder=trffirst and localcc$\mathrm{\ge}$2016, if the loop for the extended domains is started but not finished, the calculation can be restarted with the lccrest= domain option. In this case all the intermediate files located in the folder of execution have to be kept intact, and only the updated MINP file (containing lccrest=domain) 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.
In the case of restricted open-shell local correlation computations, the FOCK.OD file is also needed. For qro=on, the MOCOEF.QRO file should also be preserved.
restart, on, off, or domain.
lccrest=off
to restart a local CCSD(T) calculation set
lccrest=restart if localcc$\mathrm{\ge}$2018 (and lccoporder=lccfirst).
In the case of localcc<2018 and/or
lccoporder=trffirst set
lccrest=on or lccrest=domain depending on the point of interruption.
Controls the accuracy of local correlation calculations by setting the relevant thresholds: bp*, lnoepso, lnoepsv, laptol, naf_cor, osveps, wpairtol, spairtol (see also Refs. 72, 107, and 108 for details).
Very loose thresholds will be used, see Table 2.
Very Tight thresholds will be used, see Table 1.
Even tighter or very very Tight thresholds will be used, see Table 1.
The truncation thresholds will be set so that the canonical energy be reproduced, it is only useful for testing.
The values of the thresholds controlled by lcorthr are summarized in Table 1 for the localcc$\mathrm{\ge}$2018 case and in Table 2 for the previous schemes. The performance of lcorthr=Loose to vTight settings for localcc=2018 is extensively benchmarked in Ref. 106. The vLoose and vvTight settings may be employed for convergence studies but their accuracy/cost performance is mostly not practical for production calculations.
Expected accuracy. Using the Normal settings for local MP2 and
CCSD(T), if localcc$\mathrm{\ge}$2018, 1 kJ/mol (1 kcal/mol) average (maximum)
errors are expected in energy differences even for relatively complicated or sizable systems
(see Refs. 107, 108, and 106).
Note that molecules with complicated electronic structure can exhibit
slower convergence with the local approximations and higher local errors [106, L7CCDMC].
Thus, converge tests with increasingly tighter lcorthr settings
are suggested at least for some representative systems
to determine appropriate local correlation settings and for
the robust estimation of the local approximation errors [106].
In the case of monotonic convergence with respect to lcorthr,
the corresponding LNO approximations of Ref. 106 are also recommended
to further increase the accuracy.
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. 72).
The restricted open-shell local MP2 approach is available via localcc$\mathrm{\ge}$2021 [ROLocalMP2]. Currently, only MP2 and the corresponding double-hybrid DFT computations are available with our local approximations. ROHF, ROKS, and in combination with unrestricted references, also QRO references can be employed. With UHF (UKS) references, qro=on is set for localcc$\mathrm{\ge}$2021.
For localcc>2018, the NAF approximation is not employed for the local MP2 correlation energy of calc=LMP2, thus this local MP2 result is exactly the same as the byproduct local MP2 energy of an LNO-CCSD(T) computation with localcc>2018. The NAF approximation was employed in Ref. 107, thus naf_cor=2e-3 is set in combination with localcc=2018 and calc=LMP2, leading to a slightly different local MP2 correlation energy than the one obtained within a localcc=2018 type LNO-CCSD(T) computation (where the NAF approach is only employed at the CC level but not for MP2). For local MP2 naf_cor=off, is set if localcc=2015. Keywords lnoepso, lnoepsv, and laptol are irrelevant for local MP2 calculations.
The closed-shell localcc=2018 and localcc=2021 variants are very close but not always identical. In the case of off-center atomic orbitals employed via the ghost keyword, the domain construction scheme was improved, and the localcc=2021 version is considerably more accurate than localcc=2018.
calc | LMP2${}^{a}$ or LNO-CC | all | |||||
---|---|---|---|---|---|---|---|
localcc | 2018 or 2021 | ||||||
lcorthr | vLoose | Loose | Normal | Tight | vTight | vvTight | 0 |
bpedo | 0.9999 | 0.9999 | 0.9999 | 0.99995 | 0.99999 | 0.99999 | 1.0 |
wpairtol | 1e-4 | 3e-5 | 1e-5 | 3e-6 | 1e-6 | 3e-7 | 0.0 |
lnoepso | 1e-4 | 3e-5 | 1e-5 | 3e-6 | 1e-6 | 3e-7 | 0.0 |
lnoepsv | 1e-5 | 3e-6 | 1e-6 | 3e-7 | 1e-7 | 3e-8 | 0.0 |
laptol | 1e-1 | 1e-1 | 1e-2 | 1e-3 | 1e-3 | 1e-3 | 1e-6 |
naf_cor | 1e-2 | 1e-2 | 1e-2 | 1e-2 | 1e-3 | 1e-3 | 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 |
^{a} See note 4 for lcorthr for additional details of LMP2.
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 |
lcorthr=Normal
to use tight thresholds set lcorthr=Tight
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).
The MP2 density matrix fragments are calculated using the correct, non-canonical expressions.
The MP2 density matrix fragments are calculated using the approximate canonical expressions (as defined in Ref. 135).
To reproduce the method described in Ref. 135 use lmp2dens=off.
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.
lmp2dens=on
to use the canonical expressions give lmp2dens=off
Threshold for the occupation numbers of occupied local natural orbitals (LNOs) in the case of local correlation calculations, for state-averaged MP2/CIS(D) occupied natural orbitals for reduced-cost excited-state calculations, or for MP2 occupied natural orbitals for reduced-cost correlation calculations, see also keyword lnoepsv. See Ref. 138 as well as Refs. 93 and 94 for more details.
Orbitals with occupation numbers greater than 1-lnoepso will be frozen.
For default settings with options other than lcorthr=Normal, see the description of lcorthr.
For FNO-CC calculations, occupied NOs are only implemented for closed-shell references with density fitting.
Frozen occupied NOs are not implemented for explicitly correlated methods.
lnoepso=1e-5 for local correlation, lnoepso=0 for excited-states
to set a threshold of $5\cdot {10}^{-6}$ type lnoepso=5e-6
Threshold for the occupation numbers of virtual local natural orbitals (LNOs) in the case of local correlation calculations, for state-averaged MP2/CIS(D) virtual natural orbitals for reduced-cost excited-state calculations, or for MP2 virtual natural orbitals for reduced-cost correlation calculations, see also keyword lnoepso. See Ref. 138 as well as Refs. 93, 94, FNONAF, and RedCostF12 for more details.
Orbitals with occupation numbers smaller than this number will be dropped.
For default settings with options other than lcorthr=Normal, see the description of lcorthr.
lnoepsv=1e-6 for local correlation; lnoepsv=7.5e-5 for excited states if no spin scaling is applied; lnoepsv=3e-5 for spin-scaled excited-state models (such as SCS-CC2, SOS-CC2, …); lnoepsv=5e-5 for reduced-cost (FNO) ground-state calculations if no explicit correlation is used. For explicitly correlated CC calculations, lnoepsv=1e-4 with the cc-pVDZ-F12 and aug-cc-pVDZ bases, lnoepsv=1e-6 with cc-pVQZ-F12 and aug-cc-pVQZ, and lnoepsv=1e-5 otherwise.
to set a threshold of $5\cdot {10}^{-7}$ type lnoepsv=5e-7
Specifies if local correlation calculation is performed. See Refs. 135, 138, 72, 107, 105, and 108 for more details.
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).
localcc=off
for local correlation calculations give localcc=on
Controls the point group symmetry recognition feature of ground-state local correlation methods with localcc$\mathrm{\ge}$2018 as discussed in Sect. 3.4 of Ref. 108.
Point group symmetry is not exploited for the local correlation computation.
The symmetry equivalent domains are recognized by analyzing both the diagonal and off-diagonal elements of the Fock matrix transformed to the localized MO basis.
The symmetry equivalent domains are determined only based on the diagonal Fock matrix elements. This is a less strict criterion than localcorrsymm=on.
The degeneracy levels of each localized MO can be read from the MINP input file. The domain correlation energy contribution of each localized MO will be multiplied by the positive integer degeneracy level in the MO order determined by the MO localization program. Orbitals without symmetry equivalent MOs should have a degeneracy level of 1. Negative integers denote that the correlation energy contribution of the corresponding MO is already accounted for in the domain of a symmetry equivalent MO. The correctness of these degeneracy levels is not checked by the program. The orblocguess=read or restart setting is recommended with the option to retain the order of the localized MOs.
localcorrsymm=off
While Abelian point group symmetry is recognized and exploited by the HF/KS-DFT programs, the orbital localization algorithms are not yet symmetry adapted. Hence the localized MOs may not exhibit (all of) the symmetry properties of the molecular structure and canonical MOs. Because of that, the localcorrsymm feature is turned off by default even if localcc$\mathrm{\ne}$off. If localcorrsymm$\mathrm{\ne}$off is set, it is recommended to manually check the degeneracy level determination in the output section “Determine symmetry equivalent LMO lists”.
This point group symmetry identification approach [108] can also take advantage of non-Abelian symmetry if the canonical and localized MOs retain the corresponding symmetry properties even if (non-Abelian) symmetry is not enforced by the SCF and localization algorithms.
to read the list of symmetry-equivalent MOs
for ${C}_{2v}$ water:
localcorrsymm=read
2 -2 2 -2
Specifies the number of active orbitals in an MCSCF calculation. See also the description of keyword docc.
The number of active orbitals per irrep should be given in the
following format:
mact=$$
where $$ is the number of active orbitals in irrep $i$,
and ${N}_{ir}$ is the number of irreps.
There is no default, the active orbitals must be given in the case of an MCSCF calculation.
See the description of keyword docc for examples.
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.
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.
maxact=off
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
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.
$$any positive integer$\mathrm{>}$
maxdim=max(200, 3*max(ntrip, nsing))
for a maximum of 50 expansion vectors set maxdim=50
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. 70 and 60 for more details.
The excitation manifold is not truncated.
The excitation manifold is truncated at $n$-fold excitations, see above.
maxex=0
to truncate the excitation manifold at triple excitations set maxex=3
Maximum number of microiterations in a Newton step in the case quadratic SCF calculations (see keyword qscf).
$$any positive integer$\mathrm{>}$
maxmicroit=100
to increase the maximum number of microiterations to 200 give maxmicroit=200
Initial guess for the MCSCF calculation.
The sad, ao, core, mo, and restart options can be used as described for keyword scfiguess. In addition, there is one further option:
An RHF or ROHF calculation is run before the MCSCF calculation, and the HF MOs will be used as initial guess in the MCSCF calculation.
To restart the MCSCF calculation using MOs from a previous MCSCF run, e.g., when calculating potential energy surfaces, use mcscfiguess=mo
mcscfiguess=HF
set mcscfiguess=sad to start the MCSCF calculation from SAD initial guess.
Specifies the core memory used.
The amount of memory to allocate is specified in megabytes
The amount of memory to allocate is specified in gigabytes
mem=256MB
to allocate 8 GB core memory the user should set mem=8GB
Specifies the molecular mechanics program to be used for ONIOM calculations. See the description below for the supported methods and programs. Currently, it is not necessary to specify the keyword as Mrcc automatically checks the program availability for the specified molecular mechanics method.
molecular mechanics method is not requested.
The xTB package of Grimme [xtbprog].
The default for GFN-FF is xtb, otherwise mmprog=off.
to select the xTB program for a calculation, set mmprog=xtb.
Specifies whether input file for the Molden program and an xyz-file containing the Cartesian coordinates are written (see also Sect. 14).
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.
The construction of the MOLDEN input and the COORD.xyz file is turned off.
molden=on
if you do not need Molden input and the COORD.xyz file, add molden=off
Specifies the number of MPI tasks spawned by dmrcc. Currently only executables scf, mrcc, and ccsd can be spawned for MPI parallel execution. All of the current MPI implementations replicate all arrays stored in memory for each MPI task. The mem keyword sets the available memory for a single MPI task. Hence, mem should be set to the total available memory of the compute node divided by the number of MPI tasks used on a single compute node.
$$any positive integer$\mathrm{>}$
mpitasks=1
For optimal performance, please set mpitasks to the total number of available CPUs/NUMA nodes/nodes/cores, etc., as the additional number of processes on top of mpitasks are driver processes running mostly in the background and do not require dedicated resources.
to run an MPI-parallel CC calculation using $2$ MPI tasks set mpitasks=2
Specifies the multipole approximation used for the evaluation of pair energies of distant pairs.
The simplified dipole-dipole estimate of Ref. 134 will be used.
Full dipole-dipole estimate [125].
All the terms are included in the multipole expansion up to the contributions of the quadrupole moment [107].
All the terms are included in the multipole expansion up to the contributions of the octapole moment [107].
mulmet=3 if localcc$\mathrm{\ge}$2016, mulmet=0 otherwise
to use the octapole approximation set mulmet=3
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.
$$any positive integer$\mathrm{>}$
for atoms the corresponding experimental multiplicity is set, for molecules mult=1 (singlet) for an even number of electrons, mult=2 (doublet) otherwise
for a triplet state one should give mult=3
Specifies whether the natural auxiliary basis (NAB) approximation will be used for explicitly correlated CC calculations and also specifies the threshold for the occupation numbers of NAB functions (see Ref. RedCostF12).
the NAB approach will not be used.
A NAB will be constructed and auxiliary functions with occupation numbers smaller than this number will be dropped.
Sets nab=0.1
nab=0.1 for reduced-cost explicitly correlated calculations (that is, if the FNO- prefix is added to the acronym of the selected method), nab=off otherwise.
to use the NAB approximation and set a threshold of $0.2$, type nab=0.2
Number of active occupied spinorbitals for multi-reference (active-space) CI/CC calculations. 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. 70 for more details. See also keywords nactv, maxex, and maxact. Note also that this keyword only sets the active orbitals for the post-SCF calculation, the MCSCF active orbitals can be specified by keyword mact.
$$any positive integer$\mathrm{>}$ or 0
If MCSCF reference is used, nacto is set so that the MCSCF active orbitals (given by keyword mact) will also be the active orbitals in the post-SCF calculation. nacto=0 in the general case.
for two active occupied spin-orbitals give nacto=2
Number of active virtual spinorbitals for multi-reference (active-space) CI/CC calculations. 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.
$$any positive integer$\mathrm{>}$ or 0
If MCSCF reference is used, nactv is set so that the MCSCF active orbitals (given by keyword mact) will also be the active orbitals in the post-SCF calculation. nactv=0 in the general case.
for two active virtual spin-orbitals give nactv=2
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. 71 for the definitions). The latter option is somewhat more efficient but can be dangerous for processes involving atoms.
NAFs are constructed from (${\mathbf{W}}^{\alpha}$ $+$ ${\mathbf{W}}^{\beta}$)/2.
NAFs are constructed from ${\mathbf{W}}^{\alpha}$.
nafalg=albe
to use ${\mathbf{W}}^{\alpha}$, set nafalg=alpha
This keyword can be used to control the state-averaged density construction for local ADC(2) calculations. See Ref. 92 for further details.
The NAF approximation is not used for the density construction.
A fast, local domain-based algorithm is used for the density construction invoking the NAF approximation.
nafdens=off
If redcost_exc=on and localcc=2016 settings are used, nafdens=on is set by default.
If local CVS-ADC(2) calculations are carried out, all the NAFs are retrieved during the density construction.
to switch on the NAF approximation at the construction of density matrices, set nafdens=on
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. 71).
NAFs will not be constructed.
A NAF basis will be constructed and NAFs with occupation numbers smaller than this number will be dropped.
Sets naf_cor=5e-2 when ccprog=ccsd is set or naf_cor=5e-3 in other cases if the use of NAFs is not specified as default. If NAFs are used as default, e.g., in the cases of certain local correlation and reduced-cost excited-state approaches, naf_cor=on should not be written to the input file.
If ccprog=mrcc or dfbasis_cor=none, no correction for the dropped NAFs is calculated.
according to the value of lcorthr for local correlation methods; naf_cor=0.1 for the reduced-cost excited-state approaches of Refs. 93 and 94 (see keyword redcost_exc) without spin scaling; naf_cor=0.075 for the spin-scaled version of the latter methods (such as SCS-CC2, SOS-CC2, …); naf_cor=0.05 for reduced-cost ground-state calculations (that is, if theFNO- prefix is added to the acronym of the selected method); naf_cor=off otherwise.
to use NAFs and set a threshold of ${10}^{-2}$, type naf_cor=1e-2
Specifies whether natural auxiliary functions (NAFs) will be used for explicitly correlated CC calculations at the calculation of the F12-dependent intermediates and also specifies the threshold for the occupation numbers of NAFs (see Ref. RedCostF12).
NAFs will not be constructed.
A NAF basis will be constructed and NAFs with occupation numbers smaller than this number will be dropped.
Sets naf_f12=5e-2
naf_f12=5e-2 for reduced-cost explicitly correlated calculations (that is, if the FNO- prefix is added to the acronym of the selected method), naf_f12=off otherwise.
to use NAFs and set a threshold of ${10}^{-2}$, type naf_f12=1e-2
Specifies whether NAFs will be used for density-fitting SCF calculations and also specifies the threshold for the occupation numbers of NAFs (see Ref. 71). The syntax is analogous with that for keyword naf_cor.
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. 71 and 108 for details).
NAFs are constructed from ${J}_{ai}^{P}$
NAFs are constructed from ${J}_{pi}^{P}$
NAFs are constructed from ${J}_{pq}^{P}$ using the $p\ge q$ values.
NAFs are constructed from all ${J}_{pq}^{P}$ integrals. This option generates NAFs invariant to orbital rotations among the MOs but it overweighs the $p\ne q$ integrals compared to the diagonal $p=q$ integrals, leading to a minor numerical deviation compared to the naftyp=Jpq option.
NAFs are constructed from ${J}_{\mu i}^{P}$
to construct NAFs using ${J}_{pq}^{P}$ set naftyp=Jpq
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.
The number of Cholesky vectors/quadrature points will be automatically determined to achieve the required precision.
The number of Cholesky vectors/quadrature points will also be automatically determined but the maximum number of the vectors cannot exceed this number.
nchol=auto
to use ten Cholesky vectors/quadrature points give nchol=10
Step size for the numerical differentiation in atomic units.
$$any positive real number$\mathrm{>}$
ndeps=1e-3
for a step size of $5\cdot {10}^{-4}$ a.u. for numerical Hessian evaluation set ndeps=5e-4
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, ntrip, and optex.
$$any positive integer$\mathrm{>}$
nstate=max(1, nsing+ntrip)
for three states give nstate=3
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. 69. This enables us to search for singlet and triplet roots separately. See also keywords nstate and ntrip.
$$any positive integer$\mathrm{>}$
nsing=1 for closed-shell reference determinants, nsing=0 otherwise
for two singlet states give nsing=2
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.
$$any positive integer$\mathrm{>}$
ntrip=0
for two triplet states give ntrip=2
If set, natural transition orbitals (NTOs) or state-averaged natural orbitals (NOs) are written (see also Sect. 14).
NTOs or NOs are saved to the corresponding files (see Notes). These files can be handled using the same programs as file MOLDEN can (see keyword molden).
The construction of the corresponding NTO or NO files is turned off.
For the NTOs, in the case of CIS/TDHF/TDDFT/TDA-TDDFT calculations, the MOLDEN_CIS_NTO.$\mathrm{n}$ files are generated for each state, where $n$ is the serial number of the excited state according to energy order.
For the NTOs, in the case of ADC(2)/CC2/CIS(D${}_{\mathrm{\infty}}$) calculations, the
MOLDEN_CORR_NTO.$\mathrm{n}$ files are generated. For a canonical
calculation, the files are labeled according to the order of the
excitation energies of the second-order method. In contrast, if a
state-specific, for example, reduced-cost calculation is performed, the
files are labeled according to the order of the CIS excitation
energies. Note that the ADC(2) and CIS roots can be swapped. You can
check at the end of the output in the summary if such rearrangement
occurred.
For the NTOs, in the case of double hybrid calculations, the MOLDEN_DH_NTO.$\mathrm{n}$
files are generated. In the case of the
ADC(2)-based ansatz, the files are
labeled according to the order of the excitation energies. In contrast,
if the CIS(D)-based ansatz is used, the files are labeled according
to the order of the TDA excitation
energies.
For the second-order methods, NTOs are calculated using only single excitation coefficients.
nto=off
if you need the NTO or NO files, add nto=on
This keyword controls the calculation of the numerical gradient of the energy with respect to nuclear positions.
Calculates the numerical gradient with a four point central difference scheme with step size of $0.01$ bohr. Please note that these settings may not be suitable for all types of calculations and molecules. Currently there is no way to control the step size or the numerical stencil from the input file.
Mrcc does not calculate numerical gradients.
numgrad=off
Specifies the occupation of the Hartree–Fock determinant.
If this keyword is not given, the occupation is automatically determined in the SCF calculations.
For RHF calculations the occupation should be given in the following
format:
occ=$$
where $$ is the number of occupied orbitals in irrep $i$, and
${N}_{ir}$ is the number of irreps.
For ROHF and UHF calculations the occupation should be given as
occ=$$/$$
where $$ is the number of occupied $\sigma $ spinorbitals in irrep $i$.
occ is not specified, that is, the occupation is set by the SCF program.
Water, RHF calculation:
occ=3,1,1,0
Water, UHF calculation:
occ=3,1,1,0/3,1,1,0
Carbon atom, ROHF or UHF calculation:
occ=2,0,0,0,0,1,0,1/2,0,0,0,0,0,0,0
Turns on the OCC-RI-K [occri, DualDF] algorithm. It can be used together with local and non-local direct DF-SCF calculations. It is available for both Coulomb- and overlap-metric DF-SCF methods. The desired fitting metric can be selected with the fitting and dual_df keywords.
Turns on the OCC-RI-K algorithm.
Turns off the OCC-RI-K algorithm.
It is turned off except for overlap-metric SCF calculations.
to turn of the OCC-RI-K algorithm, set occri=on
Initiates an ONIOM calculation [86, 57] and specifies the number of layers. In an ONIOM calculation, the system is split up into several layers described at different levels of theory. Both mechanical and electronic embeddings are supported. In the case of the mechanical embedding, the presence of the lower layers are omitted in upper-layer calculations, while a point charge representation of the lower layers is used in electrostatic embedding. The top layer (model system, active subsystem) is treated at the highest level, and the corresponding method, basis set, …are set by the other keywords in the MINP file, calc, basis, …, as usual. The atoms and the boundary handling schemes of the upper layers, and the methods of the lower layers must be specified in the subsequent lines as follows, starting from the lowest-level layer. For each layer three lines are required.
The atoms of the layer 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 atoms $$ through $$ are part of the layer. 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. Atoms of an upper-level layer have to be a subset of the lower-level layers.
The method for the lower-level layer must be specified in the second line using the corresponding options of keyword calc, or the supported methods of external programs (see keywords sqmprog and mmprog). The highest-level method is specified by keyword calc (or keywords calc and dft if a KS reference is used in a correlation calculation). The specification of the basis set (for the low-level method), charge, and multiplicity of the lower-level layer is optional but can be given in the second line separated by forward slashes and spaces, e.g., as PBE / cc-pVDZ / -1 / 1. See keywords basis, charge, and mult for the available options, respectively. Note that keywords basis, charge, and mult specify the basis of the last calculation, the charge, and the multiplicity of the top layer, respectively. If the basis set is not specified under oniom, the basis set of the upper-level layer will be used, while the default options will be applied for the charge and the multiplicity if they are not set here.
In the third line, an integer should be given which specifies the border handling scheme: the automatic (0), the semi-automatic (1), and the manual (2) border handling schemes are supported.
In the case of the automatic handling of layer borders, the following procedure is executed if the method of the lowest layer is not calculated by an external program. The MOs are localized after the first, full-system calculation (see orbloce for the available options for the localization scheme), and the localized MOs are assigned to maximum of two atoms by the Boughton–Pulay algorithm [9]. If the completeness criteria of a truncated LMO is greater than 0.8, the bond order between the corresponding two atoms is increased by one. The borders of the layers are then determined by those bonds which are located on a layer and a non-layer atom. Depending on their integer bond order, a hydrogen, an oxygen, or a nitrogen atom is selected as link atom. The bond distance is set based on the type of the link atom and the layer atom (see the average bond distances at https://cccbdb.nist.gov/expbondlengths1.asp), and the link atom is placed in the direction of the original bond. Note that the implementation follows the convention of the Amber program, where the distances between the link atoms and the layer atoms are fixed independently of the geometry. Note also that quadruple bonds are not separable automatically, and the separation of aromatic structures should be avoided. If the method of the lowest-level layer is handled by an external program, the same procedure is executed except that the the bonding table is not built by using the Boughton–Pulay algorithm. Instead, the bond table is based on the bond order table of the external program.
In the case of the manual border handling, the user has to
specify the atom pairs of the borders, the link atoms, and the
distance between the link atoms and the layer atoms.
These specifications can be given in
three input sections separated by forward slashes and spaces
on the same line as follows:
$$ ,
where ${B}_{n}$, ${x}_{n}$, ${y}_{n}$, ${L}_{n}$, and ${R}_{n}$
are the serial number of the border, the
layer atom of the border bond, the non-layer atom of the border bond,
the atomic symbol of the link atom, and the distance
between the layer atom
and the link atom, respectively. Note that the serial numbers
of the borders are always followed by a colon,
different specifications of the same section are separated by
a comma, and different specifications of the same
border hold the same border serial number in all input sections.
All the border parameters
have to be specified in the case of manual border handling,
and ${B}_{n}$ and $n$ have to be positive integers
smaller than the number of the atoms in the system. Note that
the unit of ${R}_{n}$ has to match that set by keyword unit as
there are no consistency checks in this handling scheme.
Note also that, currently, the manual border handling scheme
is the only option if one considers using the GFN-FF method
via the Mrcc-xtb interface.
The semi-automatic border handling (1) is the combination of the automatic and manual schemes: the borders are selected and set by the automatic approach, but these selections can be modified manually. Any automatically selected ${B}_{m}$ border can be deleted by writing a colon and an $X$ after the serial number of the border in the border atom pair section: $$. New borders can also be defined in the same manner as in the case of the manual border handling, however, all the parameters of the newly defined borders have to be specified. Please be aware that the simple overwriting of single parameters of the automatic approach is not supported by this scheme.
The ONIOM implementation allows arbitrary number of layers, thus the input requirement of an ONIOM calculation with $N$ layers is $(N-1)*3$ lines: the atoms of the second layer are specified in the first line, while atoms of the third layer are specified in the fourth line after the keyword oniom, and so on. To allow full control over the keywords of the subcalculations, the subminp keyword is advised for use. Note that that projection-based embedding and multi-level local correlation are supported as the top-level calculations (see the keywords embed and corembed for the description). Note also that the restart of the ONIOM procedure is controlled by the comprest keyword. See also the examples below.
No ONIOM calculation is performed.
An ONIOM calculation with mechanical embedding will be carried out with this integer specifying the number of layers. The further parameters should be given in the subsequent lines as described above.
An ONIOM calculation with electronic embedding will be carried out with this integer specifying the number of layers. The further parameters should be given in the subsequent lines as described above. See also keywords oniom_eechg and oniom_qcorr for options of the point charge representation of the lower layers.
oniom=off, while mechanical embedding is used if only a positive integer is given as input.
A two-layer ONIOM calculation is requested
where the model system (layer 2) consists of atoms 8 and 9
(using the automatic border handling scheme).
The subcalculations will be executed in the following order:
PBE/STO-3G (full system), PBE/STO-3G (model system),
and CCSD/aug-cc-pVDZ (model system):
basis=aug-cc-pVDZ
calc=CCSD
oniom=2
8-9
PBE / STO-3g
0
The same settings are applied as in the previous case
except that the borders are selected manually:
layer 2 will have two borders between atoms 8-5 and 8-1.
The first (second) link atom is a hydrogen (fluorine),
which will be placed 1.09 Å (1.01 Å) away from atom 8
in the direction of atom 5 (1):
basis=aug-cc-pVDZ
calc=CCSD
oniom=2
8-9
PBE / STO-3g
2 / 1:8-5,2:8-1 / 1:H,2:F / 1:1.09,2:1.01
The same settings are applied as in the previous case
except that the borders are selected semi-automatically:
the first border of the automatic selection is deleted and
an additional third border is specified which is located
between atoms 8 and 2. In this case, chlorine will be used as
link atom, and it will be placed 1.01 Å away from atom 5
in the direction of atom 2:
basis=aug-cc-pVDZ
calc=CCSD
oniom=2
8-9
PBE / STO-3G
1 / 1:X,3:8-2 / 3:Cl / 3:1.01
Four-layer ONIOM calculation,
LNO-CCSD(T):LMP2:PBE:LDA (aug-cc-pVDZ:cc-pVDZ:6-31G*:STO-3G).
Each layer uses automatic border handling scheme and
has a different negative charge, except the top layer,
which has zero charge:
basis=aug-cc-pVDZ
calc=LNO-CCSD(T)
charge=0
mult=1
oniom=4
1-8
LDA / STO-3G / -3 / 1
0
1-5
PBE / 6-31G* / -2 / 1
0
1-2
LMP2 / cc-pVDZ / -1 / 1
0
Two-layer ONIOM calculation,
LNO-CCSD(T)-in-LMP2:PBE (aug-cc-pVTZ:6-31G*).
The automatic border handling scheme is utilized. The multilevel
local correlation approach is used as the top-level calculation
with an automatic selection of active occupied molecular
orbitals.
basis=aug-cc-pVTZ
calc=LNO-CCSD(T)
corembed=on
1-4
LMP2
0
charge=0
mult=1
oniom=2
1-8
PBE / 6-31G* / 0 / 1
0
Two-layer ONIOM calculation,
LNO-CCSD(T)-in-PBE0:PBE (aug-cc-pVTZ:6-31G*).
The automatic border handling scheme is utilized. The projection-based
embedding is used as the top-level calculation
with an automatic selection of active occupied molecular
orbitals.
basis=aug-cc-pVTZ
calc=LNO-CCSD(T)
embed=huzinaga
1-4
PBE0
0
charge=0
mult=1
oniom=2
1-8
PBE / 6-31G* / 0 / 1
0
Two-layer ONIOM calculation,
PBE0:PM6-D3H4 (aug-cc-pVTZ:minimal basis).
The automatic border handling scheme is utilized. The Mopac
semi-empirical program is used for the PM6-D3H4 calculations.
basis=aug-cc-pVTZ
calc=PBE0
charge=0
mult=1
sqmprog=mopac2016
oniom=2
1-8
PM6-D3H4 / / 0 / 1
0
Two-layer ONIOM calculation,
LNO-CCSD(T)-in-LMP2-in-PBE-D3:GFN2-xTB (aug-cc-pVTZ:minimal basis).
The automatic border handling scheme is utilized for
the ONIOM calculation, while the amul scheme
is used for the selection of occupied orbitals in the
embedding and multilevel local correlation
calculations. The xTB
semi-empirical program is used for the GFN2-xTB calculations.
basis=aug-cc-pVTZ
calc=LNO-CCSD(T)
charge=0
mult=1
orbloce=pm
sqmprog=xtb
corembed=on
1-4
0
embed=huzinaga
1-8
PBE-D3
0
oniom=2
1-12
GFN2-xTB / / 0 / 1
0
Specifies the type of atomic charges used in the ONIOM calculations with electronic embedding. Note that in the case of the mulli, lowdin, iao, chelpg , and mk options, the atomic charges are generated after the first ONIOM subcalculation. Note also that in the cases where gopt$\mathrm{\ne}$off and freq$\mathrm{\ne}$off, these charges will be used in all calculations.
No embedding charges. This is equivalent with mechanical embedding.
Mulliken atomic charges are used as embedding charges.
Löwdin atomic charges are used as embedding charges.
Intrinsic atomic orbital (IAO) atomic charges atomic charges are used as embedding charges.
CHELPG atomic charges are used as embedding charges.
Merz–Singh–Kollman atomic charges are used as embedding charges.
The charges are obtained by an external program if the same program is used in the first calculation. Currently, only the xTB program is supported.
Perturbation-independent charges are obtained from the topology file of the Amber program. See the manual of Amber how to generate such files and the example below how to use them with Mrcc. If oniom_eechg=amber and qmmm=off, then it is assumed that the serial numbers of the atoms are the same in the topology file and in MINP.
Those atomic charges will be used that are specified after the line of oniom_eechg. Note that each given atomic charge must be in a new line, the order of the charges must be the same as the order of the corresponding atoms in the specification of the geometry, and the charge of all atoms must be specified.
oniom_eechg=lowdin, if oniom$\mathrm{\ne}$off
and electronic embedding is utilized,
oniom_eechg=off otherwise.
Two-layer ONIOM calculation is requested
with electronic embedding for a system with 9 atoms.
The upper layer consists of atoms 8 and 9,
while the layer borders are handled automatically. The subcalculations
will be executed in the following order: PBE/STO-3G (full system),
PBE/STO-3G (model system, in the presence of the IAO atomic charges of
atoms 1-7), and CCSD/aug-cc-pVDZ (model system, in the presence of
the IAO atomic charges of atoms 1-7).
basis=aug-cc-pVDZ
calc=CCSD
oniom_eechg=iao
oniom=2-ee
8-9
PBE / STO-3G
0
The same settings are applied to a 2-atom system, except
that user-defined charges will be used:
basis=aug-cc-pVDZ
calc=CCSD
charge=-1
oniom_eechg=user
1.00
-1.00
oniom=2-ee
2
PBE / STO-3G
0
A two-layer ONIOM calculation is requested
with electronic embedding, PBE0:GFN2-xTB (6-31G*:minimal basis).
After the full system calculation, a GFN2-xTB
and a PBE0 calculation are executed using the point charges
obtained from the xTB program.
basis=6-31G*
calc=PBE0
oniom_eechg=external
oniom=2-ee
2
GFN2-xTB / / 0 / 1
0
A two-layer ONIOM calculation is requested
with electronic embedding, PBE0:GFN2-xTB (6-31G*:minimal basis).
After the full system calculation, a GFN2-xTB
and a PBE0 calculation are executed using the point charges
obtained from the topology file example.prmtop of the Amber program.
basis=6-31G*
calc=PBE0
oniom_eechg=amber
example.prmtop
oniom=2-ee
2
GFN2-xTB / / 0 / 1
0
Specifies the PCM scheme for ONIOM calculations [101, 149].
ONIOM subcalculations will not use PCM.
PCM is utilized only in the first calculation, i.e., in the lowest-level calculation performed for the entire system.
PCM will be used in all subcalculations. The cavity of the upper layer calculation will be the same as in the first calculation. The wall of the cavity is repolarized by the density of the subcalculation.
If the lowest-level calculation is handed by an external program, the implicit solvent model of the external program can be utilized with oniom_pcm=c.
If the Mopac program is utilized in the first calculation and pcm=$$solvent$\mathrm{>}$, then the eps and rsolv keywords are set for Mopac. The corresponding static permittivity and probe radius values of the given solvent are taken from the PCMSolver library.
oniom_pcm=x if oniom$\mathrm{\ne}$off and pcm$\mathrm{\ne}$off and mechanical embedding is used, oniom_pcm=c if oniom$\mathrm{\ne}$off and pcm$\mathrm{\ne}$off and electronic embedding is used and the lowest-level calculation is handled by an external program, oniom_pcm=off otherwise.
Two-layer ONIOM calculation is requested using the
IEFPCM solvent model utilizing water for the medium.
The upper layer consists of atoms 8 and 9,
while the layer borders are handled automatically. The subcalculations
will be executed in the following order: PBE/STO-3G (full system, in the
cavity of the full system), PBE/STO-3G (model system, in the cavity
of the full system), and CCSD/aug-cc-pVDZ (model system, in the
cavity of the full system). Note that the interface is repolarized in
all subcalculations.
basis=aug-cc-pVDZ
calc=CCSD
pcm=water
oniom_pcm=x
oniom=2
8-9
PBE / STO-3G
0
The same settings are applied as in example 1, except that
PCM is only utilized in the first calculation:
basis=aug-cc-pVDZ
calc=CCSD
pcm=water
oniom_pcm=c
oniom=2
8-9
PBE / STO-3G
0
The same settings are applied as in example 2 except that
the implicit solvent model is only utilized in the first calculation
by an external program:
basis=aug-cc-pVDZ
calc=CCSD
pcm=water
oniom_pcm=c
oniom=2
8-9
PM6-D3H4 / / 0 / 1
0
Specifies the point charge correction scheme for ONIOM calculations to avoid the overpolarization of the electron density by those classical point charges that are saptially close to the quantum mechanical system. This feature also guarantees integer net charges for the subsystems.
The correct integer charge of the subsystems is guaranteed. The point charges of the upper-layer atoms are summed, the integer charge of the upper-layer is subtracted from the sum, and the residual charge is distributed equally among those atoms that were not counted in the summation.
Same as option 0, but moreover, the point charges of the lower-layer host atoms are added to the point charge summation, and then they are zeroed.
Same as option 1, but moreover, user-defined point charges of the lower-layer host atoms are added to the point charge summation, and then they are zeroed. The serial numbers of the additional atoms have to be specified in the subsequent lines after the keyword oniom_qcorr for each layer. See the examples below.
No charge correction scheme is utilized. The initial input point charges will be used throughout the subcalculations, thus the sum of the point charges may differ from the integer charge of the environment subsystem.
In the case of options 0, 1, and special, the integer charges are guaranteed for each sublayer, thus a unique charge correction is utilized for each layer.
Option 1 only differs from option 0 if the subsystems are covalently bound.
oniom_qcorr=off if oniom=off, and oniom_qcorr=1 otherwise.
A CCSD:PBE (aug-cc-pVDZ:STO-3G) calculation is performed for the system that consists of 6 atoms. The manual border handling scheme is employed: a hydrogen link atom is positioned to the layer border, which is between atoms 2 and 5. No charge correction is utilized, thus the input point charges will be employed for the model system calculations. Note that the sum of point charges of the environmental subsystem ($-$1.1) differs from the integer charge of the subsystem ($-$1), that is, the charge of the total system (0) minus the charge of the model system ($+$1)
basis=aug-cc-pVDZ
calc=CCSD
charge=1
oniom_qcorr=off
oniom_eechg=user
0.85
-0.60
0.85
0.20
-1.5
0.20
oniom=2-ee
1-3
PBE / STO-3G / 0 / 1
2 / 1:2-5 / 1:H / 1:1.09
The same settings are applied as in the previous example except that the integer charge of the subsystem is guaranteed with oniom_qcorr=0. In this case, the sum of point charges of atoms 1-3 is summed ($+$1.1), the charge of the model system ($+$1) is subtracted, and the difference (0.1) is distributed equally among the rest of the charges, which results in point charges $+$0.233, $-$1.467, and $+$0.233 for atoms 4, 5, and 6, respectively. Note that the point charge representation of the environmental subsystem matches the integer charge of the subsystem.
The same settings are applied as in the previous example except that oniom_qcorr=1. In addition to the previous example, the point charge of the lower-layer host atom, which is atom 5, is included in the summation. After the integer charge of the subsystem ($+$1) is subtracted from the resulting point charge sum ($-$0.4), the residual charge ($-$1.4) is distributed equally among the rest of the point charges, which results in $-$0.5 and $-$0.5 for atoms 4 and 6, respectively. Note that the point charge representation of the environmental subsystem matches the integer charge of the subsystem.
A CCSD:PBE0:PBE (aug-cc-pVDZ:6-31G*:STO-3G) calculation
is carried out for a system which consists of 12 atoms, and manual border handling scheme is employed. The same charge correction procedure is
utilized as in the previous example, except that additional point charges are used in the summation procedure and then they are zeroed.
In the case of the second layer calculations,
atoms 7 and 8 are added to the charge correction procedure. In the case of the last layer calculations, the point charge of atoms 5-9 will also be used in the summation and then zeroed.
basis=aug-cc-pVDZ
calc=CCSD
charge=1
oniom_qcorr=special
7,8
5-9
oniom_eechg=lowdin
oniom=2-ee
1-6
PBE / STO-3G / 0 / 1
2 / 1:6-9 / 1:H / 1:1.09
1-3
PBE0 / 6-31G / 0 /1
2 / 1:1-4 / 1:H / 1:1.09
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 [111] is the only available option. A geometry optimization can be carried out by either the Broyden–Fletcher–Goldfarb–Shanno (BFGS) or the simplex algorithm.
the simplex method of Nelder and Mead.
the BFGS algorithm.
optalg=simplex for basis set optimization, optalg=BFGS otherwise.
To run a geometry optimization with the simplex algorithm set optalg=simplex
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.
$$any positive integer$\mathrm{>}$
optmaxit=50
to allow 60 iteration steps set optmaxit=60
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.
$$any positive real number$\mathrm{>}$
optetol=1e-6
for a convergence threshold of $5\cdot {10}^{-6}$ E${}_{h}$ set optetol=5e-6
Serial number of the excited state (excluding the ground state) for which the geometry optimization or property calculation is performed. See also keywords nsing, ntrip, and nstate.
$$any positive integer$\mathrm{>}$ that is smaller than the total number of states set by keywords nsing, ntrip, and nstate.
optex=nstate-1, that is, the geometry optimization will be performed for the highest excited state considered.
if you set
calc=CIS
nsing=3
optex=1
for a closed-shell molecule, the ground-state and two excited-state
roots will be determined in each geometry optimization step, but the
gradients will be calculated for the lowest singlet excited state.
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.
$$any positive real number$\mathrm{>}$
optgtol=1e-4
for a convergence threshold of $3\cdot {10}^{-4}$ E${}_{h}$/bohr set optgtol=3e-4
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.
$$any positive real number$\mathrm{>}$
optstol=1e-3 for a basis set optimization, optstol=1e-4 otherwise.
to set a threshold of ${10}^{-5}$ bohr for a geometry optimization type optstol=1e-5
Specifies what type of orbital localization is performed for the core molecular orbitals.
All the options introduced for keyword orbloco also work for orblocc, see the description of keyword orbloco for details.
orblocc=orbloco if localcc=2013, or core=corr and
localcc$\mathrm{\ne}$off;
orblocc=off otherwise
to localize of core orbitals with the Pipek–Mezey algorithm specify orblocc=on
Specifies what type of orbital localization is performed for the molecular orbitals in the case of localization-based embedding calculations.
All the options introduced for keyword orbloco also work for orbloce, see the description of keyword orbloco for details. An additional option is orbloce=special, where one must specify orblocc, orbloco, and orblocv, respectively, in the three subsequent lines after orbloce.
Except for orbloce=special, the same option is set for the core and valance orbitals.
The keyword also controls the orbital localization in the case of automated link atom handling of ONIOM calculations.
orbloce=special is the only option to allow the selection of virtual orbitals in the subsequent projector constructions, where all options of orblocv can be used, except pao-subsys. Note that it is allowed to utilize the full virtual subspace in the high-level SCF run but to use different virtual orbitals during the correlation calculations. See examples 10 and 11 for keyword embed.
If orbloce=spade is utilized, the atoms specified under the keyword embed are used to construct the projector.
orbloce=spade if embed$\mathrm{\ne}$off and
oniom=off,
orbloce=pm if oniom$\mathrm{\ne}$off,
orbloce=off otherwise
to localize the orbitals with the Pipek–Mezey algorithm specify orbloce=pm
to use Pipek–Mezey localization for the occupied orbitals and
SPADE localization for the virtual orbitals, type
orbloce=special
pm
pm
spade
note that this setup allows the selection of active virtual orbitals, thus
the projector of the embedding calculations will contain frozen virtual
orbitals during the high-level SCF run.
Initial guess for the orbital localization.
Guess orbitals are calculated by the Cholesky decomposition of the one-particle density matrix [4].
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.
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.
orblocguess=cholesky
to speed up the orbital localization by using the localized orbitals of a previous calculation as guess specify orblocguess=restart.
Specifies what type of orbital localization is performed for occupied molecular orbitals.
No orbital localization.
Boys localization is performed [27].
Pipek–Mezey localization is performed [126].
intrinsic bond orbitals of Knizia are constructed [76].
generalized Boys localization by Jansík et al. [58], where $m$ is a positive integer defining the power of orbital variance. Note that for $m=1$, the Boys localization is retrieved.
localized orbitals are calculated by the Cholesky decomposition of the one-particle density matrix [4].
the SPADE fragment localization of Claudino and Mayhall [12].
orbloco=off in the general case, orbloco=boys for local correlation calculations
In the case of orbloco=spade, a list of atoms is read after the line orbloco in the following format: $$,$$,…,$$-$$,…, where