I have a question regarding energy calculations. I understand that force fields are more for the purpose of cleaning up the geometry of a molecule, and are not particularly accurate when it comes to energy calculations. However, I have compiled a series of adefovir analogues, and am collecting data and identifying differences and similarities relative to Adefovir, such as bond lengths for example. I would like to do the same with the energy optimisation tool, but am wondering what the best way to do it is, and which force field to use. I have input adefovir to begin with and optimised its energy using the Universal Force Field, and got an answer of 563.7. I then did the same thing with MMFF94 and got the answer -237.6. I am wondering why there is such a drastic difference between answers, and which one would be better to use for the other molecules. I completely understand that quantum mechanical methods would be much more accurate, however, for the sake of this project, the goal is just to make very broad conclusions about this group of molecules in comparison to each other by comparing these energy values. As long as each figure is calculated in the same way the accuracy is not of huge importance. It is just where molecules rank compared to each other. I would really appreciate any help.
One technical problem with adefovir (and analogues) is conformational flexibility. This is to say in presence of e.g., alkyl chains and other groups which may rotate (at least a little) you have to check if the energetic minimum of the optimization not only is a local one (as for the two terminal methyl groups in gauche relationship in n-butane), but the global one (the two terminal methyl groups of butane in anti relationship) when comparing molecule A with molecule B. Openbabel can achieved this by a systematic rotor search, and Avogadro2 equally provides this. For each of these conformers, an energy can be assigned (hence multiple 3D conformer on pubchem, too.) This however is for an isolated molecule in vacuum – once it binds to a receptor, an other conformation can be more favourable; and you have to compare the same conformation of molecule A (while bound) with molecule B (while bound).
A second point to consider is the force field used. For different structure motifs the authors considered useful, the authors established parameters. The U in UFF is for universal in the authors aimed for a parametrization to cover every atom in every geometry (known at time of publication). MMFF94 however was parametrized for small organic molecules (as then relevant to the researchers at Merck Rahway/NJ). I recommend to check if P in general, and organic phosphates (as motif P is embedded) in particular are covered, or not.
In general, MMFF94 or GAFF are parameterized for small drug-like molecules, while UFF as indicated is intended as a universal force field (e.g., across the periodic table).
Force fields will give very different energies because they tabulate different things. For example, MMFF94 includes electrostatics and hydrogen-bond terms, but UFF does not.
But I’m not sure what you want to do with the energies? Or bond lengths?
My project focuses on a series of adefovir analogues where the main structural variation is ortho, meta, and para substitution patterns relative to the core scaffold. Preliminary results suggest that ortho isomers tend to have higher minimized energies than para isomers, likely due to greater steric strain and less favorable spatial arrangement.
Regarding bond lenghts, for each molecule, I’m measuring the distance between the phosphorus atom in the phosphonate group and the N9 nitrogen in the adenine base, as an indicator of how ring connectivity and substitution influence the overall molecular dimensions. The idea is not just to obtain the “perfect” structure or absolute energy, but to illustrate what trends and geometric data can be rapidly generated using 3D modeling software like Avogadro, and how these can inform structure–property relationships or guide future analogue design—especially in a first-pass, in silico setting like mine.
I’m not synthesizing these compounds, my goal is to demonstrate what can be learned about a chemical series using computational tools, understanding the limitations of force field-based energies. I’ll be reporting trends—not absolute values—and making it very clear that this is just a first-pass computational analysis. If you have recommendations for workflow steps, force field choice, or any other advice I’d greatly appreciate it.
I would highly recommend optimizing geometries and calculating energy trends with GFN2-XTB. For example, you can download the xtb plugin with Avogadro and this should help install xtb for optimizing geometries, calculating energies and make it easier to run.
In general, GFN2 energies and geometries are significantly more accurate than traditional force fields.