Dihedral Angle for Ethylene Glycol

Good Morning,
I am new in using Psi4 and I have some problems with the potential energies for the dihedral angle of the Ethylene Glycol.
I ran different codes (changing some parameters) but the Potential Energies are always the same for all the angles (more or less). My purpose is to plot these energies to find the different conformers of this molecule.

Here is one of the codes that I used.
It will be a pleasure for me if someone can help me and explain me why the energies are all the same.


molecule Ethylene_Glycol{
0 1
H
C  1  r2
H  2  r3  1  a3
O  2  r4  1  a4  3  d4
H  4  r5  2  a5  1  d5
C  2  r6  1  a6  3  d6
H  6  r7  2  a7  1  d7
H  6  r8  2  a8  1  d8
O  6  r9  2  a9  1  d9
H  9  r10  6  a10  2  d10

r2= 1.1107
r3= 1.1139
a3= 107.32
r4= 1.4014
a4= 108.51
d4= 119.25
r5= 0.9929
a5= 106.65
d5= 198.10
r6= 1.5276
a6= 109.46
d6= 239.75
r7= 1.1107
a7= 109.46
d7= 299.05
r8= 1.1139
a8= 110.75
d8= 180.92
r9= 1.4014
a9= 110.37
d9=  58.41
r10= 0.9929
a10= 106.65
d10=  78.17
}
#optimize molecule at different dihedral angles
steps = 37
step_size = 10
start_angle = 0
ecp = {}
for counter in range (0, steps):
    # set smaller basis set for geometry HT geometry optimization
    set {
        basis sto-3g
        geom_maxiter 100
        scf_type df
        guess sad
        ints_tolerance 1.0E-8
    }

    dihedral = 1.0 * counter * step_size + start_angle
    dihedral_string =  "2 4 6 9%.1f" % dihedral
    set optking {
        frozen_dihedral = ""
        fixed_dihedral = $dihedral_string
    }
    optimize('scf')   # geometry optimization with HF

#change basis set for single point energy calculation
    set {
        basis cc-pVDZ 
        scf_type df
        guess sad
    }

    ecp[counter] = energy('scf')
psi4.print_out("\n")
psi4.print_out("Energies for different dihedral angles with sto-3g \n\n")
psi4.print_out("        Dihedral [degree]         E_int [kcal/mol]             \n")
psi4.print_out("-----------------------------------------------------\n")
for counter in range (0, steps):
   dihedral = 1.0 * counter * step_size + start_angle
   e = ecp[counter] * psi_hartree2kcalmol
   psi4.print_out("        %3.1f            %10.6f\n" % (dihedral, e)) 

Thank you so much!

Could you edit your post so the input file is all in backticks? The forum obliterates all your code formatting unless you do that, and I need that formatting to know what code is in your for loop and what isn’t.

Also, you can easily tell from the output file whether the problem is

  • The optimizations converge to the same thing despite you fixing a dihedral
  • The energy call isn’t using the geometry from the optimizations
  • Despite each energy call starting from a different geometry, they all return compute approximately the same energy
  • The energy computed by the energy call is different, but the energies printed out are the same

Posting information like that is extremely helpful for anybody who takes a look at this.

This the final part of my output and although the angle changes, the energy converges to the same value.
Also the geometry in all the steps does not change.

Thank you for your help.

        Dihedral [degree]         E_int [kcal/mol]             
-----------------------------------------------------
        0.0            -143662.446708
        10.0            -143662.460359
        20.0            -143662.452544
        30.0            -143662.452544
        40.0            -143662.452544
        50.0            -143662.452544
        60.0            -143662.452544
        70.0            -143662.452544
        80.0            -143662.452544
        90.0            -143662.452544
        100.0            -143662.452544
        110.0            -143662.452544
        120.0            -143662.452544
        130.0            -143662.452544
        140.0            -143662.452544
        150.0            -143662.452544
        160.0            -143662.452544
        170.0            -143662.452544
        180.0            -143662.452544
        190.0            -143662.452544
        200.0            -143662.452544
        210.0            -143662.452544
        220.0            -143662.452544
        230.0            -143662.452544
        240.0            -143662.452544
        250.0            -143662.452544
        260.0            -143662.452544
        270.0            -143662.452544
        280.0            -143662.452544
        290.0            -143662.452544
        300.0            -143662.452544
        310.0            -143662.452544
        320.0            -143662.452544
        330.0            -143662.452544
        340.0            -143662.452544
        350.0            -143662.452544
        360.0            -143662.452544

There’s a tiny typo in your input file. Change "2 4 6 9%.1f" to "2 4 6 9 %.1f".

I tried the input file with that change, and the first 10 optimizations went smoothly (with different energies!), but the 11th optimization doesn’t converge. My first piece of trouble-shooting advice is to take hessians every few steps. Because this is RHF, that can be done analytically, so the additional computation cost should be relatively minor.