Assignment 1, due date: February 10, 2010 (Dr. Piotr Piecuch)

Dipole moment function for the HF molecule (Gaussian and GAMESS calculations).
------------------------------------------------------------------------------
The main purpose of this assignment is the calculation of the dipole
moment function for the HF molecule (this function is
important in calculations of the infrared absorption line strengths).
In general, the dipole moment function describes the dependence
of the electronic dipole moment of a molecule on the nuclear
geometry. We would like to compare the results of calculations
performed with the following three methods:
the restricted Hartree-Fock approach (RHF), the coupled-cluster
singles and doubles method (CCSD, in two different ways),
and the full configuration interaction (FCI) approach. The FCI
approach gives the exact dipole moment function for a given basis set.
This will be your reference for the analysis of the results.
The RHF approach does not describe many-electron correlation
effects. These are described by the CCSD method (and FCI, of course).
You will be able to learn what the role of the electron correlation
effects in describing the dipole moment function of the HF molecule is,
and experience the significance of electron correlation in describing
bond breaking.

-------------------------------------------------------------------------------

Instructions:

GAUSSIAN 98 RUNS

1. Prepare inputs for the CCSD (coupled-cluster singles and doubles)
calculations for the HF molecule, as described by the DZ basis set.
The calculations are for the dipole moment of HF as a function
of the internuclear distance R. Here is the set of R values
that we want to obtain the dipole moments for:
1.12632, 1.2996, 1.7328 (the approximate equilibrium bond length),
2.166, 2.5992, 3.0324, 3.4656, 4.332, and 5.1984
(all in bohr).

The Gaussian input decks should have the following lines:

# rccsd=(rw,conver=9)/gen UNITS=Bohr SCF=(QC,Conver=9)
# Polar=EnOnly

for R < 4.332 bohr, and

# rccsd=(rw,conver=8)/gen UNITS=Bohr SCF=(QC,Conver=8)
# Polar=EnOnly

for R = 4.332 and 5.1984 bohr.

In describing the nuclear geometry, use the Cartesian coordinates
rather than the Z matrix, orienting the HF molecule along
the z axis in the following way:

F   0.0  0.0  0.0
H   0.0  0.0  R  

for example,

F   0.0  0.0  0.0
H   0.0  0.0  1.7328

for R = 1.7328 bohr. This is important, since dipole moment is a vector
whose components depend on the choice of the molecular frame.

The DZ basis set of Dunning which we will use has the following form:

H 0
    S    3
    1.924060D+01   3.282800D-02
    2.899200D+00   2.312080D-01
    6.534000D-01   8.172380D-01
    S    1
    1.776000D-01   1.000000D+00
****
F 0
    S    6
    9.994790D+03   2.017000D-03
    1.506030D+03   1.529500D-02
    3.502690D+02   7.311000D-02
    1.040530D+02   2.464200D-01
    3.484320D+01   6.125930D-01
    4.368800D+00   2.424890D-01
    S    1
    1.221640D+01   1.000000D+00
    S    1
    1.207800D+00   1.000000D+00
    S    1
    3.634000D-01   1.000000D+00
    P    4
    4.435550D+01   2.086800D-02
    1.008200D+01   1.300920D-01
    2.995900D+00   3.962190D-01
    9.383000D-01   6.203680D-01
    P    1
    2.733000D-01   1.000000D+00
****

Make sure that all electrons are correlated in the CCSD runs, which
means that you should add the following window at the end of your
input file (below the last ****):

(...)
****

1,0

An example of an input file for the CCSD calculation for
R = 1.7328 bohr will be provided to you during the
lab on Wednesday, February 3, 2010 or you can copy it
yourself to your directory by typing

cp ~piecuch/cem888_SS10/students1/hfccsd1.7328.com ./

2. Prepare similar input files for the RHF calculations.
The Gaussian input deck should begin with the line:

# rhf/gen UNITS=Bohr SCF=(QC,Conver=10)

An example of an input file for the RHF calculation for
R = 1.7328 bohr will be provided to you during the
lab on Wednesday, February 3, 2010 or you can copy it
yourself to your directory by typing

cp ~piecuch/cem888_SS10/students1/hfrhf1.7328.com ./

3. Once you are done with the calculations, read the values of the
dipole moment from the output files. In each case,
go to the end of the output file and read the value
of the dipole reported as \Dipole= ...
You will see three numbers there corresponding to the x, y, and
z components of the dipole (z is the internuclear axis).
Use the z component for plotting and reporting.
The dipole moment values are in atomic units (a.u. or bohr-electron;
1 a.u. = 2.541766 Debye). Use atomic units for reporting and plotting.
In these units, the experimental value of the dipole moment
(the vibrationless value) is 0.707 a.u.

GAMESS RUNS:

1. Prepare inputs for the CCSD property calculations
for the HF molecule, as described by the DZ basis set.
Again, the calculations are for the dipole moment of HF as a function
of the internuclear distance R using the same set of R values
as used in the Gaussian calculations.

The GAMESS input files should have the following format

 $CONTRL SCFTYP=RHF CCTYP=CCSD
         RUNTYP=ENERGY UNITS=BOHR $END
 $SYSTEM TIMLIM=9999999 MWORDS=20 $END
 $BASIS GBASIS=DZDUN EXTFIL=.TRUE. $END
 $CCINP CCPRP=.TRUE. NCORE=0 $END
 $DATA
HF
Cnv  2

F   9.0   0.0  0.0  0.0
H   1.0   0.0  0.0  R
 $END

for R = 1.12632, 1.2996, 1.7328, 2.166, 2.5992, 3.0324, 3.4656, and 4.332 bohr.

For R = 5.1984 bohr, add the following additional line:

 $SCF DIIS=.TRUE. SOSCF=.FALSE. $END

Note that as in the case of Gaussian calculations, all electrons are
correlated (NCORE=0; the default in GAMESS is to freeze core electrons).
Use the same DZ basis set as in the Gaussian calculations (called DZDUN
in the above input example). The easiest thing to do (the above input
format reflects on this) is to place the basis set in a separate file
in the same directory as your GAMESS inputs are (in this way, you can build
your own library of basis sets, in addition to the basis sets stored
in the GAMESS source code, and keep inputs short).
In order for GAMESS to read the basis set informaton from the external
file (EXTFIL=.TRUE. in the $BASIS input group), you must submit
the Gamess runs with the -b flag in gmssub, for example,
'gmsub -b basis hfccsd1.7328', where 'basis' is the name
of your file containing the basis set information.
You must make sure that the basis set is formatted using
the GAMESS style, which is slightly different than the Gaussian style.
In the case of the DZ basis set that interests us here,
you must use something like this:

H DZDUN
 S   3
  1        19.24060000         0.3282800000E-01
  2        2.899200000         0.2312080000
  3       0.6534000000         0.8172380000
 S   1
  1       0.1776000000          1.000000000

for the hydrogen atom and

F DZDUN
 S   6
  1        9994.790000         0.2017000000E-02
  2        1506.030000         0.1529500000E-01
  3        350.2690000         0.7311000000E-01
  4        104.0530000         0.2464200000
  5        34.84320000         0.6125930000
  6        4.368800000         0.2424890000
 S   1
  1        12.21640000          1.000000000
 S   1
  1        1.207800000          1.000000000
 S   1
  1       0.3634000000          1.000000000
 P   4
  1        44.35550000         0.2086800000E-01
  2        10.08200000         0.1300920000
  3        2.995900000         0.3962190000
  4       0.9383000000         0.6203680000
 P   1
  1       0.2733000000          1.000000000

for fluorine.

An example of an input file for the CCSD calculation for
R = 1.7328 bohr will be provided to you during the
lab on Wednesday, February 3, 2010 or you can copy it
yourself to your directory by typing

cp ~piecuch/cem888_SS10/students1/hfccsd1.7328.inp ./

2. Once you are done with the calculations, read the values of the
dipole moment from the output files. The RHF dipoles, in Debye, are
listed a few lines above

     ---------------------------
     COUPLED CLUSTER CALCULATION
     ---------------------------

The CCSD dipoles can be found in two places,
a few lines below the line GROUND STATE CCSD PROPERTIES in a.u., and
and at the end of the output file in Debye.
Use atomic units for reporting and plotting.
As mentioned earlier, the unit conversion
factor is 1 a.u. = 2.541766 Debye.

THINGS TO DO:

1. Report the results by preparing the following table:

R      RHF(Gaussian) RHF(Gamess) CCSD(Gaussian) CCSD(Gamess)  FCI
--------------------------------------------

1.12632      ...         ...          ...          ...        ...
...

etc.

where, as explained earlier, FCI stands for the full configuration
interaction approach (which is the exact result for the DZ basis set).

The FCI/DZ values of the HF dipole moment (the z component; the
x and y components are zero) are:

R (bohr)    Dipole moment (a.u.)

1.12632     0.75852
1.2996      0.80163
1.7328      0.89789
2.166       0.94752
2.5992      0.92262
3.0324      0.80517
3.4656      0.60494
4.332       0.20937
5.1984      0.05006

2. Prepare a plot with the dipole moment functions obtained with
the RHF, CCSD, and FCI methods (use the FCI results provided
above). Do it twice, for the GAMESS and Gaussian calculations.
You can use gnuplot or any other plotting program you like.
A legend and labels for axes are important.
Do not forget about them.

3. Comment on the results. Why is the RHF curve so different
than FCI or CCSD? What may, in your view, cause the wrong behavior
of the RHF dipole moment function, whose values do not approach zero as R
approaches infinity (as they should, see the FCI results)? The simple
answer "lack of correlation effects in RHF" will not suffice. Try to think
about the unphysical components of the RHF wave function that are
no longer present when the exact, FCI approach is employed.
The dipole moment of HF should vanish at R=infinity, since
in this limit the dipole moment of HF becomes a sum of dipole moments
of the H and F atoms.

4. Are the dipole moments values obtained in the RHF Gaussian and GAMESS
calculations identical? Are the dipole moments values obtained in the
CCSD Gaussian and GAMESS
calculations identical? Do not worry about the differences at the level of
10^-7 a.u. that are related to convergence criteria used in the calculations,
which were not the tightest possible, and focus on the first six decimal
places. This is a little tricky, but can you
try to explain the differences, if any, based on
the contents of your output files and the way the Gaussian and GAMESS
dipole moment calculations are executed (outputs are rather transparent
about it, but you can always consult the relevant manuals if need be).

5. Examine the Gamess outputs and look for the largest T1 and largest T2 cluster
amplitudes printed after the CCSD iterations, which reflect on
the role of the leading singly and doubly excited configurations
in the wave function expansion relative to the RHF reference determinant.
Focus on the largest cluster amplitude in each category
and create a table of the values of these
amplitudes as a function of the H-F internuclear distance R. Can you
determine what orbital excitations are involved in these amplitudes
(using the language of sigma, pi, etc. molecular orbitals of the HF molecule)?
The outputs have all the information. Can you describe and explain the observed
dependence of the largest T1 and T2 amplitudes on the H-F distance R?
Please elaborate.