Those who are well-versed in the computational arts won't get much out of this post. On the other hand, happy amateurs like myself might find it a bit more useful as it clarifies something that's been bothering me.
Short version: Hcorr in Nwchem and Gaussian include PV.
As usual, I am not an expert in either computational or theoretical chemistry. I try to use it as a tool, and I try to use it as well as I can. But I am not an authority. Also, if you consult older posts on the blog you'll find plenty of examples of me misunderstanding basic computational concepts (with 550+ posts it's difficult to go back and erase all the embarrassing little gems).
The background
I had a bit of a fright the other day. I'm currently working on computationally characterising a system which undergoes a reaction that can lead to a large number of isomers. Only one of them is experimentally observed, and so it was of interest to see whether this is the thermodynamically favoured product as predict by reasonably cheap computational methods.
Because DFT calculations like these are based on gas phase reactions (even if you use a solvation model) the free energy that you get is based on the standard conditions for gas phase corrections i.e. 1 bar of partial pressure of each reactant.
If you want to calculate the free energy of reaction in solution you need to use concentrations instead. As you will see, you'll only have to worry about this for reactions that involve an unequal number of reactants and products. Normally your best results will be obtained by using isodesmic reaction schemes anyway, which is a great way of avoiding this. However, if you do have unequal numbers of reactants and products you /must/ correct for it when making solution phase predictions.
A gas at 1 bar of pressure is ca
So for an example reaction like this:
Using
In other words: you can't ignore the standard state change when doing solution phase calculations. This is obviously of extra importance in pH calculations, which are notoriously tricky.
Enthalpy
So knowing the need for standard state corrections I was a bit paranoid about how to treat the reaction enthalpy and came across this document: http://chemistry.illinoisstate.edu/standard/che38037/handouts/380.37assign3.pdf
On page, equation 2 states that for Gaussian
(where Ee is the electronic energy) which is an expression that doesn't include Δ(PV). In that case
Consulting the excellent Gaussian thermochemistry whitepaper (http://www.gaussian.com/g_whitepap/thermo.htm) offers the following:
At a first glance, kBT isn't equivalent to RT, but in fact kB=R/NA -- in the words of Wikipedia: "[The gas constant R] is equivalent to the Boltzmann constant, but expressed in units of energy".
In other words, Δ(PV) is already accounted for in Hcorr in gaussian.
Somewhat clearer: the Freq() page, http://www.gaussian.com/g_tech/g_ur/k_freq.htm, on the gaussian website now states
We test that claim:
E=-527.489751 Hartree
RT=(1.987*298.15/1000)/627.503=9.4410e-04 Hartree
E+RT=-527.489751+9.4410e-04=-527.488807
It holds up.
A similar example from an nwchem calculation:
Nwchem also includes Δ(PV) in the thermal correction to enthalpy.
Short version: Hcorr in Nwchem and Gaussian include PV.
As usual, I am not an expert in either computational or theoretical chemistry. I try to use it as a tool, and I try to use it as well as I can. But I am not an authority. Also, if you consult older posts on the blog you'll find plenty of examples of me misunderstanding basic computational concepts (with 550+ posts it's difficult to go back and erase all the embarrassing little gems).
The background
I had a bit of a fright the other day. I'm currently working on computationally characterising a system which undergoes a reaction that can lead to a large number of isomers. Only one of them is experimentally observed, and so it was of interest to see whether this is the thermodynamically favoured product as predict by reasonably cheap computational methods.
Because DFT calculations like these are based on gas phase reactions (even if you use a solvation model) the free energy that you get is based on the standard conditions for gas phase corrections i.e. 1 bar of partial pressure of each reactant.
If you want to calculate the free energy of reaction in solution you need to use concentrations instead. As you will see, you'll only have to worry about this for reactions that involve an unequal number of reactants and products. Normally your best results will be obtained by using isodesmic reaction schemes anyway, which is a great way of avoiding this. However, if you do have unequal numbers of reactants and products you /must/ correct for it when making solution phase predictions.
A gas at 1 bar of pressure is ca
V=n*R*T/P= 1*8.314*298.15/101350=0.024458 m^3= 24.46 litres
So for an example reaction like this:
A+B -> C
Using
G_SATP=G(gas phase) + R*T*ln Q
=G(gas phase)+ R*T*ln([C]/([A]*[B])
=G(gas phase)+ R*T*ln((1/24.46)/((1/24.46)*(1/24.46)))
=G(gas phase)+1.89 kcal/mol
In other words: you can't ignore the standard state change when doing solution phase calculations. This is obviously of extra importance in pH calculations, which are notoriously tricky.
Enthalpy
So knowing the need for standard state corrections I was a bit paranoid about how to treat the reaction enthalpy and came across this document: http://chemistry.illinoisstate.edu/standard/che38037/handouts/380.37assign3.pdf
On page, equation 2 states that for Gaussian
Ee+Hcorr=Ee+Hvib+Hrot+Htrans
(where Ee is the electronic energy) which is an expression that doesn't include Δ(PV). In that case
Δ(H)=Σ(Hproducts)-Σ(Hreactants)+Δn*R*Twhere Δn is the difference in number of moles of products and reactants.
Consulting the excellent Gaussian thermochemistry whitepaper (http://www.gaussian.com/g_whitepap/thermo.htm) offers the following:
Hcorr=Etot+kBTand
The Gibbs free energy includes Δ(PV) , so when it's applied to calculating ΔG for a reaction, ΔNRT=ΔPV is already included. This means that ΔG will be computed correctly when the number of moles of gas changes during the course of a reaction.[Note that H=Ee+Hcorr=Ee+Etot+kBT]
At a first glance, kBT isn't equivalent to RT, but in fact kB=R/NA -- in the words of Wikipedia: "[The gas constant R] is equivalent to the Boltzmann constant, but expressed in units of energy".
In other words, Δ(PV) is already accounted for in Hcorr in gaussian.
Somewhat clearer: the Freq() page, http://www.gaussian.com/g_tech/g_ur/k_freq.htm, on the gaussian website now states
Sum of electronic and zero-point Energies=-527.492585 E0=Eelec+ZPE
Sum of electronic and thermal Energies= -527.489751 E= E0+ Evib+ Erot+Etrans
Sum of electronic and thermal Enthalpies=-527.488807 H=E+RT
We test that claim:
E=-527.489751 Hartree
RT=(1.987*298.15/1000)/627.503=9.4410e-04 Hartree
E+RT=-527.489751+9.4410e-04=-527.488807
It holds up.
A similar example from an nwchem calculation:
Zero-Point correction to Energy = 120.416 kcal/mol ( 0.191895 au)(0.203513-0.202569)=9.4400e-04 Hartree(1.987*298.15/1000)/627.503
Thermal correction to Energy = 127.114 kcal/mol ( 0.202569 au)
Thermal correction to Enthalpy = 127.706 kcal/mol ( 0.203513 au)
Nwchem also includes Δ(PV) in the thermal correction to enthalpy.
This is dummy text. It is not meant to be read. Accordingly, it is difficult to figure out when to end it. But then, this is dummy text. It is not meant to be read. Period.
ConversionConversion EmoticonEmoticon