Hi Everybody,
I am new to Psi4. Would somebody be so kind and check my implementation of the MP2.X method.
I implemented a test run for a H2O molecule and think I got it right but would like somebody to take a look just to be sure.
The basic idea is that I run:
- a MP3/6-31+g* single point calculation
-> ec_mp2 mp2_correlation energy
-> ec_mp3 mp3 correlation energy - a MP2/CBS calculation with the cc-pv[dt]z basis set
-> e_mp2_cbs
And then compute:
E(MP2.X/CBS) = e_mp2_cbs + C * (ce_mp3 - ce_mp2)
With C taken from http://onlinelibrary.wiley.com/doi/10.1002/cphc.201200850/abstract
The input file is below.
Thanks,
Alberto Gobbi
#! MP2.X calculation of H2O
memory 500 mb
molecule h2o {
O
H 1 1.0
H 1 1.0 2 104.5
}
set freeze_core true
# MP3/6-31+g* MP3/6-31+g* MP3/6-31+g* MP3/6-31+g* MP3/6-31+g*
set basis 6-31+g*
energy(‘mp3’)
e_scf = psi4.get_variable(‘SCF TOTAL ENERGY’)
ce_mp2 = psi4.get_variable(‘MP2 CORRELATION ENERGY’)
ce_mp3 = psi4.get_variable(‘MP3 CORRELATION ENERGY’)
mp3_basis = psi4.get_global_option(“BASIS”)
clean()
# MP2/CBS MP2/CBS MP2/CBS MP2/CBS MP2/CBS MP2/CBS MP2/CBS
cbs_basis = ‘cc-pv[dt]z’
set guess read
cbs(name=‘mp2’, scf_wfn=‘scf’, corl_wfn=‘mp2’, corl_basis=cbs_basis, corl_scheme=corl_xtpl_helgaker_2)
e_mp2_cbs = psi4.get_variable(‘CBS TOTAL ENERGY’)
# compute MP2.X
# coefficients for mp2.x from doi: 10.1002/cphc.201200850
mp2_x_c = { ‘3-21++g*’ : 1.12,
‘6-31g*’ : 0.86,
‘6-31+g*’ : 0.73,
‘6-311++g**’ : 0.65,
‘def2-tzvp’ : 0.67,
‘def2-tzvpp’ : 0.53,
‘cc-pvdz’ : 0.72,
‘aug-cc-pvdz’: 0.52,
‘aug-cc-pvtz’: 0.50 }
c = mp2_x_c[mp3_basis.lower()]
cor_mp2_x = c * (ce_mp3 - ce_mp2)
e_mp2_x = e_mp2_cbs + cor_mp2_x
psi4.print_out("\n\nMP2.X MP2.X MP2.X MP2.X MP2.X MP2.X MP2.X MP2.X \n")
psi4.print_out(“MP2.X (MP3(%s) CBS(%s)) C=%4.2f\n” %
(mp3_basis, cbs_basis, c))
psi4.print_out(" SCF TOTAL ENERGY = %12.7f\n" % e_scf)
psi4.print_out(" MP2 TOTAL ENERGY = %12.7f\n" % (e_scf + ce_mp2))
psi4.print_out(" MP3 TOTAL ENERGY = %12.7f\n" % (e_scf + ce_mp3))
psi4.print_out(" CBS TOTAL ENERGY = %12.7f\n" % e_mp2_cbs)
psi4.print_out(" MP2.X CORRECTION = %12.7f\n" % cor_mp2_x)
psi4.print_out(" MP2.X TOTAL ENERGY = %12.7f\n" % e_mp2_x)