How to setup a dihedral scan

I am a first time Psi4 user and I am trying to setup an input file for performing a dihedral scan around a rotatable bond. I am having problems assigning Python variables inside a “set optking” command. Any help will be greatly appreciated. I would like the following code to work with Cartesian input.

# Constrained minimization with H-O-O-H dihedral fixed to 120.0 degrees.
molecule HOOH {
    H
    O 1 0.90
    O 2 1.40 1 100.0
    H 3 0.90 2 100.0 1 115.0
}
# optimize molecule at different dihedral angles
# Psithon wraps strings in quotes by default
# use $ to tell Psithon that $var is a Python variable
steps = 36
step_size = 10
start_angle = 0
for counter in range (0, steps):
    # set smaller basis set for geometry HT geometry optimization
    set {
        basis 6-31G*
        scf_type df
        guess sad
        ints_tolerance 1.0E-8
    }

    dihedral = 1.0 * counter * step_size + start_angle
    dihedral_string =  "(\" 1 2 3 4 %.1f \")" % dihedral
    # how to set optking so that dihedral angle changes
    # for each step?
    set optking {
        frozen_dihedral = ""
        fixed_dihedral = (" 1 2 3 4 $dihedral ")
    }
    optimize('scf')   # geometry optimization with HF

# change basis set for single point energy calculation
    set {
        basis 6-311+G**
        scf_type df
        guess sad
    }
    DFT_energy = energy ('b3lyp')  # single point energy
    MP2_energy = energy ('df-mp2')  # calculate single point energy with MP2

    print "Dihedral=%.1f" % dihedral
    print "DFT_Energy(kcal/mol)=%.4f" % (DFT_energy  * psi_hartree2kcalmol)
    print "MP2_Energy(kcal/mol)=%.4f" % (MP2_energy * psi_hartree2kcalmol)

Thanks,

JW

Hey JW,

Changing the following:

dihedral_string =  "(\" 1 2 3 4 %.1f \")" % dihedral

to

dihedral_string =  "1 2 3 4 %.1f" % dihedral

and

fixed_dihedral = (" 1 2 3 4 $dihedral ")

to

fixed_dihedral = $dihedral_string

will resolve the issue you were running into. However, I’ve run into a number of problems trying to use optking for this purpose. Some are easily fixable (bad true/false in some if statements), but others are not. For example, it seems that if psi4 can’t determine the point group of the molecule at some rotations, it will simply not change the geometry and the optimization will fail. I don’t think this part of the code is ready for production use yet (but I would happy to be told I’m wrong).

Best,
Nate

I did PES scan for dihedrals in psi4. For this, first I rotated the dihedral which I want to calculate then I calculated single point energy for each step. I compared “relative energies” and I saw that Relative energies was equal in gaussian09 and psi4. I tested many times and I got the same results.