Hello,
I am very newbie user of psi4. I want to calculate ESP charges with HF/6-31G basis set. I tried this input file which I found on github (https://github.com/psi4/psi4/blob/master/samples/props4/test.in) and I found ESP charges in grid_esp.dat file. But in gaussian09 I obtained slightly different charges. Also I want to change this file to calculate ESP charges on DFT B3LYP 6,311G++(2d,2p) how can I do this?
Thank you.
#! Electrostatic potential and electric field evaluated on a grid around water.
molecule h2o {
noreorient
nocom
O 0.250254404867 0.126248114412 0.000000000000
H 0.428893090449 1.055731838795 0.000000000000
H 1.104987458381 -0.280303532167 0.000000000000
}
set basis cc-pvdz
with open('grid.dat', 'w') as fp:
for x in range(3):
xval = (x-1.0)*2.0
for y in range(3):
yval = (y-1.0)*2.0
fp.write("%16.10f%16.10f%16.10f\n" % (xval, yval, 1.0))
set basis cc-pvdz
E, wfn = prop('scf', properties=["GRID_ESP", "GRID_FIELD"], return_wfn=True)
Vvals = wfn.oeprop.Vvals()
Exvals = wfn.oeprop.Exvals()
Eyvals = wfn.oeprop.Eyvals()
Ezvals = wfn.oeprop.Ezvals()
set basis cc-pvtz
E, wfn = prop('scf/cc-pvdz', properties=["GRID_ESP", "GRID_FIELD"], return_wfn=True)
Vvals_2 = wfn.oeprop.Vvals()
Hello, thank you for your answer, it really helped. I understood that when I want to use DFT/6-31G++(2d,2p) I must use b3lyp/6-311++G(2d,2p).For example the input file in the question above I must change like this : wfn = prop(‘b3lyp/6-311++G(2d,2p)’, properties=[“GRID_ESP”, “GRID_FIELD”], return_wfn=True) and also I should change set_basis lines.
#! Electrostatic potential and electric field evaluated on a grid around XX.
molecule XX {
0 1
C 0.00000000 0.00000000 0.77000000
H -0.50355600 0.87414000 1.12666600
H 1.00880500 -0.00097800 1.12666600
H -0.50525000 -0.87316200 1.12666600
C 0.00000000 0.00000000 -0.77000000
H 0.50525000 -0.87316200 -1.12666600
H -1.00880500 -0.00097800 -1.12666600
H 0.50355600 0.87414000 -1.12666600
}
set basis 6-311++G(2d,2p)
# Make a regular grid to evaluate properties on
with open('grid.dat', 'w') as fp:
for x in range(3):
xval = (x-1.0)*2.0
for y in range(3):
yval = (y-1.0)*2.0
fp.write("%16.10f%16.10f%16.10f\n" % (xval, yval, 1.0))
set basis 6-311++G(2d,2p)
E, wfn = prop('b3lyp', properties=["GRID_ESP", "GRID_FIELD"], return_wfn=True)
Vvals = wfn.oeprop.Vvals()
Exvals = wfn.oeprop.Exvals()
Eyvals = wfn.oeprop.Eyvals()
Ezvals = wfn.oeprop.Ezvals()
set basis 6-311++G(2d,2p)
E, wfn = prop('b3lyp/6-311++G(2d,2p)', properties=["GRID_ESP", "GRID_FIELD"], return_wfn=True)
Vvals_2 = wfn.oeprop.Vvals()
Now, am I doing right so far ?
if I am doing right I will go the second question (about ESP)
That’s right! It’s more complicated than what you need, but it’s correct.
You aren’t doing anything with the Vvals, Exvals, and so forth variables, so you can remove them. You also don’t use the basis set when writing the grid.dat file, so you can move that over. This gives you the following simplified input file:
#! Electrostatic potential and electric field evaluated on a grid around XX.
molecule XX {
0 1
C 0.00000000 0.00000000 0.77000000
H -0.50355600 0.87414000 1.12666600
H 1.00880500 -0.00097800 1.12666600
H -0.50525000 -0.87316200 1.12666600
C 0.00000000 0.00000000 -0.77000000
H 0.50525000 -0.87316200 -1.12666600
H -1.00880500 -0.00097800 -1.12666600
H 0.50355600 0.87414000 -1.12666600
}
# Make a regular grid to evaluate properties on
with open('grid.dat', 'w') as fp:
for x in range(3):
xval = (x-1.0)*2.0
for y in range(3):
yval = (y-1.0)*2.0
fp.write("%16.10f%16.10f%16.10f\n" % (xval, yval, 1.0))
set basis 6-311++G(2d,2p)
E, wfn = prop('b3lyp', properties=["GRID_ESP", "GRID_FIELD"], return_wfn=True)
As for your other question, why Psi’s ESP charges disagrees with those from Gaussian09. First, double-check that you’re using the same method both. The input file you posted originally will compute with Hartree-Fock and the cc-pVDZ basis set. Is that the same in your Gaussian input file? Second, how large of a disagreement are you seeing? Is this in the tenth decimal place, the third…?
What do you want to do with the ESP charges?
RESP fit? Or plotting an ESP surface? For these we have special tools.
6-31G in Gaussian uses cartesian D functions, i.e. you have a slightly different basis set. For 6-311G… it should be the same (both use spherical basis functions)
Hello, thank you again for your valuable answers. I want to calculate ESP charges to perform molecular dynamic simulation. My molecules are boron compounds so RESP fitting is a little bit tricky because boron is non-standart atom for AMBER. Thus, I don’t want RESP-fitting charges ESP charges is enough for me.
I checked in gaussian09 basis sets, basis sets are the same but still I got different results. I tested with and without ReadRadii keyword ( when I used ReadRadii, I assigned 1.70 radii as for boron atom as antechamber suggested) but I got different results in psi4 and gaussian09. The Difference was between 4%-7%. Still I couldn’t figure out.
The input posted for PSI4 calculates only ESP values on a grid. These are not atomic ESP charges.
I don’t think PSI4 has atomic ESP charges out of the box. There is a RESP plugin, which may be able to print you ESP charges as well. However, to get agreement with Gaussian one needs to construct a grid that is sufficiently similar.
Dear hokru, after your reply I found that I can use psi4 files in multiwfn and I calculated ESP and even RESP charges with multiwfn after I obtained grid. Thank you so much, man.
