Recent changes

Several well parallelized new methods have been implemented in Molpro. A review of recent developments in Molpro can be found in J. Chem. Phys. 152, 144107 (2020).

We recommend always to use the most recent version, since developments are ongoing and problems reported by users are always fixed as quickly as possible. In particular, before reporting bugs, please check if these still occur in the latest version.

A window appears below the brown viewer window, to display the xyz file of the structure presently in the viewer; continuously updated even during vibrations. It can be opened/closed with the menu item Window/Show xyz file

Better placement of the windows on startup.

Various bug fixes and improvements of existing methods.

A new graphical user interface called iMolpro is now included without additional charge. This can be used for generating Molpro inputs for DFT and all single reference methods, and submit these for execution to a chosen machine. This functionality is very much like in gMolpro. In addition, structures, orbitals, and vibrations can be viewed using Jmol.

Unfortunately this has been wrong (usually too small) for pseudo-potential calculations. This bug has been fixed.

Many alias names for functionals computed with libxc have been added. These are compatible with those used in the D4 program.

This program has been very much speeded up, with and without density fitting. Various functionals can now be used in it.

A CAS(2,2) version of icMRCCSD is now available, as published in J. Chem. Phys. 158, 134801 (2023).

New variables EMP2F12_PNO and EMP2F12_OSV are set for MP2-F12 energies computed with PNO or PAO/OSV projector, respectively. This happens in calculations with option projector=PAO or projector=MIXED.

A new automated procedure for the optimisation of quantum nuclei positions has been introduced (adaptive-NEO). The latter works by updating the nuclei position (meaning the nuclear basis centers, together with the respective electronic basis functions) to the nuclear orbital centroids during the SCF cycles. See Adaptive NEO for further information.

The short names [awc]vnz-dk (n=d,t,q,5) now refer to the [aug]-cc-p[wC]V(n+d)Z-DK sets for second-row elements Al-Cl. To avoid using the +d sets (not recommended!) use the full basis set name.

Version 2.3.0 is bundled with Molpro 2024.1 .

There is now a switch “Display/Retina Screen” in the builder and viewer window to deal with those Mac's which have retina screens. It switches from a normal display (retina factor 1) to a display with doubled number of pixels in both directions (retina factor 2). When this factor is changed, then it gets auto-saved in a configuration file $HOME/.pqsmol/gmolpro.conf . If necessary, this configuration file can be edited and retina factors different from 1 or 2 may be used.

LIBGL_ALWAYS_SOFTWARE is set to be true. This avoids a crash on some Linux distributions, which is due to a bug in the Mesa library.

Molpro 2023.2 is a bug-fix release with the same features as Molpro 2023.1

The Linux MPI-PR binary now supports parallel calculations on multiple computer nodes based on UCX (which must be installed on the cluster). Currently, multi-node calculations are restricted to closed-and open-shell DF-HF, DF-KS, and PNO-LCCSD(T). It is generally recommended to use the MPI-PR version for parallel calculations, since it is more stable than the sockets version and does not require GA preallocation (-G option). Note, however, that multi-node calculations can be very slow unless a fast network (e.g. Infiniband) is available.

Unrestricted UCCSD(T) and related methods based on UHF orbitals are now available, see section Open-shell coupled cluster theories.

An improved geometry optimizer using natural internal coordinates is now available and used by default, see section Geometry optimization (OPTG). New special options for this optimizer are described in section Options for the PQS optimizer.

Rotationally invariant DFT grids have been implemented and is used by default. This is achieved by rotating the grid of each atom dependent of its local environment (see section Atomic grid orientation (ORIENT)).

The algorithm for automatically sorting inner and outer core orbitals has been improved and new a new option LOC_OUTCORE=MIX has been added. See Advanced options in section Local correlation methods with pair natural orbitals (PNOs) for details.

The multicomponent restricted NEO Hartree-Fock method has been implemented in Molpro. This allows for the quantum mechanical treatment of a selected number of protons concurrent with the electronic SCF. For more information, please consult the section Nuclear-electron orbital (NEO) method.

The basis set incompleteness error correction based on a DFT model by E. Giner and J. Toulouse (J. Phys. Chem. Lett. 10, 2931 (2019)) has been implemented, see section Basis set incompleteness error correction for details.

The σ-functionals developed by Görling and coworkers (J. Chem. Phys. 154, 014104 (2021), J. Chem. Phys. 155, 134111 (2021)) are implemented in the RIRPA code. See the corresponding section in the manual for details: RIRPA program. σ-functionals are as fast as RPA, but provide more accurate results for many chemical properties.

A new variant of multi-reference perturbation theory, denoted cispt2, has been added. This adds a perturbative double correction to a MRCI-singles (MRCIS) calculation. The initial MRCIS relaxes the reference coefficients before generating the internally contracted doubles. Valence excitation energies are more accurate than with CASPT2 or NEVPT2.

Some problems in NEVPT2 have been fixed, and tighter default integral screening thresholds for DF-NEVPT2 have been set.

Version 2.2.0 of the graphical user interface gmolpro is bundled with Molpro2023.2.0 . When installing gmolpro, then it is installed together with this Molpro version, and uses this version on the local machine.

Due to moving to gtk2.24, the menu bar on Mac is better synchronised.

Important: on some platforms such as openSUSE15.5, export LIBGL_ALWAYS_SOFTWARE=1 or export MESA_LOADER_DRIVER_OVERRIDE=i965 may be necessary before starting gmolpro , to avoid a crash (because of a bug in the Mesa library)

Frequently asked question:

Q: An optimisation (or frequency calculation) is performed, but the icon to open the optimsation (or frequency) window is greyed out, why?

A: The GUI searches for orbitals, and generates a pulldown menu with a set of orbitals found. If there is more than one set of orbitals, then it may be necessary to load a different set of orbitals. If a corresponding optimisation (or frequency) calculation is found, then the icon to open the window will become clickable (and is not greyed out any more).

The –mppx option on the molpro command line has been disabled. mppx (i.e. independent calculations for different geometries on each processor) is used automatically in numerical gradient and hessian calculations, as well as in surface calculations for VCI. The mppx mode can be disabled setting option MPPX=0 on the FREQENCIES, HESSIAN, OPTG, or FORCE command lines, e.g. FREQ,MPPX=0.

It is now possible to run CCSD or CCSD(T) calculations using orbitals from a previous Kohn-Sham (KS) calculation. New options H0=MP|KS and KSFOCK have been added. With H0=MP (default, Moeller-Plesset zeroth-order Hamiltonian) the Fock matrix is recomputed using the KS orbitals and block-diagonalised, while with H0=KS the KS Fock matrix is used as zeroth-order Hamiltonian, and thus the KS orbital energies in the denominators in the perturbative MP2 and triples steps. In both cases single excitations are included in the MP2 step (this is also possible with DF-MP2). Single excitations can be excluded using option NOSINGLES.

Option KSFOCK forces the use of the KS Fock matrix in the subsequent calculation. This is required for double-hybrid DFT and RPA calculations.

Note that the (T) energies are still slightly different in closed-shell CCSD and UCCSD calculations. This is due to missing singles terms in the closed-shell triples program. Work is in progress to fix this problem.

  • The printed LMP2 energies in PNO-LMP2-F12 and PNO-LCCSD-F12 calculations now include the CABS singles correction (if available). The domain-corrected LMP2 energy is stored in a variable EMP2C (alias EMP2_DC), the MP2-F12 energy in variable EMP2F12 (alias EMP2_F12).
  • (T) in PNO-LCCSD(T)-F12 calculations can be restarted from a dump file, see Restarting (T) for details.
  • If SAVE_LMO=record is given the localised orbitals are written to the given record as well as to the xml file for visualisation with gmolpro. If record=1 the current orbital dump record is used.

A bug in the MRCI restart with stored reference vectors has been fixed.

The maximum number of active orbitals in NEVPT2 has been increased from 14 to 32.

An parallelisation problem occurring with large reference spaces has been fixed. The performance and parallelisation of some parts of the gradient program have been improved.

The existing VCI program has been extended by a new directive (ROVIB), which allows for the automated calculation of high-resolution rotational and rovibrational spectra based on RVCI theory. Rotational and Coriolis coupling can be included up to high accuracy and several rotational bases are provided. As the final line lists can be very long, they are stored in an external file, which can be further processed by the DAT2GR program in order to simulate the corresponding spectra. To allow for maximal flexibility with respect to the vibrational bases, the VIBSTATE program has been significantly extended and controls now for all parts of the VSCF/VCI programs the list of states to be computed - including non-Abelian point groups (see the VSCF manual).

NMR shielding tensors can be computed with LDA, GGA and hybrid GGA functionals using density fitting approximations of two-electron integrals using the new dftshield program. The program is restricted to closed-shell systems and currently does not support symmetry.

The documentation of default basis sets in Molpro has been extended, see section Default basis sets

The input specifications for cartesian basis functions have been generalised. It is no possible to specify option CARTESIAN or ATTRIBUTES=CARTESIAN on basis set SET directives, e.g. SET,JK,CONTEXT=JKFIT,CARTESIAN. It should be noted, however, that density fitting with cartesian basis functions is much slower than with spherical ones. F12 calculations with cartesian basis functions are not possible.

The pymolpro Python package has been released, and is available via conda-forge. This supports the same Molpro project bundles as the gmolpro GUI, and allows the programmatic construction, running and analysis of Molpro jobs. For details, see and

Various problems concerning the use of cartesian basis functions have been fixed, including the starting guess and AVAS. Density fitting now works for cartesian basis sets up to h-functions.

The default of DFT grid option orient has been change to 1; this makes KS energies rotationally invariant, but may slightly change energies (mostly in the microhartree range. To recover the old behaviour add grid option orient=1.

A new F12 projector option is added to the PNO-LCCSD program. Different F12 geminal exponents for valence, core-valence, and core-core pairs are now supported. See GEM_BETA for details.

The memory usage and performance of the RS2C expectation value calculations are significantly improved. The NOPROP option is added to skip the expectation value calculations.

DFT options for controlling grid accuracy have been added to gmolpro.

OPTG options for choosing optimization methods and coordinates have been added to gmolpro.

Version 2.1.0 of the graphical user interface gmolpro is bundled with Molpro2022.2.3 . When installing gmolpro, then it is installed together with Molpro version 2022.2.3 , and uses this version on the local machine.

The token may be installed in $HOME/.molpro/token . Alternatively (useful for group licences), on Linux, it may be installed in /usr/share/gmolpro/molpro/lib/.token , and on Mac in /Applications/

A slider to change the vibrational amplitudes has been added. As before, the animation speed can be changed with Options/Graphics/Animation Speed

The xyz-geometry in the builder window is now displayed in a little extra window.

Various bug fixes: z-matrix containing dummy atoms now properly written by Molpro, so that gmolpro correctly displays it; small window (segment window) in builder: hydrogens + dummies properly displayed etc.

Force-field parameters autodetected when importing xyz.

Version 2.0.0 of the graphical user interface gmolpro is bundled with Molpro2022.2.2 . When installing gmolpro, then it is installed together with Molpro version 2022.2.2 , and uses this version on the local machine.

The token may be installed in $HOME/.molpro/token . Alternatively (useful for group licences), on Linux, it may be installed in /usr/share/gmolpro/molpro/lib/.token , and on Mac in /Applications/

A bug affecting KS and TDDFT calculations with cartesian basis functions has been fixed.

The perturbation-adapted zero-order hamiltonian described in J. Chem. Phys. 156, 011101 (2022); is available for closed-shell systems.

With CORE,MIXED correlation of the 2s and 2p electrons for the second-row elements Al - Ar is now excluded. Global CORE command is now supported. See defining orbital subspaces for details. The old behaviour, where these electrons were correlated, can be recovered using CORE,MIXED_OLD.

  • Quadratic optimization or mixed quadratic/first-order SO-SCI optimization is available for closed-shell and open-shell RHF by using {RHF,SO;…} or {RHF,SO-SCI;…}.
  • The CAHF optimization with the quadratic or the SO-SCI optimization is implemented into the CAHF program, and it can be called by {CAHF,SO;…} or {CAHF,SO-SCI;…}.
  • A new variant of the SO-SCI optimization is available for the RHF case. It optimizes all AVAS orbitals quadratically and is called by {RHF,SO-SCI-ACTIVE}.

The SO-SCI option significantly improves convergence and is recommended in difficult cases. The cost may even be lower than for standard ROHF or CAHF calculations, since the number of iterations is reduced. Furthermore, the probability of convergence to local minima or saddle points is much reduced.

See section options for the SO-SCI optimization for more details and options.

  • The entire export of CI vectors to subsequent Molpro programs has been rewritten. CI vectors can be exported via SAVE,CIREC=record or NATORB,orbrecord,CIREC=record; more details can be found here. All exported CI vectors are stored in the CSFs basis. Determinants are automatically transformed into the CSFs basis, which also allows an export with state-averaging over different spins. Warning: MULTI based CI records from earlier Molpro versions are not compatible anymore, and the wfu file has to be recalculated with Molpro 2022!
  • {MULTI,CASCI;…} can now be simply called by {CASCI;…}, and the input orbitals are not transformed after the CASCI optimization (also for dont,orbital).
  • Gradients of the state-averaged energy are now available.
  • The program is now more robust against spin-contamination in determinant based calculations.

The CI vectors stored in multi using the directive SAVE,CIREC=record can now be used in the MRCI and CASPT2 programs as reference vectors using option CIREC=record.

A new much more efficient algorithm for computing the spin-orbit Hamiltonian using MCSCF wavefunctions and ECP spin-orbit operators. The program to compute and diagonalize spin-orbit matrices can generally be called using

HLSMAT,type, record1, record2, record3, …

where record1, record2 etc are records saved with MULTI or MRCI. The new program is used if type is ECPLS and if no external excitations are present (so either MCSCF/CASSCF wavefunctions or MRCI ones with NOEXC).

For more details see section calculation and diagonalization of the entire SO-matrix.

The IBO localization convergence threshold has been tightened from $10^{-8}$ to $10^{-9}$. The PAO redundancy threshold THRLOC has been tightened from $10^{-6}$ to $10^{-7}$.

The localisation of core orbitals has been improved; it now uses an AO-based algorithm, which minimises mixings of orbitals of different types (e.g. s, p$_x$, p$_y$, p$_z$). It has been found that this stabilises the F12(PNO) energy contribution if core orbitals are correlated. The new default for core-orbital localization is loc_method_core=IBO(AO).

The previous behaviour of the PNO program can be recovered by putting before the first PNO command:

local, thrpip=1.d-8, thrloc=1.d-6, loc_method_core=ibo

The corresponding new options are

local, thrpip=1.d-9, thrloc=1.d-7, loc_method_core=ibo(ao)

The options can also be given on the PNO command, lines.

Version 1.4.0 of the graphical user interface gmolpro is compatible with Molpro2022.1 .

Main improvement over earlier versions is the ability to handle output files (xml files) with several optimisations, frequency calculations, or even different molecules within one xml file. Start from the brown window, open the orbital window, and then select and load an orbital set. A corresponding optimisation or vibrational modes may now be displayed.

When creating an input in guided mode, more properties may be selected by buttons.

An option has been added to use MPI files instead of large GlobalArrays to avoid potential GA allocation problems. From Molpro2020.2 the disk option (--ga-impl disk) is the default for single-node calculations. Thus, GA preallocation is not necessary unless --ga-impl ga is used. See disk option for details.

For molpro2020.2 and Linux systems, an additional binary compiled with the GA mpi-pr option is provided. It is recommended to use this binary if --ga-impl ga is chosen, since it avoids potential GA allocation problems as well (the -G option is not needed). Note, however, that this requires an extra helper process. For example, if -n 20 is given, only 19 processes will be used for program execution. See Running Molpro on parallel computers for details.

Multi-state CASPT2 is now available in the rs2c program, see Multi-State CASPT2.

NEVPT2 can now be done with density fitting, the command is DF-NEVPT2

Additional options allow a more flexible definition of the target space.

A new implementation of the HF/KS program is now used by default. In difficult cases this shows more robust convergence.

  • The old version can still be used by specifying option “old”, e.g. hf,old.
  • Optionally, the hybrid second-order/super-ci method (so-sci) can be used in difficult cases (hf,so-sci). In open-shell cases this converges more robustly than the standard ROHF procedure.
  • CAHF optimization with the multi program for the state-averaged optimization including all irreducible representations and spin manifolds. Convergence is more robust than in the equivalent Hartree-Fock optimization, especially for multiple CAHF shells.
  • Export of the raw orbitals used in the MCSCF optimization. Input orbitals stay untouched in combination with the command dont,orbital.

Analytical energy gradients are now available for RASPT2 (using rs2).

  • Rotatory strengths for simulating CD spectra of optically active molecules can be computed for spin-restricted TDDFT excitations
  • First-order nonadiabatic coupling matrix elements (NACME) are available for closed and open-shell TDHF and TDDFT (LDA or GGA-type functionals)
  • All fundamental constants used in the program have been updated to the most recent CODATA 2018 standard
  • The sparse orbital print has been modified. Now it shows for each printed coefficient the center number, the function type, and the sequence number for this type at the given atom (denoted MU).
  • DF-MCSCF (single state), DF-HF, DF-UHF, DF-KS, and DF-UKS gradients are now available with symmetry.
  • TD-DFT gradients for closed-shell singlet excited states are available with and without density fitting.
  • The core orbitals are no longer excluded from the CCSD-F12 projector by default. The approximation, which is applied by default in previous versions, saves a few percent of the computation time, but can have a significant effect in calculations involving transition metals. The old default can be recovered by MODOMC=1.
  • Correlated core orbitals are now localized separately from the valence orbitals by default, which mostly improves the results and convergence in calculations with core orbital correlation. The old default can be recovered by LOC_OUTCORE=ON.
  • The ENERGY variable is now set to a scalar representing the final energy in F12 calculations. PNO-LCCSD-F12 is now equivalent to PNO-LCCSD-F12b, and F12a and F12b total energies can be accessed from the variables ENERGYA and ENERGYB. Similar changes are made to the ENERGC variable.
  • The scaling factor for F12-scaled triples now excludes the contribution from distant pairs, for which F12 energies are neglected.

The minimum requirements for using the Linux binaries has increased to glibc 2.17 or newer.

A windows beta version of Molpro is now available.

A powerful graphical user interface (based on PQSMol) for building and preoptimising molecular structures, preparing and running Molpro inputs, and visualisation of results (will soon be available). gMolpro also supports remote job submission and is available for Linux and Mac.

The CASSCF program multi has been largely rewritten and well parallelized. Improved second-order algorithms as well as first-order methods applicable to large molecules have been implemented. The new program combines robust convergence with excellent efficiency. Details are described in J. Chem. Phys. 150, 194106 (2019) and J. Chem. Phys. 152, 074102 (2020).

The PNO-CASPT2 program for large molecules has been extended to allow multi-state MS-CASPT2 calculations. For details see J. Chem. Phys. 150, 214107 (2019).

The time-dependent DFT program has been completely rewritten to support molecular symmetry, open-shell systems (for spin-unrestricted wave functions), various integral modes (including a very fast parallelised density-fitting mode) and standard LDA, GGA, hybrid-GGA and range-separated hybrid GGA functionals and kernels. Calculations using the exact Kohn-Sham exchange (TDEXX) method can be done both by using the adiabatic and non-adiabatic EXX kernel. Linear response properties can be calculated for any one-electron operators available in Molpro. Isotropic and anisotropic $C_6$, $C_8$ and $C_{10}$ dispersion coefficients can be computed along with the calculation of frequency dependent (dipole,quadrupole,otocpole) polarisabilities.

An interface to the external LibXC library (see of density functionals has been added. Therefore virtually any density functional of the LDA, GGA, hybrid-GGA, meta-GGA and range-separated hybrid-GGA can be used in Molpro calculations in the DFT and Kohn-Sham programs. Range-separated DFT calculations with both long-range and/or short-range exact exchange admixtures can be done now also using an efficient density-fitting implementation of the Fock matrix computation.

An implementation of the nonlocal van-der-Waals DF VV10 by Vydrov & Van Voorhis is available for spin-restriced wave functions to improve the description of long-range correlation interactions of standard DFT methods.

The local coupled-cluster program PNO-LCCSD(T)-F12 has been revised to support open-shell molecules. The new revision also includes some minor changes in the local approximations and is also used for closed-shell molecules by default. The closed-shell program in Molpro 2019.2 and earlier can be executed with the option version=2019.2. For more details see Local correlation methods with pair natural orbitals (PNOs).

A new program, XSURF, has been implemented for the efficient calculation of potential energy surfaces (PES) as needed for the calculation of anharmonic spectra. In contrast to the “old” SURF program, XSURF can handle any type of symmetry, is not restricted to the length of the expansion of the PES and is restartable at any point.

Due to a complete redesign of the programs for the fitting (POLY) and the transformation of potential energy surfaces (PESTRANS), these programs could be accelerated by several orders of magnitude (only operative in combination of the new XSURF program).

Similar to the Franck-Condon program, the new EVSPEC program allows for the calculation of anharmonic electronic-vibrational absorption spectra with the inclusion of Duschinsky effects. Two different imlementations are provided, the contracted invariant Krylov subspace (CIKS) approach and a Raman wavefunction (RWF) formalism.

A new parallel density fitted Fock matrix generation has been implemented. For large molecules this scales well across computer nodes (with a fast network such as Infiniband) up to about half as many computing cores as occupied orbitals.

Symmetry can now be used in DF-HF and DF-KS (restricted and unrestricted), but not yet for gradients.

Improved efficient implementation of quasi-variational coupled-cluster (QVCCD).

Analytical energy gradients for closed shell DF-MP2-F12, DF-CCSD(T)-F12, CCSD(T) and open-shell [DF-]RMP2 have been added. The gradients for explicitly correlated (F12) methods are restricted to certain Ansätze; please refer to the users manual for details and W. Györffy, G. Knizia and H.-J. Werner, J. Chem. Phys. 147, 214101 (2017); W. Györffy and H.-J. Werner, J. Chem. Phys. 148, 114104 (2018).

Analytical nuclear gradients have been implemented for projection-based wavefunction-in-DFT (WF-in-DFT) embedding with and without atomic orbital (AO) truncation. The current available methods that can be used for the WF method on subsystem A are CCSD(T), CCSD, MP2, and HF. The current available methods that can be used for the low-level SCF method are LDA, LDAX (LDA with any amount of exact exchange) and HF. The support for using GGAs as the low-level method will be coming soon.

Analytic energy gradients for local MP2, RMP2, CC2 and ADC(2) (with density fitting) are available. These methods can be used for geometry optimisations and property calculations of larger molecules.

A hierarchy of local coupled cluster models for ionization potentials employing LMP2 or LCC2 ground state amplitudes and the Jacobian formally reaching the IP-CCSD level: IP-CCSD$_{\rm CC2}$ is now available. For details see G. Wälz, D. Usvyat, T. Korona, M. Schütz, J. Chem. Phys. 144, 084117 (2016).

Local MP2 NMR shielding, magnetisability, and rotational g tensors.

The internally contracted multi-reference program of A. Köhn et al. has been interfaced to Molpro. An embedded Multireference Coupled Cluster method as described in D. J. Coughtrie, R. Giereth, D. Kats, H.-J. Werner, and A. Köhn, J. Chem. Theory Comput. 14, 693 (2018) is also available.

A new parallel implementation of Hartree Fock with local density fitting (closed and open-shell), as described in C. Köppl, H.-J. Werner, J. Chem. Theory Comput. 12, 3122 (2016). This program can run in parallel across several nodes and is for large systems significantly faster than canonical DF-HF, DF-KS.

An accurate and efficient local PNO-CASPT2 for large molecules as described in F. Menezes, D. Kats, H.-J. Werner, J. Chem. Phys. 145, 124115 (2016). For extended systems the computational effort scales linearly with molecular size. A multi-state version will be made available soon.

Explicitly correlated well parallelized PNO methods are now available up to the PNO-LCCSD(T)-F12 level. This program yields results that are very close to the corresponding canonical CCSD(T)-F12 ones. For medium size cases where canonical calculations are still feasible, the PNO-LCCSD(T)-F12 is up to an order of magnitude faster than the canonical CCSD(T)-F12 program, while relative energies typically differ by only 0.2 kcal/mol. The largest applications so far include molecules with up to about 300 atoms and 10000 basis functions. The methods are described in H.-J. Werner, J. Chem. Phys. 145, 201101 (2016). M. Schwilk, Q. Ma, C. Köppl, H.-J. Werner, J. Chem. Theory Comput. 13, 3650 (2017); Q. Ma, M. Schwilk, C. Köppl, H.-J. Werner, J. Chem. Theory Comput. 13, 4871 (2017); Q. Ma and H.-J. Werner, J. Chem. Theory Comput. 14, 198 (2018). A review can be found in Q. Ma and H.-J. Werner, WIREs Comput. Mol. Sci. 2018;e1371. References for PNO-LMP2 and PNO-LMP2-F12 are given under “New features for Molpro2015.1”.

The new D4 dispersion model improves the description of atomic dispersion coefficients through charge-dependent local polarisabilities obtained by a self-consistent tight-binding method, see E. Caldeweyher, C. Bannwarth, and S. Grimme J. Chem. Phys. 147, 034112 (2017)

The nonlocal DFT method (NLDFT) improves the description of long-range correlation interactions of standard DFT functionals through a double-Hirshfeld partitioning of the correlation energy density, see A. Heßelmann, J. Chem. Theory Comput. 9, 273 (2013).

Random-phase approximation electron correlation methods based on Kohn-Sham reference determinants, see Refs. A. Heßelmann, Phys. Rev. A 85, 012517 (2012); A. Heßelmann and A. Görling, J. Chem. Theory Comput. 9, 4382 (2013).

Several new features have been implemented in the DFT-SAPT program, including:

  • Exchange interaction energy contributions without the single-exchange approximation (R. Schäffer and G. Jansen, Theor. Chem. Acc. 131 (2012) 1235; Mol. Phys. 111, 2570 (2013))
  • DFT-SAPT employing exact-exchange response kernels (A. Heßelmann, J. Chem. Theory Comput. 14, 1943 (2018))
  • Regularised SAPT for estimating short-range polarisation interactions (A. J. Misquitta, J. Chem. Theory Comput. 9, 5313 (2013))

This method helps to find good starting orbitals for CASSCF calculations. It is based on the work of Knizia et al., J. Chem. Theory Comput. 13, 4063 (2017).

  • Improved parallelization of coupled-cluster codes
  • Improved DFT quadrature
  • New plugin interfaces for other programs, e.g. MR-CCSD(T), DMRG and FCIQMC codes, supporting directed parallel execution
  • Implementation of the eXact-2-Component (X2C) scalar relativistic Hamilitonian
  • Support for a Gaussian finite nucleus model
  • New correlation consistent basis sets for heavy alkali and alkaline earth elements (both relativistic all-electron and with ECPs), as well as all-electron (DK3 and X2C) basis sets for the lanthanide and actinide elements.
  • X2C-contracted versions for all the standard correlation consistent basis sets for H, He, B-Ne, and Al-Cl, e.g., cc-pVDZ-X2C

Molpro 2015 contains many improvements, bug-fixes and new features. Only the most important ones are listed in the following.

The WF-in-DFT impementation in Molpro [F. R. Manby et al., J. Chem. Theory Comput. 8, 2564 (2012)] permits embedding of almost any ground-state wavefunction method in Molpro in an environment described by DFT.

IBOs provide a reliable and efficient way to generate localized orbitals and to analyze wavefunctions, see G. Knizia, J. Chem. Theory Comput. 9, 4834 (2013)

Molpro 2015 supports NBO6 Natural Bond Orbital Analysis (NBO6) via an interface. NBO6 must be licensed separately.

Usually close (and weak) pairs are treated only at the LMP2 level in LCCSD(T) calculations. This not satisfying since LMP2 often treats van der Waals interactions poorly, hence compromising the high accuracy of the actual LCCSD(T) calculation. Fortunately, close pairs can now be treated in an inexpensive way at a level beyond LMP2, ranging from local direct RPA via ring CCD to LCCD[S]-R$^{-6}$. This is described in

O. Masur, D. Usvyat and M. Schütz, J. Chem. Phys. 139, 164116 (2013);
M. Schütz, O. Masur and D. Usvyat, J. Chem. Phys. 140, 244107 (2014).

This approximation yields improved results at the cost of CCSD. In particular, DCSD equilibrium structures are as accurate as for CCSD(T), see

D. Kats and F. R. Manby, J. Chem. Phys. 139, 021102 (2013);
D. Kats, J. Chem. Phys. 141i, 061101 (2014);
D. Kats, D. Kreplin, H.-J. Werner, F. R. Manby, J. Chem. Phys. 142, 064111 (2015).

Analytical gradients are now available for DCSD. Analytical gradients for CCSD(T), as well as for MP2-F12, CCSD(T)-F12, and DCSD-F12 are under development and will be made available free of charge at a later stage.

Orbital relaxed properties, and analytical nuclear gradients for excited states are now available for local CC2 response and local ADC(2). This allows for geometry optimizations of excited states of large chromophores at these levels of theory. Furthermore, also local ADC(2) transition moments at the level of strict ADC(2) (which require second-order ground state singles and doubles) are available. For further reference see

K. Ledermüller, D.Kats and M. Schütz, J. Chem. Phys. 139, 084111 (2013);
K. Ledermüller and M. Schütz, J. Chem. Phys. 140, 164113 (2014);
M. Schütz, J. Chem. Phys. 142, 214103 (2015).

New parallel linear scaling PNO-LMP2 and PNO-LMP2-F12 programs as described in

H.-J. Werner, G. Knizia, C. Krause, and M. Schwilk, J. Chem. Theory Comput. 11, 484 (2015);
Q. Ma and H.-J. Werner, J. Chem. Theory Comput., DOI: 10.1021/acs.jctc.5b00843 (2015);
C. Köppl and H.-J. Werner, J. Chem. Phys. 142, 164108 (2015).

are available in version 2015. PNO-LCCSD(T) is currently under development and will be made freely available to licensees of Molpro 2015 at a later stage.

The program for NMR chemical shifts has been extended to include also magnetizability and rotational g-tensor. This is described in

S. Loibl and M. Schütz, J. Chem. Phys. 141, 024108 (2014).

DFT and F12 methods with density fitting are significantly speeded up using improved integration routines written by G. Knizia.

A VPT2 program based on a semi-quartic force field (QFF) has been implemented, which allows for low-cost evaluations of anharmonic vibrational frequencies and associated vibrational constants. The determination of the QFF does not require analytical gradients (see the SURF program) and can thus been used in combination with any electronic structure method.

A suite of vibrational multi-reference methods is available now. Modals can be optimized at the VMCSCF level [S. Heislbetz and G. Rauhut, J. Chem. Phys. 132, 124102 (2012)] and correlation effects can be accounted for either variationally (VMRCI [F. Pfeiffer and G. Rauhut, J. Chem. Phys. 140, 064110 (2014)].

A Franck-Condon program based on anharmonic vibrational wavefunctions has been implemented. Franck-Condon factors can either be computed by rotating the vibrational wavefunction or by transforming the potential energy surface in order to account for Duschinsky effects. This program, which allwos for the accurate calculation of photoelectron spectra (absorption and fluorescence) relies on the newly developed transformation program SURFTRANS [P. Meier and G. Rauhut, J. Chem. Phys. (2015)].

Multi-dimensional potential energy surfaces spanned in terms of normal coordinates - as computed with the SURF program - can now be rotated or transformed into an arbitrary set of normal coordinates. This allows for the efficient calculation of Franck-Condon factors or the calculation of vibrational spectra of isotopologues from just one potential energy surface [P. Meier, D. Oschetzki, R. Berger and G. Rauhut, J. Chem. Phys. 140, 184111 (2014)].

Within the SURF program, multi-dimensional potential energy surfaces can now alternatively be spanned in terms of localized normal coordinates. Three localization schemes are offered: (a) localization of the CH-stretching modes only and (b) localization of all normal coordinates to the subunits in molecular cluster calculations or © both. All localization schemes are symmetry sensitive.

In the old Molpro release the calculation of vibrational states was limited to fundamentals, first vibrational overtones and simple combination bands. The new VIBSTATE program allows to specify arbitrary vibrational states, but the user must take care, that the parameters in the SURF program and the vibrational correlation programs are adjusted accordingly.

Raman scattering activities can now be computed within all vibrational SCF and vibration correlation programs. As polarizabilities can only be determined at the Hartree-Fock level, the accuracy is currently still limited.

A facility is provided to store and interrogate sets of molecules, together with information about how they are to be combined in balanced chemical equations. This database can be generated manually, or partially by running appropriate Molpro calculations. Analysis of the database can give a summary of the energy changes associated with each described reaction, and two or more similar databases can be compared reaction by reaction, to give a statistical analysis of the differences between them. Several common databases are included.

A new DFT-SAPT program has been implemented in Molpro which can be used in conjunction with monomer-centered (MCBS), monomer-centered plus (MC+BS) or dimer-centered (DCBS) basis sets. Even in DCBS mode fairly large complexes with about 800 electrons can be studied with this program, see A. Heßelmann and T. Korona. J. Chem. Phys. 141, 094107 (2014).

A summary of the features in Molpro can be found in a recent review: H.-J. Werner, P. J. Knowles, G. Knizia, F. R. Manby, and M. Schütz, Molpro – a general purpose quantum chemistry program package, WIRES Comput. Mol. Sci. 2, 242 (2012), 10.1002/wcms.82.

The new features of Molpro version 2012.1 include the following.

This modification of standard CCSD is capable of robustly describing chemical bond breaking with a single Hartree-Fock reference determinant (see J. B. Robinson and P. J. Knowles, J. Chem. Phys. 136, 054114 (2012), 10.1063/1.3680560). It is implemented for closed-shell systems, with either the Brueckner condition or energy optimisation for the determination of orbitals. Triple excitations can be included perturbatively, BQVCCD(T) or OQVCCD(T) (J. B. Robinson and P. J. Knowles, Phys. Chem. Chem. Phys. 14, 6729-6732 (2012), 10.1039/C2CP40698E; J. Chem. Theor. Comput. 8, 2653-2660 (2012), 10.1021/ct300416b).

A new internally contracted MRCI code [see K. R. Shamasundar, G. Knizia, and H.-J. Werner, A new internally contracted multi-reference configuration interaction (MRCI) method, J. Chem. Phys. 135, 054101 (2011)] is now available. As compared to the old MRCI code it has the following advantages: (i) inactive orbitals (correlated in MRCI but closed-shell in the reference function) are treated explicitly, i.e., no density matrices and coupling coefficients need to be computed that involve these orbitals. Thus, in principle any number of inactive orbitals can be correlated, without the previously existing limitation to 32 correlated orbitals. Furthermore, additional configuration spaces are internally contracted (as in the RS2C code), resulting in a much improved efficiency, particularly for cases with many inactive orbitals. Currently, the method is implemented only for single-state calculations, but an extension to multi-state cases is under development and will be provided in the near future.

Explicitly correlated multireference theories (CASPT2-F12, MRCI-F12) as described in T. Shiozaki and H.-J. Werner, Communication: Second-order multireference perturbation theory with explicit correlation: CASPT2-F12, J. Chem. Phys. 133, 141103 (2010); T. Shiozaki, G. Knizia, and H.-J. Werner, Explicitly correlated multireference configuration interaction: MRCI-F12, J. Chem. Phys. 134, 034113 (2011); T. Shiozaki and H.-J. Werner, Explicitly correlated multireference configuration interaction with multiple reference functions: Avoided crossings and conical intersections, J. Chem. Phys. 134, 184104 (2011) lead to much improved convergence of the correlation energies with the basis set size. Explicitly correlated version of the RS2C and MRCIC codes are not yet available but will be made available as soon as possible.

A review of the explicitly correlated multireference methods can be found in T. Shiozaki and H.-J. Werner, Mol. Phys. 111, 607 (2013).

Extended multistate CASPT2 (XMS-CASPT2) [see T. Shiozaki, W. Győrffy, P. Celani, and H.-J. Werner, Communication: The extended multi-state CASPT2 method: Energy and nuclear gradients, J. Chem. Phys. 135, 081106 (2011)] provides a better description of near degenerate situations and avoided crossings. Currently this option is available only with the RS2 code, not for RS2C.

CASSCF and CASPT2 as well as the corresponding analytical gradient theories are now available with density fitting (DF-CASSCF, DF-RS2). See details in W. Győrffy, T. Shiozaki, G. Knizia, and H.-J. Werner, Analytical energy gradients for second-order multireference perturbation theory using density fitting, J. Chem. Phys. 138, 104104 (2013).

The (F12*) approximation proposed by Hättig et al. J. Chem. Phys. 132, 231102 (2010) has been implemented for closed-shell cases (in Molpro, this is denoted F12c). Other recent work (partly already included in Molpro2010.1) is described in H.-J. Werner, G. Knizia, and F. R. Manby, Explicitly correlated coupled-cluster methods with pair-specific geminals, Mol. Phys. 109, 407 (2011); K. A. Peterson, C. Krause, H. Stoll, J. G. Hill, and H.-J. Werner, Application of explicitly correlated coupled-cluster methods to molecules containing post-3$d$ main group elements, Mol. Phys. 109, 2607 (2011) for alternative approximations.

Further basis sets of K. A. Peterson and J.G. Hill for explicitly correlated methods have been included. In particular these include the aug-cc-pVnZ-PP/OptRI and aug-cc-pwCVnZ-PP/OptRI sets for the group 11 and 12 transition metals. A full set of F12 basis sets for the $p$-block elements Ga-Rn will be added in the very near future. For recent work on molecules containing transition metals see D. H. Bross, J. G. Hill, H.-J. Werner, and Kirk A. Peterson, Explicitly correlated composite thermochemistry of transition metal species, J. Chem. Phys. 139, 094302 (2013).

The local coupled cluster methods have been further improved. See H.-J. Werner and M. Schütz, An efficient local coupled-cluster method for accurate thermochemistry of large systems, J. Chem. Phys. 135, 144116 (2011) for a description of the current implementation. An open-shell implementation [DF-LUCCSD(T)] is now also available (Y. Liu and H.-J. Werner, to be published.) Furthermore, the direct random phase approximation (RPA) has been implemented as a local method and using density fitting. Direct RPA can be considered as a CCD method reduced to only ring diagrams. In contrast tp MP2 it also contains higher-order diagrams which e.g. cover Axilrod-Teller terms in intermolecular calculations. It is therefore especially attractive to replace in LCC calculations the LMP2 method for weak pairs by LRPA. This is now possible in our LCC code (O. Masur and M. Schütz, to be published.)

Local coupled cluster methods can optionally use orbital specific virtual orbitals (OSVs), see J. Yang, G. K. L. Chan, F. R. Manby, M. Schütz, and H.-J. Werner, The orbital-specific virtual local coupled-cluster singles and doubles method: OSV-LCCSD, J. Chem. Phys. 136, 144105 (2012). The main advantage of this method is that the accuracy of the domain approximation can be controlled by a single parameter. Further developments that also use pair-natural orbitals (PNOs) as descibed in C. Krause and H.-J. Werner, Perspective: Comparison of explicitly correlated local coupled-cluster methods with various choices of virtual orbitals, Phys. Chem. Chem. Phys. 14, 7591-7604 (2012) are in progress. A preliminary parallel PNO-LMP2-F12 method is already available.

The explicitly correlated coupled cluster methods as described in H.-J. Werner, Eliminating the domain error in local explicitly correlated second-order Møller-Plesset perturbation theory, J. Chem. Phys. 129, 101103 (2008); T. B. Adler, F. R. Manby, and H.-J. Werner, Local explicitly correlated second-order perturbation theory for the accurate treatment of large molecules, J. Chem. Phys. 130, 054106 (2009); T. B. Adler and H.-J. Werner, Local explicitly correlated coupled-cluster methods: Efficient removal of the basis set incompleteness and domain errors, J. Chem. Phys. 130, 241101 (2009); T. B. Adler and H.-J. Werner, An explicitly correlated local coupled-cluster method for calculations of large molecules close to the basis set limit, J. Chem. Phys. 135, 144117 (2011) are available for closed-shell cases (an open-shell implementation will follow soon). In these methods the errors of the local domain approximation are eliminated to a large extend, and the same accuracy as with the corresponding canonical methods can be achieved, even for molecules with 50-100 atoms.

The efficiency of density functional theory has been much improved. In particular, the density fitting (DF-RKS, DF-UKS) codes for analytical energy gradients are now very much faster, due to a new integral code (adaptive integral core, AIC) written by G. Knizia. Some benchmarks can be found in H.-J. Werner, P. J. Knowles, G. Knizia, F. R. Manby, and M. Schütz, Molpro – a general purpose quantum chemistry program package, WIRES Comput. Mol. Sci. 2, 242 (2012).

A large number of additional density functionals has been added, including PBE0, PBEREV, M05, M05-2X, M06, M06-2X, M06-L, M06-HF, M08-HX, M08-SO, M11-L, SOGGA, SOGGA11, SOGGA11-X. (M11 is currently not available, but will likely be added in the near future).

New analytical gradient codes are now available for DF-MP2, DF-CASSCF, DF-RS2 (including MS and XMS options, also without density fitting), and for CCSD.

First-order properties of excited states via time-dependent Local CC2 linear response theory and ADC(2) are now also available for triplet states. This is described in K. Freundorfer, D. Kats, T. Korona, and M. Schütz, Local CC2 response method for triplet states based on Laplace transform: Excitation energies and first-order properties, J. Chem. Phys. 133, 244110 (2010);

An efficient method for calculating NMR shielding tensors at the local MP2 level has been implemented. Gauge including atomic orbitals (GIAOs) are used to eliminate the gauge origin dependence. Density fitting is employed to factorize the relevant electron repulsion integrals and their derivatives w.r. to the magnetic field. So far, the method has been already applied to systems with more than 2500 contracted basis functions and 300 correlated electrons. Relevant publications are S. Loibl, F. R. Manby, and M. Schütz, Density fitted, local Hartree-Fock treatment of NMR chemical shifts using London atomic orbitals, Mol. Phys. 108, 477 (2010), and S. Loibl and M. Schütz, NMR shielding tensors for density fitted local second-order Møller-Plesset perturbation theory using gauge including atomic orbitals, J. Chem. Phys. 137, 084107 (2012).

Symmetry-adapted perturbation theory of intermolecular interactions with monomers described in the CCSD level.

For open-shell systems, RHF and RKS now use a two-step diagonalization process by default: Here the beta orbitals are found by a second diagonalization in the subspace of occupied alpha orbitals. This process usually leads to more stable convergence in difficult cases, compared to the standard diagonalization of a single open-shell Fock matrix (the latter behavior is recovered by {rhf,algo=0}). An additional section with suggestions for dealing with difficult cases has been added to be manual under “The SCF program”. Additionally, various limitations of the SCF program have been lifted. In particular, {rhf; maxit,1} now always works, and just calculates a Fock matrix and energy from the input orbitals, without updating said orbitals.

The FCIQMC program exists through an interface to the NECI codebase, which is actively developed in the group of A. Alavi, and has been integrated in to the Molpro code. The FCIQMC method is a recently introduced stochastic method which can calculate in principle FCI-quality energies for small to medium-sized molecules. See G. H. Booth, A. J. W. Thom, and A. Alavi, J. Chem. Phys. 131, 054106 (2009); D. M. Cleland, G. H. Booth, and A. Alavi, J. Chem. Phys. 134, 024112 (2011); G. H. Booth, D. M. Cleland, A. J. W. Thom, and A. Alavi, J. Chem. Phys. 135, 084104 (2011).

The AIMS module implements the Ab Initio Multiple Spawning method to perform dynamics calculations on multiple electronic states. It can also be used quite generally for first principles molecular dynamics on a single electronic surface, provided that nuclear gradients are available. Currently, non-adiabatic dynamics is limited to CASSCF wavefunctions; however, MS-CASPT2 non-adiabatic dynamics (with an implementation of analytical MS-CASPT2 non-adiabatic couplings) will be provided in the very near future. See M. Ben-Nun and T. J. Martinez, Chem. Phys. Lett. 298, 57 (1998); B. G. Levine, J. D. Coe, A. M. Virshup and T. J. Martinez, Chem. Phys. 347, 3 (2008); T. Mori, W. J. Glover, M. S. Schuurman and T. J. Martinez, J. Phys. Chem. A 116, 2808 (2012).

The partially augmented Turbomole basis sets def2-SVPD, def2-TZVPD, and def2-QZVPPD (Rappoport, Furche: J. Chem. Phys. 133, 134105 (2010)) have been added to the library. Additionally, the recently developed dhf- variants of the def2- basis sets have been included (dhf-SVP, dhf-TZVPP, …, Weigend, Baldes: J. Chem. Phys. 133, 174102 (2010)); these are reoptimized versions for the more modern MDF effective core potentials (most elements used MDF ECPs already in def2- and are unchanged). Short names like TZVPP (without prefix) now refer to dhf- basis sets.

The MPP and MPPX builds of Molpro have been merged and the decision to run in MPP or MPPX mode made a run-time option. To build parallel Molpro the -mpp flag to configure should be used; the -mppx flag is no longer required or valid. When running parallel Molpro, the --mpp option (default) can be used to specify MPP mode, and the --mppx option for MPPX mode.

New options for configuring parallel Molpro to run on a single node or workstation have been implemented. These options (see manual) are prefixed with -auto, and in conjunction with -mpp configure will download, build and install a variety of different prerequisites for parallel Molpro.

The functionality is essentially the same as in 2009.1, but many bug fixes and small improvements have been added. Please note the following major changes, in particular of the default RI basis sets in explicitly correlated methods as described below.

A faster integral program for density fitting, written by Gerald Knizia, has been added. In particular this speeds up the integral evaluation in F12 calculations by up to a factor of about 10 (depending on the basis set). This program is now used by default, but can be disabled by setting


in the beginning of the input.

Different exponents for the Slater-type geminals can be used for valence-valence, core-valence, and core-core pairs. See manual for details.

For explicitly correlated F12 calculations that use the VnZ-F12 orbital basis sets (OBS), it is now the default to use the corresponding VnZ-F12/OPTRI basis sets to construct the complementary auxiliary orbital basis (CABS). In case that CABS is not used (e.g., in LMP2-F12), the OBS and OPTRI sets are merged automatically. This yields exactly the same results as would be obtained with the CABS approach. In order to use the default RI sets of 2009.1, please specify option RI_BASIS=JKFIT on the command line, or


For compatibility reasons, it is still the default to use the JKFIT sets as RI basis for the AVnZ orbital basis sets. In order to use the corresponding OPTRI sets (where available) please specify option RI_BASIS=OPTRI.

A number of new basis sets have been added to the Molpro library since version 2009.1. The references for these sets can be found in the headers of the respective libmol files.

  • Li, Be, Na, Mg: a) New official versions of the correlation consistent basis sets for these elements have been added, both non-relativistic and those contracted for Douglas-Kroll relativistic calculations. Specifically these are: cc-pVnZ (n=D-5) cc-pwCVnZ (n=D-5) aug-cc-pVnZ (n=D-5) aug-cc-pwCVnZ (n=D-5) and the above with a -DK extension. The older cc-pVnZ basis sets for these elements can still be accessed via the keywords vdz-old, etc. b) New basis sets, including RI and MP2 auxiliary sets, have been added for F12 explicit correlation calculations: cc-pVnZ-F12 (n=D-Q) cc-pCVnZ-F12 (n=D-Q) The optimized CABS auxiliary sets have the same name but with a /OptRI context For MP2 and CCSD auxiliary sets, the cc-pVnZ/MP2FIT sets of Hättig can be used, but a new cc-pV5Z/MP2FIT set has been added that is optimal for the new cc-pV5Z basis set; new aug-cc-pVnZ/MP2FIT (n=D-5) sets have been added as well. The original cc-pV5Z/MP2FIT sets of Hättig have been renamed v5z-old/mp2fit.
  • Cu-Zn, Y-Cd, Hf-Hg: a) While the aug-cc-pVnZ-PP (n=D-5) and cc-pwCVnZ-PP (n=D-5) sets were already available, the combination aug-cc-pwCVnZ-PP was not yet defined. These have now been added for these elements. b) Triple-zeta DK sets have been included now for all of these elements. Unless otherwise noted, these were optimized for 2nd-order DKH. In the cases of Hf-Hg, sets contracted for 3rd-order DKH are also now included: cc-pVTZ-DK cc-pwCVTZ-DK aug-cc-pVTZ-DK aug-cc-pwCVTZ-DK and the above with -DK replaced by -DK3 for DKH3 calculations in the case of Hf-Hg.
  • H-He, B-Ne, Al-Ar, Ga-Kr: a) A variety of DK contracted basis sets have been added for these elements: aug-cc-pVnZ-DK (n=D-5) cc-pCVnZ-DK (n=D-5) cc-pwCVnZ-DK (n=D-5) aug-cc-pCVnZ-DK (n=D-5) aug-cc-pwCVnZ-DK (n=D-5) b) Official cc-pCV6Z and aug-cc-pCV6Z are now also available for Al-Ar c) For explicitly correlated calculations, the core-valence sets have been added: cc-pCVnZ-F12 (n=D-Q) for B-Ne, Al-Ar cc-pCVnZ-F12/OptRI (n=D-Q) for B-Ne, Al-Ar d) cc-pVnZ-F12/OptRI (n=D-Q) as also been added for He
  • Turbomole def2 basis sets: The complete Turbomole def2 basis set family has been added to the Molpro basis library (for all elements H to Rn, except Lanthanides). The def2-orbital basis sets can now be accessed as SV(P), SVP, TZVP, TZVPP, QZVP and QZVPP. In this nomenclature SVP, TZVPP, and QZVPP correspond to valence double-zeta (VDZ), valence triple-zeta (VTZ) and valence quadruple-zeta (VQZ) basis sets, respectively. Auxiliary density fitting basis sets for all elements are available as well (e.g., TZVPP/JFIT, TZVPP/JKFIT, TZVPP/MP2FIT) and are chosen automatically in density-fitted calculations. Supposedly, the JKFIT sets are universal and also applicable in combination with the AVnZ basis sets. Initial results indicate that they also work well with the cc-pVnZ-PP and aug-cc-pVnZ-PP series of basis sets. The orbital basis sets can also be accessed in singly and doubly augmented versions (carrying A or DA prefixes, respectively, e.g., ASVP, DASVP), and the auxiliary fitting sets in singly augmented versions (e.g., ATZVPP/MP2FIT). The old Turbomole basis sets have been renamed; if required, they can be accessed with a def1-prefix (e..g, def1-SVP, def1-TZVPP, etc.).

The ppidd harness that manages interprocess communication has been improved. The performance of the implementation based on pure MPI, as an alternative to use of the Global Arrays toolkit, is considerably improved, through the use of dedicated helper processes that service one-sided remote memory accesses.

The previous way to compute the Davidson correction in multi-state MRCI could lead to non-continuous cluster corrected energies. This is now avoided by ordering the MRCI eigenstates according to increasing energy (previously they were ordered according to maximum overlap with the reference wavefunctions). Furthermore, additional options for computing the Davidson correction in multi-state calculations are added (for details see manual). The old behavior can be recovered using options SWAP,ROTREF=-1.

A variant of the IPEA shift of G. Ghigo, B. O. Roos, and P.A. Malmqvist, Chem. Phys. Lett. 396, 142 (2004) has been added. The implementation is not exactly identical to the one in MOLCAS, since in our program the singly external configurations are not (RS2) or only partially (RS2C) contracted. The shift is invoked by giving option IPEA=shift on the RS2 or RS2C commands; the recommended value for shift is 0.25. For details of the implementation see manual.

The two-point formula for extrapolating the HF reference energy, as proposed by A. Karton and J. M. L. Martin, Theor. Chem. Acc. 115, 330 (2006) has been added: $E_{\rm HF,n}=E_{\rm HF,CBS} +A (n+1)\cdot \exp(-9 \sqrt{n})$. Use METHOD_R=KM for this.

Correlation consistent basis sets for Li, Be, Na, and Mg have been updated to their official versions as reported in Prascher et al., Theor. Chem. Acc. (2010). These now also include core-valence, diffuse augmented, and Douglas-Kroll relativistically contracted versions. The previous sets are still available but have been renamed vdz-old, vtz-old, etc.

Due to new findings, the default behavior of the F12 programs was changed in the following points:

  1. For open-shell systems the default wave function ansatz for was modified. This affects RMP2-F12 and open-shell CCSD-F12 calculations. The new default generally improves open-shell treatments and leads to more consistent behavior. The previous behavior can be restored by explicit,extgen=0 (for more details see manual).
  2. The procedure for the construction of complementary auxiliary basis sets (CABS) and the thresholds were changed. This affects all non-local F12 calculations. The previous behavior can be restored by explicit,ortho_cabs=0,thrcabs=1-7,thrcabs_rel=1e-8
  3. In numeric frequency calculations, the freezing of auxiliary basis sets was improved. This can affect calculations where many redundant functions are deleted.
  4. Pair energies for the explicitly correlated methods can be printed using print,pairs If inner-shell orbitals are correlated, the cc, cv, and vv contributions to the correlation energies are also printed.

The atomic density guess in Hartree-Fock has been improved and extended. Guess basis sets are now available for most atoms and for all pseudopotentials. Most pseudopotentials have been linked to the appropriate basis sets, so that it is sufficient to specify, e.g.


which will select the correlation consistent triple zeta basis sets and the associated (small core) pseudopotential. Similarly, it is mostly sufficient to specify the basis set for other pseudopotential/basis set combinations.

If the wavefunction symmetry is not given in the Hartree-Fock input and not known from a previous calculation, the HF program attempts to determine it automatically from the aufbau pricniple (previously, symmetry 1 was assumed in all cases). For example,


automatically finds that the wavefunction symmetry is 8.

  1. Rationalisation of options for molecular geometry. It is now illegal to specify symmetry and orientation options (eg x;noorient;angstrom) inside a geometry block, which now contains just the geometry specification (Z-matrix or XYZ). Options have to be specified using the new ORIENT and SYMMETRY commands, and/or existing commands such ANGSTROM. This change will, unfortunately, render many inputs incompatible with 2008.1 and earlier versions of Molpro, but has been introduced to allow correct and clean parsing of geometries containing, for example, yttrium atoms, which previously conflicted with the Y symmetry option.
  2. Simplification of geometry input. The program now detects automatically whether the geometry is specified as a Z-matrix, or using cartesian coordinates, and so there is no need any more to set the geomtyp variable. The standard XYZ format is still accepted for cartesian coordinates, but the first two lines (number of atoms, and a comment) can be omitted if desired.

The program now can be built from the source files with the Global Arrays toolkit or the MPI-2 library for parallel execution.

The new features of Molpro version 2008.1 include the following.

  1. Efficient closed-shell and open-shell MP2-F12 and CCSD(T)-F12 methods which dramatically improve the basis set convergence, as described in J. Chem. Phys. 126, 164102 (2007); ibid. 127, 221106 (2007); ibid. 128, 154103 (2008).
  2. Natural bond order (NBO) and natural population analysis (NPA) as described in Mol. Phys. 105, 2753 (2007) and references therein.
  3. Correlation regions within a localized molecular orbital approach as described in J. Chem. Phys. 128, 144106 (2008).
  4. Automated calculation of anharmonic vibrational frequencies and zero-point energies using VCI methods as described in J. Chem. Phys. 126, 134108 (2007) and references therein.
  5. Coupling of DFT and coupled cluster methods as described in Phys. Chem. Chem. Phys. 10, 3353 (2008) and references therein.
  6. Enhanced connections to other programs, including graphical display of output and 3-dimensional structures.
  7. Support for latest operating systems and compilers, including Mac OS X.

Features and enhancements in Molpro version 2006.1 most notably included efficient density fitting methods, explicitly correlated methods, local coupled cluster methods, and several new gradient programs: following:

  1. More consistent input language and input pre-checking.
  2. More flexible basis input, allowing to handle multiple basis sets
  3. New more efficient density functional implementation, additional density functionals.
  4. Low-order scaling local coupled cluster methods with perturbative treatment of triples excitations (LCCSD(T) and variants like LQCISD(T))
  5. Efficient density fitting (DF) programs for Hartree-Fock (DF-HF), Density functional Kohn-Sham theory (DF-KS), Second-order Møller-Plesset perturbation theory (DF-MP2), as well as for all local methods (DF-LMP2, DF-LMP4, DF-LQCISD(T), DF-LCCSD(T))
  6. Analytical QCISD(T) gradients
  7. Analytical MRPT2 (CASPT2) and MS-CASPT2 gradients, using state averaged MCSCF reference functions
  8. Analytical DF-HF, DF-KS, DF-LMP2, and DF-SCS-LMP2 gradients
  9. Explicitly correlated methods with density fitting: DF-MP2-R12/2A’, DF-MP2-F12/2A’ as well as the local variants DF-LMP2-R12/2*A(loc) and DF-LMP2-F12/2*A(loc).
  10. Coupling of multi-reference perturbation theory and configuration interaction (CIPT2)
  11. DFT-SAPT
  12. Transition moments and transition Hamiltonian between CASSCF and MRCI wavefunctions with different orbitals.
  13. A new spin-orbit integral program for generally contracted basis sets.
  14. Douglas-Kroll-Hess Hamiltonian up to arbitrary order.
  15. Improved procedures for geometry optimization and numerical Hessian calculations, including constrained optimization.
  16. Improved facilities to treat large lattices of point charges for QM/MM calculations, including lattice gradients (see section QM/MM interfaces) .
  17. An interface to the MRCC program of M. Kallay, allowing coupled-cluster calculations with arbitrary excitation level.
  18. Automatic embarrassingly parallel computation of numerical gradients and Hessians (mppx Version).
  19. Additional parallel codes, e.g. DF-HF, DF-KS, DF-LCCSD(T) (partly, including triples).
  20. Additional output formats for tables ( XHTML, LaTeX, Maple, Mathematica, Matlab and comma-separated variables), orbitals and basis sets (XML), and an optional well-formed XML output stream with important results marked up.

Relative to version 2002.1, there are the following changes and additions:

  1. Support for IA-64 Linux systems (HP and NEC) and HP-UX 11.22 for IA-64 (Itanium2).
  2. Support for NEC-SX systems.
  3. Support for IBM-power4 systems.
  4. Modified handling of Molpro system variables. The SET command has changed (see sections variables and system variables) .
  5. The total charge of the molecule can be specified in a variable CHARGE or on the WF card (see section defining the wavefunction) .
  6. Improved numerical geometry optimization using symmetrical displacement coordinates (see sections numerical gradients and geometry optimization (OPTG)) .
  7. Improved numerical frequency calculations using the symmetry (AUTO option (see section harmonic vibrational frequencies (FREQUENCIES)) .

Relative to version 2000.1, there are the following principal changes and additions:

  1. Modules direct and local are now included in the base version. This means that integral-direct procedures as described in

M. Schütz, R. Lindh, and H.-J. Werner, Mol. Phys. 96, 719 (1999),
linear-scaling local MP2, as described in
G. Hetzer, P. Pulay, and H.-J. Werner, Chem. Phys. Lett. 290, 143 (1998),
M. Schütz, G. Hetzer, and H.-J. Werner, J. Chem. Phys. 111, 5691 (1999),
G. Hetzer, M. Schütz, H. Stoll, and H.-J. Werner, J. Chem. Phys. 113, 9443 (2000),
as well as LMP2 gradients as described in
A. El Azhary, G. Rauhut, P. Pulay, and H.-J. Werner, J. Chem. Phys. 108, 5185 (1998)
are now available without special license. The linear scaling LCCSD(T) methods as described in
M. Schütz and H.-J. Werner, J. Chem. Phys. 114, 661 (2001),
M. Schütz and H.-J. Werner, Chem. Phys. Lett. 318, 370 (2000),
M. Schütz, J. Chem. Phys. 113, 9986 (2000)
will be made available at a later stage.

  1. QCISD gradients as described in Phys. Chem. Chem. Phys. 3, 4853 (2001) are now available.
  2. Additional and more flexible options for computing numerical gradients and performing geometry optimizations.
  3. A large number of additional density functionals have been added, together with support for the automated functional implementer described in Comp. Phys. Commun. 136 310–318 (2001).
  4. Multipole moments of arbitrary order can be computed.
  5. Further modules have been parallelized, in particular the CCSD(T) and direct LMP2 codes. The parallel running procedures have been improved. The parallel version is available as an optional module.
  6. The basis set library has been extended.
  7. Some subtle changes in the basis set input: it is not possible any more that several one-line basis input cards with definitions for individual atoms follow each other. Each new basis card supercedes previous ones. Either all specifications must be given on one BASIS card, or a basis input block must be used. BASIS,NAME is now entirely equivalent to BASIS=NAME, i.e. a global default basis set is defined and the variable BASIS is set in both cases.
  8. Pseudopotential energy calculations can now be performed with up to $i$-functions, gradients with up to $h$-functions.
  9. Many internal changes have been made to make Molpro more modular and stable. Support has been added for recent operating systems on Compaq, HP, SGI, SUN, and Linux. The patching system has been improved.

Relative to version 98.1, there are the following principal changes and additions:

  1. There was a fundamental error in the derivation of the spin-restricted open-shell coupled-cluster equations in J. Chem. Phys. 99, 5129 (1993) that is also reflected in the RCCSD code in Molpro version 98.1 and earlier. This error has now been corrected, and an erratum has been published in J. Chem. Phys. 112, 3106 (2000). Fortunately, the numerical implications of the error were small, and it is not anticipated that any computed properties will have been significantly in error.
  2. There was a programming error in the transformation of gradients from Cartesian to internal coordinates, which in some cases resulted in slow convergence of geometry optimizations. The error is now fixed.
  3. Vibrational frequencies formerly by default used average atomic masses, rather than those of the most common isotopes, which is now the default behaviour.
  4. MCSCF second derivatives (author Riccardo Tarroni) added (preliminary version, only without symmetry). Frequency and geometry optimization programs are modified so that they can use the analytic Hessian.
  5. New internally contracted multi-reference second-order perturbation theory code (author Paolo Celani) through command RS2C, as described in P. Celani and H.-J. Werner, J. Chem. Phys. 112, 5546 (2000).
  6. EOM-CCSD for excited states (author Tatiana Korona).
  7. QCISD dipole moments as true analytical energy derivatives (author Guntram Rauhut).
  8. Linear scaling (CPU and memory) LMP2 as described by G. Hetzer, P. Pulay, and H.-J. Werner, Chem. Phys. Lett. 290, 143 (1998).

M. Schütz, G. Hetzer, and H.-J. Werner, J. Chem. Phys. 111, 5691 (1999).

  1. Improved handling of basis and geometry records. 98.1 and 99.1 dump files can be restarted, but in case of problems with restarting old files, add RESTART,NOGEOM immediately after the file card. Also, if there are unjustified messages coming up in very large cases about “ORBITALS CORRESPOND TO DIFFERENT GEOMETRY” try ORBITAL,record,NOCHECK. (This can happen for cases with more than 100 atoms, since the old version was limited to 100).
  2. Reorganization and generalization of basis input. Increased basis library.
  3. Counterpoise geometry optimizations.
  4. Improved running procedures for MPP machines. Parallel direct scf and scf gradients are working. These features are only available with the MPP module, which is not yet being distributed.
  5. Important bugfixes for DFT grids, CCSD with paging, finite field calculations without core orbitals, spin-orbit coupling.
  6. Many other internal changes.

As an additional service to the Molpro community, an electronic mailing list has been set up to provide a forum for open discussion on all aspects of installing and using Molpro. The mailing list is intended as the primary means of disseminating hints and tips on how to use Molpro effectively. It is not a means of raising queries directly with the authors of the program. For clearly demonstrable program errors, reports should continue to be sent to; however, ‘how-to’ questions sent there will merely be redirected to this mailing list.

In order to subscribe to the list, send mail to containing the text subscribe; for help, send mail containing the text help.

Messages can be sent to the list (, but this can be done only by subscribers. Previous postings can be viewed in the archive at irrespective of whether or not you subscribe to the list. Experienced Molpro users are encouraged to post responses to queries raised. Please do contribute to make this resource mutually useful.

Molpro98 has the full functionality of Molpro96, but in order to make the code more modular and easier to use and maintain, a number of structural changes have been made. In particular, the number of different records has been significantly reduced. The information for a given wavefunction type, like orbitals, density matrices, fock matrices, occupation numbers and other information, is now stored in a single dump record. Even different orbital types, e.g., canonical, natural, or localized orbitals, are stored in the same record, and the user can subsequently access individual sets by keywords on the ORBITAL directive. New facilities allow the use of starting orbitals computed with different basis sets and/or different symmetries for SCF or MCSCF calculations. The default starting guess for SCF calculations has been much improved, which is most useful in calculations for large molecules. The use of special procedures for computing non-adiabatic couplings or diabatization of orbitals has been significantly simplified. We hope that these changes make the program easier to use and reduce the probability of input errors. However, in order to use the new facilities efficiently, even experienced Molpro users should read the sections RECORDS and SELECTING ORBITALS AND DENSITY MATRICES in the manual. It is likely that standard Molpro96 inputs still work, but changes may be required in more special cases involving particular records for orbitals, density matrices, or operators.

All one-electron operators needed to compute expectation values and transition quantities are now stored in a single record. Operators for which expectation values are requested can be selected globally for all programs of a given run using the global GEXPEC directive, or for a specific program using the EXPEC directive. All operators are computed automatically when needed, and the user does not have to give input for this any more. See section ONE-ELECTRON OPERATORES AND EXPECTATION VALUES of the manual for details.

Due to the changed structure of dump and operator records, the utility program MATROP has a new input syntax. Molpro96 inputs for MATROP do not work any more.

In addition to these organizational changes, a number of new programs have been added. Analytic energy gradients can now be evaluated for MP2 and DFT wavefunctions, and harmonic vibrational frequencies, intensities, and thermodynamic quantities can be computed automatically using finite differences of analytical gradients. Geometry optimization has been further improved, and new facilities for reaction path following have been added.

An interface to the graphics program MOLDEN has been added, which allows to visualize molecular structures, orbitals, electron densities, or vibrations.

Integral-direct calculations, in which the two-electron integrals in the AO basis are never stored on disk but always recomputed when needed, are now available for all kinds of wavefunctions, with the exception of perturbative triple excitations in MP4 and CCSD(T) calculations. This allows the use of significantly larger basis sets than was possible before. The direct option can be selected globally using the GDIRECT command, or for a specific program using the DIRECT directive. See section INTEGRAL DIRECT METHODS in the manual for details. Note that the DIRECT module is optional and not part of the basic Molpro distribution.

Local electron correlation methods have been further improved. In combination with the integral-direct modules, which implement efficient prescreening techniques, the scaling of the computational cost with molecular size is dramatically reduced, approaching now quadratic or even linear scaling for MP2 and higher correlation methods. This makes possible to perform correlated calculations for much larger molecules than were previously feasible. However, since these methods are subject of active current research and still under intense development, we decided not to include them in the current Molpro release. They will be optionally available in one of the next releases.