Unreasonable geometry of E-stillbene with DFT

The optimisation of the geometry of E-stillbene, combined with frequency calculation, was performed using different level of theory with QChem (wB97X-D3, wB97X-D and M06-2X as functionals and 6-311G* as basis set). Unfortunately, the results were not as foreseen. A planar geometry is expected for a single molecule of E-stillbene in a gaseous phase.

When using the functional wB97X-D3, a flat geometry with imaginary frequencies were obtained for a starting flat geometry. After slightly twisting the optimised geometry, a twisted optimised structure is obtained (but no imaginary frequency). The same happened with both wB97X-D and M06-2X functionals.

On the other hand, when using the software Turbomole with the functional wB97X-D, a flat optimised geometry is obtained from a flat initial geometry (and with real frequencies).

To get rid of these imaginary frequencies, different approaches were tested, such as pre-optimising with another functional and using the optimised structure as input, modifying bond length (up to .5 Ångström) and/or dihedral angles (up to 5 degrees).

Did anyone encounter these kind of problems and know (a) how to get rid of these frequencies, (b) how to get a reasonable flat geometry ?

Thank you in advance for your help,
R. Crits and M. Maraldi

Does this persist with XC_GRID = 000099000590, i.e., (99,590) EML grid?

Yes it does ! We used mostly this grid and also tried others.

Another difference vs. Turbomole (if I remember correctly) is the default setting of pure vs. Cartesian Gaussians, what Q-Chem calls PURCAR. (You should be able to confirm by seeing how many basis functions each program tells you it is using.) I would be surprised if this made a difference, but is something to check. How large is the imaginary frequency? Have you tried tightening Q-Chem optimization stopping criteria? (GEOM_OPT_TOL_* where * = ENERGY, DISPLACEMENT, or GRADIENT. Defaults are 100, 1200, 300, respectively, in some multiples of a.u., I usually divide all those by 2 if I have a question.)

1 Like

Update: I can verify the imag freq with wB97X-D/6-311G* and a lot of variations in thresholds, etc. (using Q-Chem 6.1). Not sure what the issue is.

Update #2: also confirmed with numerical frequencies (finite diff of analytic gradients). Is it possible that the Turbmole result is wrong and that the molecule isn’t quite planar?

With wB97X-D, frequencies are 37.46 i and 17.42 i and with M06-2X they are 50.00 i and 2.21 i. After slightly modifying the geometries, the frequencies become 34.39 i and 17.84 i with wB97X-D, 28.24 i and 10.00 i with M06-2X.

Regarding the optimization criteria, we tried GEOM_OPT_TOL_ENERGY=1 and GEOM_OPT_TOL_GRADIENT=1. A highly distorted geometry was obtained.

This is the geometry we got with wB97X-D/6-311G* using Turbomole.

Reducing the stopping criteria all the way to 1 is a bad idea, that’s unrealistically tight in view of other thresholds; my suggestion was to reduce by 2x. Anyway I tried that, along with tightening various thresholds. Moreover that imaginary frequency persists in a finite difference calculation. So my question is: why are you convinced that Turbomole’s planar geometry is the right answer?

By reading in the Hessian from the imaginary freq calculation,

$rem
jobtype		opt
scf_guess	read
method		wb97x-d
basis		6-311G*
xc_grid		000099000590
thresh		16
mem_total	110000
geom_opt_tol_energy 50
geom_opt_tol_displacement 600
geom_opt_tol_gradient 150
geom_opt_hessian 	read
purcar		2222
sym_ignore	true
$end

I was able to optimize a structure with all real freqs. It has a dihedral angle of about 45 deg between the phenyl rings.

We are in gas phase at 0 K, so we expect E-stillbene to be planar because it is highly conjugated. Hereafter are papers that confirm our statement:
https://www.sciencedirect.com/science/article/abs/pii/S0009261406001369?via%3Dihub
http://www.unige.ch/sciences/chifi/publis/refs_pdf/ref00861.pdf

I’m not disputing that, but the wB97X-D/6-311G* model chemistry seems to want to give it a twist. Experiments with thresholds and with finite-difference frequencies don’t suggest any problems with the minimum-energy structure.

@jherbert suggestion with using the exact Hessian for the start of the geometry is a good one to assist with finding the minimum. A more costly optimization would be to use exact Hessian at every step (algorithm = Newton) or every few steps, see Example 9.5.

In addition, you could run a PES scan around the dihedral angle to see the energy profile of the twist - planar - other twist.