# The ENERGY-SCAN module

In CamCASP it is possible to compute the first-order energies $E^{(1)}_{\rm elst}$ and $E^{(1)}_{\rm exch}$ for a large number of dimer configurations without repeating the wavefunction calculation on the monomers. This is possible only because the first-order energies can be accurately evaluated in a monomer basis set, that is, without the presence of the mid-bond and far-bond functions that are necessary for converging the second-order interaction energy components. In this tutorial we will cover how these energies can be computed using the ENERGY-SCAN module, and additionally, how the distributed density-overlap integrals can be computed.

Important

Reference Geometries The energy scan with module ENERGY-SCAN works on the principle that we computer the molecular orbitals and energies for the molecule in their reference orientations, and then rotate/translate them as needed to compute the interaction energies for the interacting dimer. Consequently, it is important that this reference orientation must be used to define all molecules in the CamCASP command files.

Important

//Integral Accuracy//: When the DF is done accurately - here when the DC-aux basis is used - then it seems that we must use 1-electron integrals computed as accurately as is possible. So we use Integral switch = 1. But if we use a relatively poor DF (as is the case with the MC-MC case when we use the MC-aux basis type), then we need to rely on an error cancellation in the terms that appear in $E^{(1)}_{\rm elst}$. So we need to use Integral switch = 0.

## Commands in ENERGY-SCAN

BEGIN SCAN/ENERGY-SCAN
{PROBE || SCAN} [molecule || atom] <molname> +++
WITH [<molname> [[and] CHARGE <charge>]]
UNITS [A || ANG || ANGSTROM || BOHR || AU || A.U.]  +++
NOTE: UNITS must come before POINTS/RANDOM
CENTRE/CENTER PROBE [on] [COM]/[XYZ <x> <y> <z>]
REGULARIZE [E1] [[and] E2] [with] ETA [=] <reg_eta>
{ENERGIES || ENERGY} [DISP/DISPERSION/E2DISP] +++
[&/and/,] [IND/INDUCTION/E2IND] [&/and/,] [ELST/E1ELST/ELECTROSTATICS] +++
[&/and/,] [E1EXCH] [&/and/,] [E2EXIND] & [E2EXDISP] & [OVERLAP] [NOTHING]
LATTICE
POINTS
[TRANSLATIONS || TRANSLATIONS-ONLY]
[SKIP-FIRST-COLUMN]
x1 y1 z1    [angle1 nx1 ny1 nz1]
x2 y2 z2    [angle2 nx2 ny2 nz2]
etc
---  <----end POINTS block with three hyphens.
RANDOM
POINTS [=] <number of points>
DRMIN  [=] <dRmin>
DRMAX  [=] <dRmax>
SKIP   [=] <number to skip> [angular] [configurations]
WRITE  <filename> [APPEND]
---
ENERGY-FILE [PREFIX] <file prefix> [[in || using] STYLE <style index>]
OVERLAP-FILE [PREFIX] <file prefix> [[in || using] STYLE <style index>]
DEBUGGING T/TRUE/ON or F/FALSE/OFF
{QUIET || VERBOSE}
END

## First-order energies: $E^{(1)}_{\rm elst}$ and $E^{(1)}_{\rm exch}$

### MC main basis and MC auxiliary basis

Important points to get right:

• REDO-DF-ON-ROTATION False : Needed to avoid re-calculation of the monomer density-fitting as the monomers are rotated.

• Auxiliary basis sets should be Spherical in this case as CamCASP can only rotate integrals that use spherical angular momenta.

• Integral switch = 0: force 1-electron integrals to use density-fitting. Results in far more accurate $E^{(1)}_{\rm elst}$ energies. But with an MC-MC type of basis this energy will always be inaccurate!

• DF-TYPE-MONOMER OO and DF-TYPE-DIMER OO: For first-order energies only integrals in the OO space are needed. This causes the DF-INTEGRAL module to change its defaults and use the smaller space. Results in much faster execution.

! This is important as it save a lot of computer time
! As both the main and auxiliary basis sets are type MC,
! all objects, including the DF solution matrices can be
! rotated with the molecules. So no re-calculations are needed.

SET DF
REDO-DF-ON-ROTATION   False
END

! Various DF parameters.

SET DF
Molecule pyridine1
Type OO
Eta = 0.0
Lambda = 0.0
Gamma = 0.0
Print only normalization constraints
Solver LU
END
SET DF
Molecule pyridine2
Type OO
Eta = 0.0
Lambda = 0.0
Gamma = 0.0
Print only normalization constraints
Solver LU
END

! Only occupied-occupied (OO) type integrals are needed.
! This saves computer time. The default is to compute
! integrals in the OV and VV spaces too. Not needed for
! first-order energies.

Set DF-INTS
DF-TYPE-MONOMER  OO
DF-TYPE-DIMER    OO
End

! Set the integral switches in the energy modules.
! Switch = 1 causes the code to compute the overlap and nuclear
! integrals without density-fitting (i.e. they are more accurate), but this introduces
! a large error in E(1)elst as there is no error cancellation.
! For some reason this error cancellation is needed when we use an MC-type
! aux basis. So here we use switch = 0 to force CamCASP to use the 1-electron
! integrals with density-fitting. This results in more accurate E(1)elst, but
! bear in mind that E(1)elst with an MC-type main basis will never be
! accurate enough.

Set E1elst
Integral switch = 0
End
Set E1exch
Integral switch = 0
End

Begin Energy-Scan
Probe pyridine1 with pyridine2
Energy  E1elst & E1exch
Units Bohr
! Debug
#include pyr2-scan.geom
End

### MC main basis and DC auxiliary basis

This is the more accurate way of calculating the $E^{(1)}_{\rm elst}$ energy. For reasons we still do not understand, the DC-type auxiliary basis results in far more accurate electrostatic energies. The $E^{(1)}_{\rm exch}$ energy seems more insensitive to this. As the density-fitting will be repeated for each dimer configuration, this scan will be much more computationally expensive than the MC/MC case. However it will also be more accurate.

## Distributed density-overlap integrals

The distributed density-overlap integrals needed for the density-overlap model can be calculated using the ENERGY-SCAN module using the keyword:

  Energy Overlap

This can be done using two algorithms:

• The density-fitting-based partitioning, and
• the ISA-based partitioning.

The latter is much more stable and may be used with extensive basis sets. The former is much faster but must be used with relatively small auxiliary basis sets.

### DF-based distribution

Basis restrictions:

• Small auxiliary basis sets: TZVPP is OK but does lead to artifacts.
• The auxiliary basis must used spherical GTOs.
• Both the main and auxiliary basis sets must be type MC.

### ISA-based distribution

Basis restrictions:

• At present, the s-functions in the auxiliary basis must be those required for the ISA. This will soon change (it may have already changed, but I need to check the code!)
• The auxiliary basis must used spherical GTOs.
• Both the main and auxiliary basis sets must be type MC.

Other important points to get right:

• REDO-DF-ON-ROTATION False : Needed to avoid re-calculation of the monomer density-fitting as the monomers are rotated.

• Auxiliary basis sets should be Spherical in this case as CamCASP can only rotate integrals that use spherical angular momenta.

• Integral switch = 0: force 1-electron integrals to use density-fitting. Results in far more accurate $E^{(1)}_{\rm elst}$ energies. But with an MC-MC type of basis this energy will always be inaccurate!

• DF-TYPE-MONOMER OO and DF-TYPE-DIMER OO: For first-order energies only integrals in the OO space are needed. This causes the DF-INTEGRAL module to change its defaults and use the smaller space. Results in much faster execution.

# Accuracy : MC-DC versus MC-MC and the choice of integral switch

Here is some data for the thiophene dimer in ten random dimer orientations:

#        Reference energies    ------------------MC main basis-------------------------
#      ---augA-Sadlej-MC+----  -----aTZ DC-aux sw=1-------   -----aTZ DC-aux sw=0------
# INDX    E1elst     E1exch    E1elst        E1exch          E1elst         E1exch
1       -0.82874     5.35489 -6.6592046E-01  4.9352442E+00  -2.1537427E-01  4.9184553E+00
2       -9.69716    29.56556 -9.6704016E+00  3.0054717E+01  -8.9839291E+00  3.0025835E+01
3       -6.73208    21.89199 -6.6285621E+00  2.1964070E+01  -5.9417571E+00  2.1931287E+01
4       -1.21089     6.55433 -1.0410530E+00  6.1161696E+00  -5.2299977E-01  6.0962712E+00
5       -2.48791     6.39493 -2.3995603E+00  6.2131489E+00  -2.2424074E+00  6.2058767E+00
6       -1.28705     3.16231 -1.2235851E+00  3.0334854E+00  -1.2353517E+00  3.0317530E+00
7      -23.70028    68.83931 -2.3378553E+01  6.9401212E+01  -2.5423431E+01  6.9406683E+01
8       -1.24095     3.04211 -1.1786025E+00  2.9168411E+00  -1.1988017E+00  2.9152978E+00
9        0.36216     1.18635  4.3647476E-01  1.2357284E+00   5.9249559E-02  1.2331486E+00
10      -1.36398    10.34951 -1.2889322E+00  1.0329523E+01  -4.4949354E-01  1.0303064E+01

Here's the same for the MC-MC case:

#        Reference energies    ------------------MC main basis-------------------------
#      ---augA-Sadlej-MC+----  -----aTZ MC-aux sw=1-------   -----aTZ MC-aux sw=0------
# INDX    E1elst     E1exch    E1elst        E1exch          E1elst         E1exch
1       -0.82874     5.35489  4.3197305E+00  4.9326346E+00  -9.0835363E-01  4.9198501E+00
2       -9.69716    29.56556 -2.6562133E+00  3.0045234E+01  -1.1141111E+01  3.0029574E+01
3       -6.73208    21.89199  4.0413476E-02  2.1955698E+01  -7.8247516E+00  2.1935417E+01
4       -1.21089     6.55433  4.1474269E+00  6.1129764E+00  -1.3724563E+00  6.0980295E+00
5       -2.48791     6.39493  2.8214742E+01  6.2149486E+00  -5.0185247E+00  6.2141753E+00
6       -1.28705     3.16231  2.6948757E+01  3.0337061E+00  -3.0910482E+00  3.0362936E+00
7      -23.70028    68.83931  1.1345649E+01  6.9423383E+01  -2.7840501E+01  6.9372187E+01
8       -1.24095     3.04211  2.6857644E+01  2.9170139E+00  -3.0089173E+00  2.9196679E+00
9        0.36216     1.18635  8.4263178E+00  1.2338500E+00   1.8745135E-01  1.2334557E+00
10      -1.36398    10.34951  5.0912047E+00  1.0308759E+01  -3.4663106E-01  1.0303721E+01

Warning

The basis sets used in these calculations are not equivalent: the augA-Sadlej* basis is nearly as good as the Dunning aug-cc-pVTZ basis, but it is not the same. Besides, it is used in the MC+ type, with mid-bonds and far-bonds. So treat the comparison with some caution.

* : augA-sadlej is a basis developed for intermolecular interactions by AJM. It is available in CamCASP 6.x.

Here's the above data analysed:

 MC-DC MC-MC E(1)elst E1exch E(1)elst E1exch sw=1 sw=0 sw=1 sw=0 sw=1 sw=0 sw=1 sw=0 rmse 0.14 0.77 0.32 0.32 20.04 1.85 0.32 0.31 mad 0.11 0.61 0.25 0.25 16.09 1.42 0.25 0.25 rel rmse 0.10 0.45 0.04 0.04 13.14 0.77 0.04 0.04

Conclusions:

• $E^{(1)}_{\rm exch}$ is relatively insensitive to the type of auxiliary basis. The small relative rms errors suggest that these energies can be used fairly reliably with either of the integral switches and either of the MC-DC or MC-MC auxiliary basis sets.

• $E^{(1)}_{\rm elst}$ is very sensitive to the switch. For the MC-DC type of basis, the switch = 1 should be used. Even so, the errors made in this energy are relatively large. However for the MC-MC type of basis, the switch = 0 must be used: errors are enormous with switch = 1.

AJMPublic/camcasp/energy-scan-module (last edited 2021-04-14 11:44:02 by apw109)