Solid → liquid: derive the transition

1. Start from the same interactions, change the sampled configurations

Definitions & inputs. N is atom count; m is atomic mass; r and p are positions and momenta; U is the full interaction energy; R denotes reference lattice sites; u denotes displacement in this section; Φ is the force-constant matrix.

1. The microscopic Hamiltonian is the starting point for both phases; a liquid is not obtained by deleting its interaction energy.

   $$H=\sum_{i=1}^N\frac{\mathbf p_i^2}{2m}+U(\mathbf r_1,\ldots,\mathbf r_N)$$

2. Choose a stationary crystal structure, so its linear displacement term vanishes.

   $$\mathbf r_i=\mathbf R_i+\mathbf u_i,\qquad\left.\frac{\partial U}{\partial\mathbf u_i}\right|_{\mathbf u=0}=0$$

3. A Taylor expansion separates the static lattice energy, quadratic vibrations, and anharmonic corrections.

   $$U=U_s^0+\frac12\mathbf u^{\mathsf T}\Phi\mathbf u+O(u^3)$$

4. Keeping only the quadratic energy produces harmonic lattice dynamics, confined to a local neighborhood of the crystal.

   $$m\ddot{\mathbf u}=-\Phi\mathbf u$$

5. For thermodynamic comparison, use phase-restricted configurational integrals: Ωs samples the solid basin and its relevant permutations; Ωl samples liquid configurations with rearranging neighbors.

   $$Q_\alpha=\int_{\Omega_\alpha}e^{-\beta U(\mathbf r^N)}\,d^{3N}r,\qquad\alpha\in\{s,l\}$$

Interpretation. Melting requires comparing phase free energies. A large vibrational amplitude is a warning that the local harmonic approximation fails; it is not itself an exact melting criterion. Phase restrictions require a consistent order parameter and normalization.

• Phonopy — force constants, dynamical matrix and thermodynamics ↗
• Tong — Statistical Physics: phase transitions ↗

2. Construct comparable solid and liquid free energies

Definitions & inputs. β=1/(kBT); kB is Boltzmann’s constant; Λ is the nuclear thermal wavelength; F is Helmholtz free energy; G is Gibbs free energy; V is volume; p is pressure. λ below is a dimensionless coupling parameter, not a wavelength.

1. Integrate the Gaussian momenta and include the indistinguishability normalization consistently for the phase-restricted configuration domain.

   $$\Lambda=\frac{h}{\sqrt{2\pi m k_{\mathrm B}T}},\qquad Z_\alpha=\frac{Q_\alpha}{N!\Lambda^{3N}}$$

2. The logarithm converts the statistical weight of each phase into Helmholtz free energy.

   $$F_\alpha=-k_{\mathrm B}T\ln Z_\alpha$$

3. For a practical free-energy calculation, define a reversible path from a reference with known free energy to the target interactions, remaining in the intended phase.

   $$U_\lambda=(1-\lambda)U_{\mathrm{ref}}+\lambda U_{\mathrm{target}}$$

4. Differentiate the partition integral under the integral sign; the derivative of its exponential produces an equilibrium ensemble average.

   $$\frac{\partial F_\lambda}{\partial\lambda}=-k_{\mathrm B}T\frac{\partial_\lambda Z_\lambda}{Z_\lambda}=\left\langle U_{\mathrm{target}}-U_{\mathrm{ref}}\right\rangle_\lambda$$

5. Integrate that exact identity. Crystal references, liquid references, constraints, and finite-size corrections must be chosen and accounted for consistently.

   $$F_{\mathrm{target}}=F_{\mathrm{ref}}+\int_0^1\left\langle U_{\mathrm{target}}-U_{\mathrm{ref}}\right\rangle_\lambda\,d\lambda$$

6. In the thermodynamic limit, minimize within each phase at the imposed pressure; this allows the solid and liquid to have different equilibrium volumes.

   $$G_\alpha(T,p,N)=\min_V\{F_\alpha(T,V,N)+pV\}$$

Interpretation. There is no universal closed-form liquid free energy derivable from a harmonic solid alone. Thermodynamic integration is a computational route; it still needs a known reference, adequate sampling, and a reversible path without an uncontrolled phase change.

• Tong — Statistical Physics: phase transitions ↗
• LAMMPS — thermodynamic-integration energy derivatives ↗
• LAMMPS — crystalline free energy via an Einstein-crystal reference ↗

3. Derive the melting criterion and the latent heat

Definitions & inputs. ν is the total amount in moles; x is the liquid mole fraction; gα, hα, sα, vα are molar Gibbs energy, enthalpy, entropy, and volume. Δ always means liquid minus solid. Tm is the equilibrium melting temperature at the specified pressure.

1. Partition the conserved amount between solid and liquid and add their bulk Gibbs energies.

   $$G_{\mathrm{bulk}}=\nu[(1-x)g_s+xg_l]$$

2. Transfer material from solid to liquid while holding total amount, temperature, and pressure fixed.

   $$\left(\frac{\partial G_{\mathrm{bulk}}}{\partial x}\right)_{T,p,\nu}=\nu(g_l-g_s)=\nu\Delta g$$

3. The sign gives the direction that lowers bulk Gibbs energy. It does not remove nucleation barriers or imply instantaneous equilibration.

   $$\Delta g<0\Rightarrow\text{liquid favored},\quad\Delta g>0\Rightarrow\text{solid favored}$$

4. At coexistence, transferring an infinitesimal amount between the bulk phases costs no Gibbs energy.

   $$\Delta g(T_m,p)=0$$

5. Subtract the two phase identities at the same coexistence temperature.

   $$g=h-Ts\quad\Rightarrow\quad0=\Delta h_{\mathrm{fus}}-T_m\Delta s_{\mathrm{fus}}$$

6. The molar latent heat is the enthalpy needed per mole melted reversibly at constant pressure. It is energy, not a temperature rise.

   $$L_m\equiv\Delta h_{\mathrm{fus}}=T_m\Delta s_{\mathrm{fus}}$$

Interpretation. A solid-to-liquid derivation has two distinct tasks: calculate comparable phase free energies, then solve their equality. The liquid fraction at coexistence needs an additional constraint, such as a supplied heat budget.

• Tong — Statistical Physics: phase transitions ↗

Worked example 1 — Predict melting from a specified free-energy model

Definitions & inputs. All inputs are illustrative: Δh=10,000 J mol⁻¹ and Δs=10 J mol⁻¹ K⁻¹ at a fixed pressure. Δg is liquid minus solid molar Gibbs energy. These values are not measurements of a named material.

1. Start with the difference between the two Gibbs-energy identities at common temperature.

   $$\Delta g(T)=\Delta h-T\Delta s$$

2. Insert the supplied enthalpy and entropy differences; T is in kelvin.

   $$\Delta g(T)=10\,000-10T\quad[\mathrm{J\,mol^{-1}}]$$

3. Set the free-energy difference to zero and solve for the coexistence temperature.

   $$0=10\,000-10T_m\quad\Rightarrow\quad T_m=1000\ \mathrm K$$

4. Below Tm, the liquid has higher Gibbs energy, so the solid is favored.

   $$\Delta g(900)=10\,000-900(10)=+1000\ \mathrm{J\,mol^{-1}}$$

5. Above Tm, the liquid has lower Gibbs energy, so the liquid is favored.

   $$\Delta g(1100)=10\,000-1100(10)=-1000\ \mathrm{J\,mol^{-1}}$$

6. Check the latent-heat/entropy relation at coexistence.

   $$\Delta s_{\mathrm{fus}}=L_m/T_m=10\,000/1000=10\ \mathrm{J\,mol^{-1}K^{-1}}$$

Interpretation. Result: Tm=1000 K for this stated model. If ΔCp is significant, integrate d(Δh)/dT=ΔCp and d(Δs)/dT=ΔCp/T before solving Δg=0. This example assumes those differences are negligible.

• Tong — Statistical Physics: phase transitions ↗

Worked example 2 — Derive how pressure shifts the melting point

Definitions & inputs. Use illustrative Tm=1000 K, Lm=10,000 J mol⁻¹, vs=10×10⁻⁶ m³ mol⁻¹, vl=11×10⁻⁶ m³ mol⁻¹. Δp=10 MPa is a pressure increase. Molar quantities are used consistently.

1. Differentiate the coexistence equality while moving along the melting curve.

   $$g_l=g_s\quad\Rightarrow\quad dg_l=dg_s$$

2. Substitute the molar Gibbs differential for each phase and subtract.

   $$dg=-s\,dT+v\,dp\quad\Rightarrow\quad-\Delta s\,dT+\Delta v\,dp=0$$

3. Rearrange and use Δs=Lm/Tm; this is the inverse Clapeyron slope.

   $$\frac{dT_m}{dp}=\frac{\Delta v}{\Delta s}=\frac{T_m\Delta v}{L_m}$$

4. Calculate the expansion on melting for this illustrative material.

   $$\Delta v=(11-10)10^{-6}=10^{-6}\ \mathrm{m^3\,mol^{-1}}$$

5. Substitute the inputs. The units follow from one joule equaling one pascal cubic meter.

   $$\frac{dT_m}{dp}=\frac{1000(10^{-6})}{10\,000}=10^{-7}\ \mathrm{K\,Pa^{-1}}=0.1\ \mathrm{K\,MPa^{-1}}$$

6. Apply the local slope to the specified pressure increment.

   $$\Delta T_m\simeq(0.1)(10)=1\ \mathrm K,\qquad T_m(p+\Delta p)\simeq1001\ \mathrm K$$

Interpretation. Result: pressure raises the melting point by approximately 1 K here. If a material contracts on melting, Δv is negative and the slope reverses; that behavior must come from material data.

• Tong — Statistical Physics: phase transitions ↗

Worked example 3 — Separate warming from the latent-heat plateau

Definitions & inputs. Use 1 mol of the same illustrative material, initially solid at 900 K. Let Tm=1000 K, solid Cp=25 J mol⁻¹ K⁻¹, Lm=10,000 J mol⁻¹, and supplied heat Q=7500 J. x is the fraction melted.

1. At constant pressure, heat first increases the solid enthalpy until it reaches its melting point.

   $$Q_{\mathrm{warm}}=\nu\int_{900}^{1000}C_{p,s}\,dT$$

2. Evaluate the sensible-heat integral with constant heat capacity.

   $$Q_{\mathrm{warm}}=(1)(25)(1000-900)=2500\ \mathrm J$$

3. The remaining heat is available for changing phase at coexistence.

   $$Q_{\mathrm{melt}}=Q-Q_{\mathrm{warm}}=7500-2500=5000\ \mathrm J$$

4. Divide the available latent heat by the heat needed to melt the entire amount.

   $$Q_{\mathrm{melt}}=\nu xL_m\quad\Rightarrow\quad x=\frac{5000}{(1)(10\,000)}=0.5$$

5. Both phases remain at the coexistence temperature while the fraction changes.

   $$T=1000\ \mathrm K,\qquad\nu_l=0.5\ \mathrm{mol},\qquad\nu_s=0.5\ \mathrm{mol}$$

6. Only after the full amount melts can additional heat raise the liquid temperature; that calculation needs the liquid heat capacity.

   $$Q_{\mathrm{complete}}=2500+10\,000=12\,500\ \mathrm J$$

Interpretation. Result: 7500 J produces a half-melted sample at 1000 K, not a hotter fully liquid sample. A heating curve depends on pressure, heat capacity, latent heat, and the supplied energy.

• Tong — Statistical Physics: phase transitions ↗

4. Replace fixed neighbors with a probability distribution

Definitions & inputs. n is number density; r is separation; g(r) is the dimensionless liquid radial pair distribution; z is a mean coordination number; S(k) is the static structure factor. Here g(r) is not the molar Gibbs energy gα used above.

1. Express the pair density relative to an uncorrelated reference of the same number density.

   $$n^{(2)}(\mathbf r_1,\mathbf r_2)=n^2g(|\mathbf r_1-\mathbf r_2|)$$

2. Divide the pair density by the density of the tagged atom.

   $$\text{conditional density around a tagged atom}=ng(r)$$

3. Multiply conditional density by the volume of a thin spherical shell.

   $$dz=ng(r)\,4\pi r^2dr$$

4. Integrate over a chosen neighbor shell to obtain its mean population.

   $$z(r_a,r_b)=4\pi n\int_{r_a}^{r_b}g(r)r^2\,dr$$

5. Fourier-transform the excess pair correlation for an isotropic liquid; this connects the statistical structure to scattering.

   $$S(k)=1+4\pi n\int_0^\infty[g(r)-1]\frac{\sin kr}{kr}r^2\,dr$$

Interpretation. A crystal has lattice-dependent correlations and Bragg order. A liquid can retain a strong first-neighbor peak while losing permanent sites and long-range periodicity. Neither g(r) nor a single coordination number alone determines the melting temperature.

• Schirmacher — Theory of liquids and polymers ↗

Worked example 4 — Count neighbors in a simplified liquid shell

Definitions & inputs. Illustrative number density n=0.05 Å⁻³. Approximate g(r)=2 only in the first-shell interval ra=2 Å to rb=3 Å. One ångström (Å) equals 10⁻¹⁰ m. z is a dimensionless mean neighbor count.

1. Start with the conditional shell-count formula derived above.

   $$z=4\pi n\int_{r_a}^{r_b}g(r)r^2\,dr$$

2. Insert the density and constant shell value while using ångströms consistently.

   $$z=4\pi(0.05)(2)\int_2^3r^2\,dr$$

3. Integrate the spherical-shell radial factor explicitly.

   $$\int_2^3r^2\,dr=\left[\frac{r^3}{3}\right]_2^3=\frac{27-8}{3}=\frac{19}{3}\ \mathrm{\mathring A^3}$$

4. Multiply the factors. Å⁻³ from density cancels ų from the integral.

   $$z=4\pi(0.05)(2)\frac{19}{3}=7.958701\simeq7.96$$

5. For comparison, an uncorrelated reference with g=1 has half as many neighbors in this same interval.

   $$z_{g=1}=4\pi(0.05)\frac{19}{3}=3.979351$$

Interpretation. Result: about 7.96 neighbors on average in the chosen shell. An average need not be an integer. Local order survives even though individual neighbors may exchange over time.

• Schirmacher — Theory of liquids and polymers ↗

Mathematical tools and foundations

Derivatives and chain rule · Linear systems and elimination · Monte Carlo estimation

Use the stated mathematical conditions together with the physical assumptions. Explore the Mathematics atlas →