Lab: structure factor by hand

The setup

A density wave e^(ik·r) has wavelength 2π/k. S(k) measures how strongly the particles line up with it:

System: the 3D LJ liquid of Labs 6.9–6.11 (500 atoms, ρ = 0.8, T = 1.0, L = 8.55), 50 frames.

Step 1

S(k) of four atoms by hand

Why: The definition is just two sums of cosines and sines. On four atoms in a 1D box of length 4, we can see exactly when the waves add up and when they cancel.

What the numbers say: On the lattice every wave cancels exactly, except the one whose wavelength equals the spacing (k = 2π/1): there all four atoms are in phase and S = N, a Bragg peak that grows with the number of atoms. Shift the atoms a little and the cancellation is no longer perfect: the peak drops to 2.37 and is shared with the other k values. A liquid is the extreme case: no Bragg peaks, only a broad maximum near k ≈ 2π/(neighbour distance).

Step 2

The LJ liquid, directly from positions

Why: For a 3D box, every allowed k vector gives one value of S; vectors of the same length are averaged. The smallest k is set by the box: the longest wave that fits is L.

What the numbers say: S(k) has a sharp main peak near k = 6.6, the signature of the shell of neighbours, damped oscillations around 1 at large k, and very small values (≈ 0.05) at small k: long density waves are strongly suppressed because the liquid is hard to compress. The box decides what we can see: nothing exists below k = 0.735, and the first few k values have only 3–6 vectors each, so they are noisy. To see structure on longer length scales (clusters, nano-segregation in ionic liquids, the "pre-peak") the box must be several times larger than that length.

Step 3

S(k) from g(r): truncation and the window

Why: S(k) can also be computed from g(r), which most analyses already have. The transform needs g(r) to infinity, but we only know it up to L/2. Cutting the integral there causes errors; a window function that smoothly brings the integrand to zero (here the Lorch window, sin(πr/R)/(πr/R)) is the usual fix.

What the numbers say: At large k the transform works (10.0: 0.717 vs 0.702). At small k the raw transform goes wrong, even negative (−0.085 at k = 1.04, impossible for |Σe^(ik·r)|²), and k = 0 comes out 13 times too large (0.636), because the cut at L/2 = 4.3 leaves a step in r²(g − 1) whose transform rings. The window removes the ringing (0.050, 0.046) but smears the main peak (1.92 instead of 2.42). So: compute S(k) directly from positions when the box allows it; use the transform of g(r) for k values between the allowed ones, with a window, and know that it trades small-k accuracy against peak height.

Step 4

S(0) and the compressibility

Why: As k → 0 the structure factor measures how much the number of particles in a large region fluctuates, which is set by how compressible the liquid is: S(0) = ρk_BTκ_T = k_BT / (∂P/∂ρ). This links structure to thermodynamics and gives a check of both.

What the numbers say: The pressure route (0.0485) and the small-k structure factor (0.045–0.054) agree: a dense liquid has density fluctuations 20 times weaker than an ideal gas. This is also why the 1/N offset of g(r) in Lab 6.9 is tiny for a dense liquid (10⁻⁴) but large for a gas (1/N). For mixtures and ionic systems, the same idea gives partial structure factors: the charge–charge structure factor of an electrolyte must go to zero as k² at small k (perfect screening), a strong test of a simulation's electrostatics. Limitation: S(0) from a few small-k points needs good statistics and a large box; the pressure route needs runs at several densities.

Step 5

Comparing with experiment

Why: To compare a simulated S(k) with X-ray or neutron data, the units must match and the atoms must be weighted the way the experiment sees them.

What the numbers say: For a one-component liquid the simulated S(k) can be compared directly once k is converted (liquid argon's main peak is indeed near 2 Å⁻¹). For molecules, mixtures and ionic liquids each pair of species gets its own partial structure factor, weighted by form factors or scattering lengths, and features such as the low-k pre-peak need boxes several nanometres wide. Tools: LAMMPS compute rdf plus a transform, or compute saed/xrd; dedicated codes such as TRAVIS compute weighted S(k) directly.

Summary