Fit Lennard Jones Parameters

The Lennard-Jones (LJ) parameters are critical for predicting condense phase properties.¹² The data must contain sufficient information in order to fit the parameters. DFF provides two functionalities to handle this problem: fit molecular clusters (dimers) and to fit liquid properties.

Due to strong correlations among the charge, valence and LJ parameters, fitting LJ parameters should be decoupled from the fittings of charge and valence parameters. After the LJ parameters are changed, the fitting of valence parameters should be repeated, and the entire process should be done in multiple iterations. If, due to the strong polarizations in condense phase, the charge parameters obtained by fitting small molecule data need to be adjusted, which can be signaled when the LJ parameters in reasonable ranges cannot fit the targeted condense phase properties. The entire process should be repeated. The good thing is, the iteration converges quickly, based on multiple research results.

Fit molecular clusters (dimer)

DFF can determine nonbond parameters by fitting molecular clusters (dimers). A number of configurations of the cluster should be generated and calculated with relatively high-level QM calculations. DFF automatically freezes all parameters except the LJ parameters to fit the potential energies using the Levenberg-Marquardt method.

It should be noted that this approach might not work well for predicting physical properties of molecules in condensed phases. Although in principle the VDW parameters can be derived by fitting to QM data calculated for molecular clusters, this approach has the following limitations. a) it is very difficult to accurately describe weak VDW interactions using QM method, high level of theory and an extensive basis set must be used. b) the dimer is too small to represent molecular interactions in condensed phases. The result energetic data only represents molecular interactions in gas phase, which in many cases can be significantly different from those found in the condensed phases.

Fit Liquid Data

DFF can be used to optimize VDW parameters by combining a molecular dynamics simulation with a least-squares fit procedure.

The statistic mechanical perturbation theory is used to enable the automation. Assume a small change to the adjustable parameters:

The change to the potential energy can be regarded as a perturbation:

According the perturbation theory, the configuration integral of canonical ensemble consisting of N particles can be written as:

where the brackets indicate unperturbed (reference state) ensemble average. For any physical property A, we have:

For example, the internal energy and pressure of ensemble average are given by:

The above equations indicate that the change in potential energy and pressure with the perturbation can be estimated by averaging the perturbation energy and virial in the reference state obtained by NVT simulation. The derivatives of an ensemble averaged property <A> (energy or pressure) with respect to the parameter is then:

If the quantity can be analytically evaluated, such as the pressure and internal energy in canonic ensemble, the derivative can be obtained in one simulation. Subsequently, the parameters can be optimized by the least-square method automatically.

In order to obtain transferable parameters, several molecules containing the same atom types should be parameterized simultaneously. For example, the following table shows molecular liquids that can be used to parameterize non-valence parameters for alkanes:

The table lists temperature, heat of vaporization and equilibrium density of the relevant molecular liquids.

The uncertainty in ensemble average has an impact to the fitting quality. Large uncertainties associated with short simulation time lead to oscillation in the fitting parameters. In such situations, extending the simulation time is necessary.

For polar molecules with relatively large atomic partial charges, the electrostatic terms dominate the nonbond interactions. If the LJ parameters cannot be optimized to fit empirical data, the atomic charges may not be appropriate, even if the charges are derived by fitting QM electrostatic potentials. The reason is that the separation of electrostatic and VDW parameters is for convenience but not scientifically solid. As in principle, all intermolecular interactions are electrostatic in nature. Comparing the charge values with those in successful force fields might be helpful. If necessary, consider using restraints in the fit of the charge parameters.

Similarly, if the valence parameters are incorrect, it could be very difficult to optimize the LJ parameters. Even if a fit seems to be successful, the LJ parameters might not be physically meaningful. Careful examination of the intramolecular properties such as bond lengths, angles and dihedral angles is a key to solve this problem.

Temperature Dependent Parameters

Strictly speaking, force field is not potential energy function because of the empirical (non-polarizable) functions deployed. The nonbond interaction is found to be dependent on temperature because of the polarization effects.³ One solution is to scale the dispersion coefficient C_6 using the following form of Lennard-Jones (LJ-12-6) function

where the scaling factor is temperature dependent, which can be written as a linear function,

The reference temperature T_ref is conveniently fixed at 298K.

Using the expression of Lennard-Jones (LJ-12-6) function, the well-depth and radius parameters are expressed as:

It is found that the scaling factor is rather generic, it takes only two values depending on the hybridization state of the atom type⁴⁵⁶:

The scaling factors are specified by a global parameter file: /data TemperatureFactor.dat, which subject to further optimizations.

References