Adding Anisotropy to the Standard Quasi-Harmonic Approximation Still Fails in Several Ways to Capture Organic Crystal Thermodynamics

We evaluate the accuracy of varying thermal expansion models for the quasi-harmonic approximation (QHA) relative to molecular dynamics (MD) for 10 sets of enantiotropic organic polymorphs. Relative to experiment we find that MD, using an off-the-shelf point charge potential, gets the sign of the enthalpic contributions correct for 6 of the 10 pairs of polymorphs and the sign of the entropic contributions correct for all pairs. We find that anisotropic QHA provides little improvement to the error in free energy differences from MD relative to isotropic QHA, but does a better job capturing the thermal expansion of the crystals. A form of entropy-enthalpy compensation allows the free energy differences of QHA to deviate less than 0.1 kcal/mol from MD for most polymorphic pairs, despite errors up to 0.4 kcal/mol in the entropy and enthalpy. Deviations in the free energy of QHA and MD do not clearly correlate with molecular flexibility, clarifying a previously published finding. Much of the error previously found between QHA and MD for these flexible molecules is reduced when QHA is run from a lattice minimum consistent with the same basin as MD, rather than the energy-minimized experimental crystal structure. Specifically, performing anisotropic QHA on lattice minimum quenched from low-temperature replica exchange simulations reduced the error previously found by 0.2 kcal/mol on average. However, these conformationally flexible molecules can have many low-temperature conformational minima, and the choice of an inconsistent minima causes free energies estimated from QHA to deviate from MD at temperatures as low as 10 K. We also find finite size errors in the polymorph free energy differences using anisotropic QHA, with free energy differences as large as 0.S kcal/mol between unit and supercells loosely correlated with differences in anisotropic thermal expansion. These larger system sizes are computationally more accessible using our cheaper 1D variant of anisotropic QHA, which gives free energies within 0.02 kcal/mol of the fully anisotropic approach at all temperatures studied. The errors between MD and experiment are 1-2 orders of magnitude larger than those seen between QHA and MD, so the quality of the force field used is still of primary concern, but this study illustrates a number of other important factors that must be considered to obtain quantitative organic crystal thermodynamics.