Lab: dielectric constant and rotation by hand

The setup

Every polar molecule carries a dipole μᵢ. Their sum is the box dipole M = Σ μᵢ. With Ewald sums and conducting ("tin-foil") boundaries, the default in most codes:

To keep the physics visible we use a model: 216 unit dipoles on a 6 × 6 × 6 lattice (periodic), each rotating by small random steps (rotational diffusion D_r = 1), with a coupling J between nearest neighbours (k_BT = 1). Real dipolar liquids also need the long-range dipole–dipole interaction, which Ewald handles.

Step 1

The box dipole and ε by hand

Why: ε is built from M, so first compute M for a few molecules, then put real numbers for water into the fluctuation formula. Comparing with independent dipoles shows how much the molecules' cooperation matters.

What the numbers say: Three dipoles pointing in different directions partly cancel (|M|² = 7.55 against Σ|μ|² = 14.81). If water molecules pointed independently, ε would be only about 20; the real value near 70–80 needs G_K ≈ 3.7: hydrogen bonds make neighbouring dipoles point the same way, so their fluctuations add up. This is why ε is a sensitive test of a water model: a small change in how molecules orient together changes ε a lot. Limitation: the formula assumes conducting boundaries; with a reaction field or vacuum boundary it changes form, and for ionic systems the translational current also contributes to M and has to be separated.

Step 2

The Kirkwood factor in the lattice model

Why: With the model we can switch the neighbour coupling on and off and watch G_K respond. For J = 0 the answer must be exactly 1.

What the numbers say: Independent dipoles give G_K = 0.97 ± 0.04, consistent with 1. A modest preference for parallel neighbours (J = 0.4 k_BT) triples G_K and would raise ε from 19 to 59; antiparallel neighbours halve it. ε measures correlations between molecules, not the dipole of one molecule. That is also why ε is hard to fit: small errors in the orientational correlations are multiplied by 18.65 (for water) in ε.

Step 3

Rotational correlation functions: P1, P2, diffusion or jumps

Why: How does one molecule change direction? C₁ (from the first Legendre polynomial) is what dielectric spectroscopy sees; C₂ is what NMR, Raman and fluorescence anisotropy see. For small random rotations theory says τ₁/τ₂ = 3; for big jumps to new directions, τ₁/τ₂ = 1. So the ratio tells the mechanism.

What the numbers say: P₂ falls faster than cos θ for small steps, because squaring makes it sensitive to small tilts (0.933 → 0.806 already at t = 1). Both mechanisms give the same τ₁ = 0.5, so C₁ alone cannot tell them apart; the ratio τ₁/τ₂ can: 3.01 for diffusion, 1.03 for jumps, exactly as theory predicts. In water the ratio is below 3, a signature of the large-angle jumps by which molecules swap hydrogen-bond partners. Limitation: real molecules are not spheres; a non-spherical molecule has different rates about different axes, and a single τ hides that.

Step 4

One molecule vs the whole box

Why: Dielectric spectroscopy measures how fast the box dipole M relaxes, not one molecule. Are they the same? Compare the single-dipole τ₁ with τ_M, the area of ⟨M(0)·M(t)⟩/⟨M²⟩.

What the numbers say: For independent dipoles the box forgets at the same rate as one molecule (0.44 vs 0.49, equal within noise). With parallel coupling the box dipole lives 2.7 times longer than any single molecule's orientation: a fluctuation of M involves many aligned molecules and decays slowly. With antiparallel coupling it decays faster. τ_M / τ₁ follows G_K (2.68 vs 3.10, 0.46 vs 0.51). So a dielectric relaxation time is not a molecular rotation time; comparing the two needs a theory that includes these correlations. For ionic liquids the translational current adds its own contribution to the measured spectrum.

Step 5

Convergence: why ε needs very long runs

Why: G_K is built from M, one vector for the whole box, relaxing on the slow time τ_M. Like the conductivity (Lab 6.7) and the viscosity (Lab 6.10) it gets no help from the number of molecules. How precise is it after a given run length?

What the numbers say: Even 65 relaxation times of M give G_K only to ±12 %, which would be ±7 in a water-like ε. Getting ε to a few percent takes hundreds of τ_M; for water that means tens of nanoseconds, and for viscous ionic liquids, where M relaxes in nanoseconds, microseconds or many independent runs. A quoted ε without a run length and an error bar is hard to interpret. Bigger boxes do not help, because M is one number per frame whatever N is.

Step 6

Practical notes

Why: A few choices in real codes change the formula or the number, and are easy to overlook.

What the numbers say: With the same ⟨δM²⟩ the three boundary conditions give different ε, so the formula must match the electrostatics of the run. A nonzero ⟨M⟩ usually signals an error (molecules broken across the periodic boundary, or a slowly drifting system) and must be subtracted. GROMACS: gmx dipoles (which also prints G_K); LAMMPS: compute dipole with your own post-processing.

Summary