crest-lab / crest

CREST - A program for the automated exploration of low-energy molecular chemical space.
https://crest-lab.github.io/crest-docs/
GNU Lesser General Public License v3.0
190 stars 42 forks source link

About using ORCA for DFT calculations in Calculation of absolute molecular entropies and heat capacities made simple paper #49

Closed xiki-tempula closed 3 years ago

xiki-tempula commented 3 years ago

I was reading the recent paper of "Calculation of absolute molecular entropies and heat capacities made simple" and the higher level DFT calculations are done using TURBOMOLE. Given that TURBOMOLE is commercial software and might not be available to anyone. I wonder if it is possible to use support ORCA? Since it also has the level of B97-3c and of course B3LYP-D3(BJ)/def2-TZVP with RIJCOSX support.

pprcht commented 3 years ago

You can use any QM program you want for the DFT part, that's independent of crest. The only thing to look out for is whether you are able to apply the frequency scaling factor and rotor cut-off in the thermostatistical calculation. xtb's and crest's thermo implementation can do that, but only in Turbomole's/xtb's vibspectrum format. You might even use other functionals if you wish to do so and determine appropriate frequency scaling factors.

xiki-tempula commented 3 years ago

@pprcht Thank you for the comment. I'm trying to use ORCA for this calculation.

Based on my understanding, the S = S_msRRHO + S_conf. The crest standalone could give the S_conf at GFN2-xTB or GFN-FF level (or other xtb supported level).

The S_RRHO is the original Total entropy correction scaled according to Truhlar (https://comp.chem.umn.edu/freqscale/version3b2.htm) and applying the Grimme's entropy interpolation between harmonic and free-rotor approximations.

I wonder if this is correct? I'm aware of the recent development of the Single Point Hessian but I cannot find a tool for applying it to orca.

pprcht commented 3 years ago

yes, crest gives you only S_conf , which is a combination of Gibbs-Boltzmann conformational entropy S'_conf and a |S|_msRRHO average calculated at some GFN level (i.e., S_conf = S'_conf + |S|_msRRHO). These terms are entirely included in the S_conf crest gives you, and no manual adjustments are required here.

For the absolute entropy (S = S_msRRHO + S_conf) you then need then S_msRRHO for the lowest conformer at DFT level, for which the frequencies are scaled and the interpolation is applied, as you already mentioned. The B97-3c ord B3LYP scaling factors we applied are in the Chem. Sci. paper rather than in Truhlars database because we newly determined them for this study. In ORCA the scaling factor can be easily modified, but I'm not completely sure sure about the interpolation. I know the interpolation is implemented in ORCA with a default threshold (35 cm-1). Apparently this can be modified in ORCA's %freq input block with the keyword CutOffFreq, but please look at the ORCA manual again. If that doesn't work, a workaround would be to take the frequencies from the ORCA calculation and write a vibspectrum file with it, which then could use xtb's or crest's thermo implementation.

This is unrelated to the Single Point Hessian procedure. SPH is used within the crest calculations, so it doesn't influence your DFT calculation.

xiki-tempula commented 3 years ago

@pprcht Thank you. I just want to check if my understanding is correct and S_msRRHO doesn't contain anything else. I usually use https://github.com/bobbypaton/GoodVibes or http://sobereva.com/soft/shermo/ for scaling and Grimme's entropy interpolation.

xiki-tempula commented 3 years ago

@pprcht Thanks for addressing this question.

Another short question. The output is:

 ------------------------------------------
 | FINAL MOLECULAR ENTROPY AT T= 298.15 K |
 ------------------------------------------
   Sconf   =          1.296712
 + δSrrho  =         -0.029233
 ----------------------------------------
 = S(total)  =          1.267479

   H(T)-H(0) =          0.126621
   G         =         -0.377899   (-T*S)
 ----------------------------------------
 = G(total)  =         -0.251278  (H-T*S)

   Cp(total) =          0.166227

So the S_conf would be = S(total) = 1.267479, I wonder what is the unit? Is it in in Hartree (T= 298.15 K) or cal/molK?

The unit is not explicitly stated and I think number seems to ba a big too big for hatree.

pprcht commented 3 years ago

The entropies and heat capacities are in cal/molK, the enthalpy and free energy are in kcal/mol. A bit further up there should be a note about that (below one of the tables with the temperature dependencies). I'll add a note for this printout again.