Wiki

:::bash
$DIR/fit.lj/fit.LJ.water.constr.py -ene E.comb.dat -ljf dimer.ljf \
  [-out FILE] [-prm FILE] [-bs NUM FILE] [-v] [-unfweight] [-max E_max]

Input arguments

:::text
-ene E.comb.dat: ab-initio and force-field energies
-ljf dimer.ljf: geometrical information of each snapshot
-out FILE: output optimized LJ coefficients to FILE
-prm FILE: Don't fit, instead measure quality of parameters 
           provided by FILE (same formatting as for [-out FILE])
-bs NUM FILE: bootstrap data NUM times and output to FILE
-v: verbose
-unfweight: Uniform weight. Don't refit to eliminate false-positives (see below).
-max E_max: don't include ab-initio energies that are larger than E_max (in kcal/mol).

The E.comb.dat and dimer.ljf are determined from the $DIR/fit.lj/run.fit.sh script (see Extract raw energies).

The fit runs twice to try to eliminate false-positives:

1st round: Eqm-Eel, where Eqm adn Eel are the ab initio and electrostatic energies of a snapshot is used as a weight. Snapshots for which this difference is larger than thermal energy at room temperature are exponentially suppressed (akin to a Boltzmann factor).
2nd round: min(Eqm-Eel,Elj), where Elj is the LJ energy, is used as the weight. If Elj is (too) low, it will have a strong weight to possibly correct for it.

The script will fit all LJ coefficients (except TIP3P's) using the Lorentz-Berthelot (LB) mixing rules.

Output example

:::text
Optimization with LB mixing rule and constraints
  type  epsilon  Rmin/2
    CR  -0.0100  1.5698
     F  -0.0100  1.0479
  HCMM  -0.0181  1.4282
   HOR  -0.0100  1.0000
    OR  -0.1894  1.9032
----------------------------
       R^2:  0.23
      RMSE:  0.8877 kcal/mol
Boltz-RMSE:  1.5848 kcal/mol
       MAE:  0.4048 kcal/mol
 Boltz-MAE:  0.7456 kcal/mol
----------------------------

where the optimized LJ coefficients are first displayed (epsilon and Rmin/2), and the quality of the fit is shown afterwards. R^2 is the cross-correlation coefficient, RMSE is the root-mean-squared error, Boltz-RMSE is the Boltzmann-weighted RMSE (according to the ab-initio energies, using room temperature), MAE is the mean-absolute error, and Boltz-MAE is the Boltzmann-weighted MAE.

Note that this fitting procedure can seriously suffer from undersampling of buried atoms (e.g., methyl carbon). Akin to any force-field fitting methodology, the solution is underdetermined (i.e., no unique solution), meaning that a large number of alternative sets tend to yield very similar quality of fit. Also, the resulting methodology does not necessarily provide accurate thermodynamic properties for liquids.