QuantumATK Forum
QuantumATK => General Questions and Answers => Topic started by: k.cerda on October 7, 2021, 13:27
-
Hello
The total energy of a molecule is negative when calculated using TotalEnergy() function: https://docs.quantumatk.com/manual/Types/TotalEnergy/TotalEnergy.html
When I calculate the free (vibrational) energy of the same molecule it is a positive number using the freeEnergy() function: https://docs.quantumatk.com/manual/Types/PhononDensityOfStates/PhononDensityOfStates.html
I want to calculate the Gibbs free energy therefore should I add the total energy obtained from TotalEnergy() and the vibrational free energy obtained from phonon_density_of_states.freeEnergy()? E.g. if the TotalEnergy() gives -20eV and freeEnergy() gives +1eV, is it true that Gibbs free energy of the molecule is -20+1=-19eV?
Thanks a lot.
-
Anyone here?
-
Hi Cerda,
For molecules, we have a new analysis object 'IdealGasThermochemistry' that takes a molecule configuration and its dynamical matrix object, similar to PhononDensityOfStates, to provide the Gibbs free energy including the electronic energy https://docs.quantumatk.com/manual/Types/IdealGasThermoChemistry/IdealGasThermoChemistry.html
Example:
In GUI: open scripter, under 'configurations' select 'Analysis from file'. Load the hdf5 file with the molecule configuration and dynamical matrix and select them. Add the 'IdealGasThermoChemistry' analysis object to the scripter and submit the script. After the calculation, choose the 'IdealGasThermoChemistry' object from the hdf5 file and choose text representation to get the list of Gibbs free energy at different temperatures.
PhononDensityOfStates is more relevant for bulk/slab systems as we do not have to compute the rotational and translational free energy contributions. For that, we have the 'CrystalThermoChemistry' analysis object similar to the case of the molecules above.
With that said, you are correct in adding the free energy from phononDOS to the total electronic energy to get the final free energy - but this will not include the rotational and translational contributions as you computed only vibrations for the molecule. Please use the 'IdealGasThermoChemistry' analysis object instead to get the Gibbs free energy.
I hope this answered your question.
-
Dear Shuresko
Thanks a lot for your reply. However the IdealGasThermo object gives the following error when applied:
+------------------------------------------------------------------------------+
| Ideal Gas Thermo Chemistry Report (at standard pressure 1 bar) |
| Configuration: C398 H579 |
| Total energy: -6224.1 eV |
| Zero point energy: 345.31 eV |
| ---------------------------------------------------------------------------- |
zipdir\NL\Analysis\IdealGasThermoChemistry.py:473: RuntimeWarning: invalid value encountered in log
| temperature (K) entropy (meV/K) free energy (eV) |
| 0.00 0.000000 -5878.83940665 |
| 100.00 nan nan |
| 200.00 nan nan |
| 250.00 nan nan |
| 298.15 nan nan |
| 300.00 nan nan |
| 350.00 nan nan |
| 400.00 nan nan |
| 450.00 nan nan |
| 500.00 nan nan |
| 600.00 nan nan |
| 700.00 nan nan |
| 800.00 nan nan |
| 900.00 nan nan |
| 1000.00 nan nan |
| 1100.00 nan nan |
| 1200.00 nan nan |
| 1300.00 nan nan |
| 1400.00 nan nan |
| 1500.00 nan nan |
| 1600.00 nan nan |
| 1700.00 nan nan |
| 1800.00 nan nan |
| 1900.00 nan nan |
| 2000.00 nan nan |
| 2100.00 nan nan |
| 2200.00 nan nan |
| 2300.00 nan nan |
| 2400.00 nan nan |
| 2500.00 nan nan |
| 2600.00 nan nan |
| 2700.00 nan nan |
| 2800.00 nan nan |
| 2900.00 nan nan |
| 3000.00 nan nan |
| 3100.00 nan nan |
| 3200.00 nan nan |
| 3300.00 nan nan |
| 3400.00 nan nan |
| 3500.00 nan nan |
| 3600.00 nan nan |
| 3700.00 nan nan |
| 3800.00 nan nan |
| 3900.00 nan nan |
| 4000.00 nan nan |
| 4100.00 nan nan |
| 4200.00 nan nan |
| 4300.00 nan nan |
| 4400.00 nan nan |
| 4500.00 nan nan |
| 4600.00 nan nan |
| 4700.00 nan nan |
| 4800.00 nan nan |
| 4900.00 nan nan |
| 5000.00 nan nan |
| 5100.00 nan nan |
| 5200.00 nan nan |
| 5300.00 nan nan |
| 5400.00 nan nan |
| 5500.00 nan nan |
| 5600.00 nan nan |
| 5700.00 nan nan |
| 5800.00 nan nan |
| 5900.00 nan nan
-
Hi Cerda,
Can you share the script used to optimize the molecule and compute the dynamical matrix so that I can check the settings? If I have to guess, something wrong with the dynamical matrix calculation - because the zero-point energy is unusually large.
-
Dear shuresko
Thanks again for your reply. I attached the script I used for the dynamic matrix calculation and the IdealGasThermoChemistry object. Looking forward for your valuable ideas. Thanks again.
-
Hi Cerda,
I see that you have used a FF calculator, which is fine, but the geometry should be optimized to reach the nearest minimum in the potential energy surface before performing the dynamical matrix calculation. Without geometry optimization, your script is generating a lot of -ve normal mode frequencies. Given the size of the configuration, this could take thousands of optimization steps before the forces go below 0.01 eV/A. A good practice before computing the free energy is to check the vibrational frequencies using the 'VibrationalMode' analysis object. For nonlinear molecules, all eigenfrequencies except the first 6 should be positive. If not, then either the FF is not good enough or the geometry is not a local minimum to perform harmonic analysis.
Please try geometry optimization until the vibrational frequencies are as mentioned above before performing IdealGasThermoChemistry - negative frequencies result in nan during vibrational entropy calculation under ideal gas approximation. When you use PhononDOS, then we only consider the positive frequencies for calculation, that is why it does not give you nan. But as I mentioned earlier, PhononDOS uses Quasi Harmonic Approximation meant for solids and does not include contributions for rotation and translation. Also, it is a constant volume approximation - meaning that it only gives you the Helmholtz free energy and not the Gibbs free energy.
Let me know if there are further questions.
-
Dear shuresko
Thanks for your detailed explanations. I'll try optimization again...