I’m trying to port over some of our gaussian teaching examples to psi4 and I am struggling with an open shell system - the methyl radical.

It will converge just fine to the correct trigonal bipyramidal geometry using gaussian.

%NProcShared=2
%Mem=1GB
#P MP2/aug-cc-pvdz Opt Symmetry=None
CH3.
0 2
C 0.749 0.000 0.000
H 1.112 0.000 1.028
H 1.112 -0.890 -0.514
H 1.112 0.890 -0.514

But not using psi4. I tried different methods and basis sets (dz/tz, non-aug), setting "g_convergence":"gau_loose", "maxiter":200, "step_type":"nr","full_hess_every":1, "geom_maxiter":200 individually and in combination but it won’t converge. The last geometry printed in the log file is also nowhere near trigonal byramidal.

The psi4 steps I am using are as follows:

ch3_radical = psi4.geometry("""
0 2
symmetry c1
C 0.749 0.000 0.000
H 1.112 0.000 1.028
H 1.112 -0.890 -0.514
H 1.112 0.890 -0.514
""")
psi4.set_options({"reference" :"uhf"})
psi4.optimize('mp2/aug-cc-pvdz', molecule=ch3_radical)

This is the result of a bug introduced between 1.4rc1 and 1.4rc2. It’s fixed on the current developer version, so won’t be there in the actual 1.4.

Options:

Depending on when you need this teaching example, use 1.4rc1 (already released) or 1.4 (not released yet, should be before fall, but I can get a sharper estimate if that’d help)

Change the theoretical method. In this case, it’s enough to set mp2_type conv. This changes MP2 from using the density-fitted approximation (default in Psi4), to exact two-electron integrals (which I believe is the default in Gaussian anyways)

Run the geometry optimization on one processor instead of in parallel

Thanks for your quick response.
I ran the same example using: mp2_type conv, psi4.set_num_threads(1) and using 1.4rc1.
The example still does not converge.

In addition to what I wrote in my previous post, you’ll want to set "opt_coordinates": "both".

If you know what to look for in the output file, you can tell that the geometry optimizer is having some difficulty converting between displacements in internal coordinates and displacements in Cartesian space. The command I posted tells the optimizer to treat Cartesians and internals on a more equal footing when calculating displacement. This is slower, but lets the system converge.

Paging @Rollin_King. This should be on your radar for the geometry optimizer rewrite, if it hasn’t been already.

@jmisiewicz , thanks for notice. FYI, this situation (in which dihedrals stink) has been fixed to use out-of-plane angles and optimize fine (without cartesians) in new optimizer. For example, this works:

ch3 = psi4.geometry("""
0 2
symmetry c1
C 0.749 0.000 0.000
H 1.112 0.000 1.028
H 1.112 -0.890 -0.514
H 1.112 0.890 -0.514
“”")
psi4options = {
‘reference’: ‘uhf’,
‘basis’: ‘aug-cc-pvdz’,
‘mp2_type’: ‘conv’
}
psi4.set_options(psi4options)
import optking
optking.optimize_psi4(‘mp2’)

Glad it worked! Using just Cartesians would probably work too, but rather on accident a few years ago I discovered that you can do pretty well just throwing all the internals you like plus all the Cartesians into your redundant set is pretty effective. Don’t know if anyone has benchmarked the approach.

In this case, it took me 17 iterations using just Cartesians, and 7 using both. I think using both is the best solution here. In any case, I’ve marked the issue as solved.