Chemistry matters. Join us to get the news you need.

If you have an ACS member number, please enter it here so we can link this account to your membership. (optional)

ACS values your privacy. By submitting your information, you are gaining access to C&EN and subscribing to our weekly newsletter. We use the information you provide to make your reading experience better, and we will never sell your data to third party members.


Computational Chemistry

Density functional theory error discovered

Free energy results can vary significantly with the computational method, but a simple fix reduces errors

by Sam Lemonick
July 29, 2019


Credit: Steven Wheeler
Density functional theory can give a range of different results when calculating the relative free energies of 1,3-butadiene and 2-butyne when users choose a coarse integration grid (top, black). Using a finer grid (bottom) reduces that variation in results.

Chemists use density functional theory (DFT) to accurately approximate the exact properties—like free energies—of molecules or materials in a reasonable amount of time. The method is widely used even outside of computational chemistry groups, often to back up or explain experimental results from the lab. But a new study suggests that some researchers, using standard settings for routine DFT computations, may be unknowingly introducing significant errors into their free-energy calculations (ChemRxiv 2019, DOI: 10.26434/chemrxiv.8864204.v4).

As a method that relies on approximations, DFT has a number of known errors and fixes for them. But Kendall N. Houk, a computational chemist at the University of California, Los Angeles, who was not involved in the new work, says, “This paper identifies a problem that many of us did not know about.”

Steven Wheeler of the University of Georgia and former group-member Andrea N. Bootsma, now a senior scientist at Pfizer, discovered the problem when Bootsma was testing an automated computer program that uses DFT to calculate transition-state energies as part of a screen for new catalysts. When she double-checked the computer’s calculations by hand, the results didn’t match. After rejecting other possible sources of error, she and Wheeler realized they had found a flaw inherent to DFT.

Bootsma and Wheeler found that DFT calculations of free energy can vary by as much as 5 kcal/mol depending on the orientation of the molecule, and these calculations shouldn’t depend on orientation, Wheeler says. Calculating free energy can be useful for predicting a reaction’s stereoselectivity or regioselectivity and a 5 kcal/mol error can mean the difference between producing one isomer versus another. Fortunately, the researchers also found that users can reduce these errors by choosing the right parameters for calculations.

The researchers determined that the errors are related to the size of the integration grid, which DFT algorithms use to divide up a molecule’s 3-D electron density into smaller, easier-to-handle densities to approximate the overall energy of the system. They say a (75,302) grid, which is a relatively coarse grid and is the default setting in some commonly used DFT programs, is insufficient. They recommend that users choose a finer grid, such as (99,590) at minimum.

Wheeler explains that the error is a result of how DFT handles molecules with low-frequency bending modes. In other words, floppy molecules. He and Bootsma used DFT to calculate the free energies of 2-butyne and several other molecules with low-frequency modes. They found free energy variances of up to 5 kcal/mol using the coarsest grids, while using finer grids reduced the differences to less than 1 kcal/mol in most cases.

Experts know that molecules with low-frequency modes can pose problems for DFT, but while DFT manuals may warn users that grid size can affect calculations on floppy molecules, “I’m not sure those warnings are always heeded,” says computational chemist John Herbert of the Ohio State University. Filipp Furche, a computational chemist at the University of California, Irvine, says users need to educate themselves about the methods they use, but he acknowledges that developers have a responsibility as well to educate users.

Chemists say they think these kinds of errors have probably made it into many papers. Daniel Singleton, an organic chemist at Texas A&M University, says since reading the preprint, he has recalculated some of his own results and found one instance in a published paper in which the free energy changed by more than 1 kcal/mol depending on its orientation, although he says the discrepancy doesn’t affect the paper’s conclusions (J. Am. Chem. Soc. 2015, DOI: 10.1021/ja5111392).

Singleton says, while many chemists may have published similar errors unknowingly, there’s also a possibility that some chemists re-ran their calculations using different settings until they got a free energy that agreed with their experimental results, even if they didn’t understand why the computed value changed. “Chemists don’t see themselves as subject to the same kinds of issues that other fields have seen, for example with the ‘replication crisis,’ but the bias of humans toward good results is universal,” Singleton says.

Bootsma wonders whether some researchers fell victim to this orientation problem and then abandoned a project because their computational and experimental results didn’t match.

She says the next step is for chemists to decide how best to address this issue, which is why she and Wheeler decided to publish this preprint. They wanted to raise awareness about the problem and start a conversation about how to handle it.


This story was updated on Aug. 6, 2019, to clarify Daniel Singleton's comment. He thinks chemists may have used different settings, not different orientations, to find DFT results that agreed with experiment.



This article has been sent to the following recipient:

Hans Senn (August 1, 2019 7:40 AM)
I would strongly disagree with the characterisation of the described problem as "flaw inherent to DFT": the problem is entirely a numerical one. As ever in computational science, the translation of a, perhaps even exact, theoretical model into a practical, computationally efficient method involves a plethora of approximations. The accuracy of these approximations is normally controlled by some parameters that need to be carefully chosen by the developer or the user, and the right choice is normally system- and application-dependent. So far so unsurprising.

For the specific case of the reliability/robustness of free energies of "floppy" molecules calculated within the harmonic-oscillator model using DFT, two quite separate issues need to be distinguished:

1) The numerical problems introduced by the integration grid affect mostly low vibrational frequencies. These low frequencies contribute substantially to the vibrational entropy, which is used to compute the free energy. However, the frequencies are obtained within the harmonic-oscillator model, which is invalid for low-frequency vibrational modes. Such modes are significantly anharmonic or even correspond to hindered rotations, rather than vibrations. So notwithstanding their numerical robustness, these frequencies have little physical meaning, but they unduly affect the result.
This is a well-understood and widely known limitation of the standard quantum-chemical approach to calculating free energies. Various solutions have been proposed to eliminate, or at least alleviate, this problem. A simple, effective approach, proposed by Truhlar et al., is to raise all vibrational frequencies below a certain threshold (e.g., 100 cm^–1) up to that threshold. This removes the dependence of the free energy on the exact numerical values of low frequencies, which are physically ill-justified anyway.
Moreover, when things are really "floppy", e.g., two or more initially separate molecules associating into a transition state, the free-energy minimum may differ significantly from the potential-energy minimum. The basic underlying assumption of "static" computational thermochemistry that thermodynamic quantities can be calculated as "single-point corrections" at the potential-energy minimum no longer holds in such cases.

2) Purely numerical accuracy and robustness is quite a separate matter. Irrespective of the physical validity of the approach, the user is entitled to expect that the values obtained (here, the low harmonic frequencies) are numerically accurate and robust. I would therefore agree to this limited extent only that the original article has highlighted an issue that developers should try to address, e.g., by appropriate warnings or automatic tests.
Klavs Hansen (August 6, 2019 12:53 AM)
Hans Senn makes some very good points. For one, application of any code that relies on a discretization must be run with different values to check for convergence. This ought to be 101 of numerical work.
Davide Donadio (September 9, 2019 12:42 PM)
I totally support Hans' comment. This article is misleading. It reports a technical flaw, which is not at inherent to DFT, but rather to default settings of certain specific codes.

Leave A Comment

*Required to comment