Efficient Determination of Force Field Parameters Using a Physically Based Equation of State

One of the most prominent problems when applying molecular simulations to represent real substances is the choice of (or the optimization of) the force field parametrization [1]. With ``parametrization'', we here refer to both the choice of a functional form (such as Lennard-Jones or Buckingham potentials for van der Waals interactions) and the values of its parameters. Usually, interactions can be divided into bonded and non-bonded interactions. Although the parametrization that describes the bonded interactions is accurately obtained from quantum chemical calculations or spectroscopic data [2, 3], it is difficult and computationally expensive to obtain accurate parameters describing the non-bonded van der Waals interactions from quantum mechanics; in practice, the non-bonded parameters are usually obtained by fitting molecular simulation results to experimental data like vapor-liquid equilibrium (VLE) data. Since the calculation of the vapor-liquid equilibria depends on the van der Waals interactions of all atom types in the system, the fitting of molecular simulation results to VLE data requires many time-consuming molecular simulations; therefore this is currently a cumbersome and computationally expensive task.

In this work, we present a novel method to significantly reduce the computational burden required for fitting force field parameters to experimental VLE data. In sharp contrast to conventional optimization approaches, we use additional physical insight on how computed VLE data depends on the force field parameters in molecular simulations. We use an analytical equation of state, the perturbed chain statistical associating fluid theory (PC-SAFT) [4], to predict the results of Gibbs Ensemble molecular simulations without actually performing these simulations. We show that all parameters in the PC-SAFT approach are directly related to the force field parameters, thereby making it possible to optimize force field parameters in a different, much faster framework (i.e. PC-SAFT). Since PC-SAFT calculates VLE data at least a factor 10⁶ faster than molecular simulation, the presented optimization procedure results in a drastic decrease of the computational burden.

As a proof of concept, the method is illustrated by optimizing transferable Lennard-Jones parameters for the n-alkane series. The target is to obtain united atom force field parameters that describe the vapor-liquid part of the phase diagrams of n-butane and n-hexane. Optimizing a vector p of 4 force field parameters p = (ε_CH2CH2, ε_CH3CH3, σ_CH2CH2, σ_CH3CH3) we obtain excellent agreement of coexisting densities, vapor pressure and caloric properties using only two Gibbs ensemble simulations. It is important to note that as we use additional physical insight in the optimization procedure, the number of molecular simulations needed hardly depends on the initial guess for the vector p, i.e. choosing initial force field parameters that deviate 40% from the optimal values leads to an almost equally efficient convergence.

References

[1] Frenkel, D.; Smit, B., Understanding Molecular Simulation: from Algorithms to Applications, 2^nd ed. Academic Press; San Diego, 2002.

[2] Allinger, N.L.; Yuh, Y.H.; Lii, J.H. J. Am. Chem. Soc. 1989, 111, 8551-8566.

[3] Martin, M.G.; Siepmann, J.I. J. Phys. Chem. B 1998, 102, 2569-2577.

[4] Gross, J.; Sadowski, G. Ind. Eng. Chem. Res. 2001, 40, 1244-1260.