What is a force field?

The idea

Every lesson so far needed a potential energy U(r₁, …, r_N): the forces are its gradient. For real molecules we cannot solve quantum mechanics at every step, so we use a force field: an approximate, classical U built from simple pieces.

• Bonded terms (Lessons 4.1–4.2): bonds, angles, proper and improper dihedrals. They act only between atoms connected through the molecule's topology.
• Non-bonded terms: Lennard-Jones (repulsion plus dispersion) and Coulomb interactions between fixed partial charges. They act between all pairs in different molecules, and between atoms in the same molecule that are more than three bonds apart (1–4 pairs scaled, 1–2 and 1–3 excluded).

Atom types. Parameters are not per element but per atom type: a carbon in CH₃, in C=O and in a benzene ring each have their own LJ parameters and bonded terms. A force field is essentially a large table of these types and the parameters for every combination.

Where the numbers come from. Bond and angle constants from vibrational spectra or quantum (QM) calculations. Torsions fitted to QM rotation profiles. Partial charges fitted to the QM electrostatic potential (e.g. RESP) or chosen for liquid properties. LJ parameters tuned so that the simulation reproduces densities, heats of vaporisation and solvation free energies. Every parameter is fitted together with all the others, the functional form, the combining rule and the water model.

Families. AMBER (GAFF for small molecules), CHARMM (CGenFF), OPLS-AA and GROMOS each have their own rules (combining rule, 1–4 scaling, cutoffs). Never mix parameters from different families. Beyond fixed charges, polarizable force fields (Drude, AMOEBA) let charges respond to their environment, and machine-learned potentials (MACE, NequIP, …) learn U directly from QM data, at a higher cost.

The algorithm

Evaluating a classical force field each step:

U = 0
for each bond (i, j):         U += ½ k_b (r_ij - r0)²
for each angle (i, j, k):     U += ½ k_θ (θ_ijk - θ0)²
for each dihedral (i,j,k,l):  U += Σ_n k_n [1 + cos(n φ - δ_n)]
for each improper:            U += ½ k_ξ (ξ - ξ0)²
for each non-excluded pair (i, j) within the cutoff:
    σ_ij, ε_ij from the combining rule
    s = 1–4 scaling factor (or 1)
    U += s_LJ * 4 ε_ij [(σ_ij/r)¹² - (σ_ij/r)⁶]
       + s_QQ * q_i q_j / (4π ε0 r)              # plus long-range electrostatics (Ch. 5)

1. Combining rules: Lorentz–Berthelot σ_ij = (σ_i + σ_j)/2, ε_ij = √(ε_iε_j) (AMBER, CHARMM), or geometric for both (OPLS). Use the one your force field was fitted with.
2. The topology lists atoms, types, charges and all bonded terms. Exclusions are generated from it automatically (GROMACS nrexcl = 3).
3. Small-molecule parameters come from tools such as antechamber/GAFF, CGenFF, LigParGen (OPLS) or the ATB (GROMOS). Always check the assigned charges and torsions.
4. Validate before use: density, heat of vaporisation, diffusion or a solvation free energy against experiment.