International
Tables for Crystallography Volume D Physical properties of crystals Edited by A. Authier © International Union of Crystallography 2013 
International Tables for Crystallography (2013). Vol. D, ch. 2.1, pp. 286294
Section 2.1.2. Fundamentals of lattice dynamics in the harmonic approximation^{a}Institut für Physikalische Chemie, Universität Göttingen, Tammannstrasse 6, D37077 Göttingen, Germany 
In order to reduce the complexity of lattice dynamical considerations, we describe the crystal's periodicity by the smallest unit needed to generate the whole (infinite) lattice by translation, i.e. the primitive cell. Each individual primitive cell may be characterized by a running index l and a vector pointing to its origin. Let there be N atoms per cell, the equilibrium positions of which are given by being the vector of the κth atom with respect to the origin of the primitive cell (see Fig. 2.1.2.1).
The set of vectors describes the structure of the perfect lattice. At a particular time t, however, the κth atom within the lth primitive cell, denoted by (κl), may be found at a position which differs slightly from the equilibrium position, the timedependent displacement being The potential energy V of the whole crystal depends on the position vectors of all atoms, and is minimal if all atoms occupy their equilibrium positions. For small displacements, it can be expanded in a Taylor series with respect to :where denotes the Cartesian coordinate of in direction α. In the harmonic approximation, third and higherorder terms are neglected. In order to simplify the formulae, we now drop the time argument, keeping in mind that we are always dealing with dynamical displacements. The expansion coefficients in equation (2.1.2.4) are the partial derivatives of the potential energy with respect to the atomic displacements taken at the equilibrium positions:Using the matrix notationand dropping the constant , equation (2.1.2.4) readsThe product is just the force acting upon atom if the atom is displaced by (Fig. 2.1.2.2). Hence, the matrix may be regarded as a force constant matrix and its elements as force constants. These parameters may be calculated with the help of specific interaction models such as pair potentials, tensorforce models or more complicated manybody interactions.
The Hamiltonian of the perfect harmonic crystal can now be written in the formif and are the momentum and the mass of atom (κl), respectively. Consequently, the equation of motion for a particular atom (κl) is given by Solutions of this set of coupled differential equations are of the form which are plane waves with wavevector q and polarization vector . If a finite crystal is considered or if periodic boundary conditions are applied, the wavevector is restricted to a sequence of discrete and equidistant values which are, however, very close to each other. Thus, for practical work q can be treated as a continuous variable. The polarization vectors are, in general, different for every atom κ. Moreover, they depend on q and for each wavevector there are 3N different modes of vibration characterized not only by different 's but also by different vibrational frequencies ω. Hence, equation (2.1.2.10) can be written more specifically aswhere the running index labels the different fundamental vibrations or phonons. The upper index distinguishes two waves with identical frequencies which are travelling in opposite directions.
Not all of the elements of the forceconstant matrix are independent. From its definition, equation (2.1.2.5), it is clear that the forceconstant matrix is symmetric: Moreover, there are general stability conditions arising from the fact that a crystal as a whole is in mechanical equilibrium: If a macroscopic crystal is rigidly translated by a vector , no interatomic interactions are affected and, hence, the force acting on any particular atom must vanish: and, consequently, This relation is known as the condition of translational invariance.
In a similar way, it is argued that no interatomic interactions are affected when the crystal is rigidly rotated by infinitesimal amounts about arbitrary axes. This condition of rotational invariance leads to the following restrictions for the force constants: for all and .
In mechanical equilibrium, there must not be any strains within the crystal. The conditions of an unstrained crystal are also known as Huang conditions and may be formulated as for all .
All these stability conditions are independent of the particular crystal structure. There are other restrictions that are due to the symmetry of the atomic arrangement. They will be considered in detail in Section 2.1.3.
If the ansatz (2.1.2.10a) is inserted into the equation of motion (2.1.2.9), the following eigenvalue equation is obtained: The summation over all primitive cells on the righthand side of equation (2.1.2.15) yields the Fouriertransformed forceconstant matrix which is independent of l for infinite crystals. contains all interactions of type atoms with type atoms. Using this notation, equation (2.1.2.15) reduces toIf for a given vibration characterized by we combine the threedimensional polarization vectors of all atoms within a primitive cell to a 3Ndimensional polarization vector ,and simultaneously the matrices to a matrix F(q)equation (2.1.2.17) can be written in matrix notation and takes the simple formwhere the diagonal matrix contains the masses of all atoms. The matrix is called the dynamical matrix. It contains all the information about the dynamical behaviour of the crystal and can be calculated on the basis of specific models for interatomic interactions. In analogy to the matrices , we introduce the submatrices of the dynamical matrix: Owing to the symmetry of the forceconstant matrix, the dynamical matrix is Hermitian:^{1}or more specifically Obviously, the squares of the vibrational frequency and the polarization vectors are eigenvalues and corresponding eigenvectors of the dynamical matrix. As a direct consequence of equation (2.1.2.20), the eigenvalues are real quantities and the following relations hold:Moreover, the eigenvectors are mutually orthogonal and can be chosen to be normalized.
The wavevector dependence of the vibrational frequencies is called phonon dispersion. For each wavevector q there are 3N fundamental frequencies yielding 3N phonon branches when is plotted versus q. In most cases, the phonon dispersion is displayed for wavevectors along highsymmetry directions. These dispersion curves are, however, only special projections of the dispersion hypersurface in the fourdimensional q–ω space. As a simple example, the phonon dispersion of b.c.c. hafnium is displayed in Fig. 2.1.2.3. The wavevectors are restricted to the first Brillouin zone (see Section 2.1.3.1) and the phonon dispersion for different directions of the wavevector are combined in one single diagram making use of the fact that different highsymmetry directions meet at the Brillouinzone boundary. Note that in Fig. 2.1.2.3, the moduli of the wavevectors are scaled by the Brillouinzone boundary values and represented by the reduced coordinates ξ. Owing to the simple b.c.c. structure of hafnium with one atom per primitive cell, there are only three phonon branches. Moreover, for all wavevectors along the directions [00ξ] and [ξξξ], two exhibit the same frequencies – they are said to be degenerate. Hence in the corresponding parts of Fig. 2.1.2.3 only two branches can be distinguished.
Whereas in this simple example the different branches can be separated quite easily, this is no longer true for more complicated crystal structures. For illustration, the phonon dispersion of the highT_{c} superconductor Nd_{2}CuO_{4} is shown in Fig. 2.1.2.4 for the main symmetry directions of the tetragonal structure (space group , seven atoms per primitive cell). Note that in many publications on lattice dynamics the frequency is used rather than the angular frequency ω.
The 21 phonon branches of Nd_{2}CuO_{4} with their more complicated dispersion reflect the details of the interatomic interactions between all atoms of the structure. The phonon frequencies ν cover a range from 0 to 18 THz. In crystals with strongly bonded molecular groups, like SiO_{4} tetrahedra in quartz or SO_{4} tetrahedra in sulfates, for example, the highest frequencies are found near 35 THz and correspond to bondstretching vibrations. Soft materials like organic molecular crystals, on the other hand, exhibit a large number of phonon branches within a rather small frequency range which cannot easily be separated. Deuterated naphthalene (C_{10}D_{8}) is a well investigated example. The lowfrequency part of its phonon dispersion is shown in Fig. 2.1.2.5.
Whereas neutron inelastic scattering is the most powerful method for the determination of phonons at arbitrary wavevectors, long wavelength phonons may also be detected by optical spectroscopy. The determination of phonon frequencies alone is, however, not sufficient for a concise determination of dispersion branches. Rather, individual phonons have to be assigned uniquely to one of the 3N branches, and this may prove a rather hard task if N is large. Here, symmetry considerations of eigenvectors are of special importance since phonons belonging to the same branch must exhibit the same symmetry properties. Moreover, inspection of Figs. 2.1.2.3 to 2.1.2.5 shows that some of the branches cross each other and others do not. It is a general statement that crossing is only allowed for branches with different symmetries – a property which yields a classification scheme for the different phonon branches. The symmetry of fundamental vibrations of a lattice will be discussed in some detail in Section 2.1.3.
In the limit of long wavelengths, there are always three particular modes with identical polarization vectors for all atoms, which will be considered in the following. At exactly (the Γ point) or infinite wavelength, the eigenvalue equation (2.1.2.15) reduces to One immediately recognizes that there are special solutions withi.e. the (massweighted) eigenvectors of all atoms are identical. There are three orthogonal eigenvectors of this kind and the displacement pattern of such phonons corresponds to rigid translations of the whole lattice along the three orthogonal coordinates in direct space. These motions do not affect any interatomic interaction. Hence, there is no change in potential energy and the condition of translational invariance (cf. Section 2.1.2.2) guarantees that the frequencies of these modes are zero:The phonon branches that lead to zero frequency at the Γ point () are called acoustic, whereas all other branches are called optic. The dispersion of acoustic branches in the vicinity of the Γ point can be investigated by expanding the phase factor in equation (2.1.2.15) in powers of q. Using (2.1.2.28) one obtainsNeglecting higherorder terms, summing up both sides of equation (2.1.2.30) over κ and multiplying by yieldsM being the total mass of all atoms within the primitive cell (). The first term on the righthand side is zero according to equation (2.1.2.29). The second term vanishes due to the symmetry property of the forceconstant matrices, equation (2.1.2.23). Hence (2.1.2.31) is simplified toThe righthand side no longer depends on the moduli of the wavevector and displacement but only on their orientation with respect to the crystal lattice. Consequently, acoustic dispersion curves always leave the Γ point as a straight line with a constant slope ().
The displacement pattern of these longwavelength modes corresponds to a continuous deformation of a rigid body. Hence, acoustic phonons near the Γ point can be regarded as sound waves and the slope of the dispersion curve is given by the corresponding sound velocity, Sound velocities, on the other hand, can be calculated from macroscopic elastic constants using the theory of macroscopic elasticity (cf. Chapter 1.3 ). Thus we are able to correlate macroscopic and microscopic dynamic properties of crystals. Using the generalized Hooke's law, the equation of motion for the dynamic deformation of a macroscopic body may be written as being the macroscopic density, the ith Cartesian component of the deformation and () the symmetric tensor of elastic stiffnesses, which is discussed in detail in Chapter 1.3 . The solution of this differential equation using plane waves,leads to the following relation: If we define the components of the propagation tensor byequation (2.1.2.36) may be written as the eigenvector equationFor a given propagation direction as defined by the Cartesian components of q, the eigenvectors of the corresponding propagation tensor yield the polarization of three mutually orthogonal deformation waves. Its eigenvalues are related to the respective sound velocities . If the tensor of elastic stiffnesses is known, the elements of and, hence, the velocity of elastic (sound) waves can be calculated for arbitrary propagation directions (see Section 1.3.4 ). These data, in turn, allow the prediction of the slopes of acoustic phonon dispersion curves near .
The planewave solutions (2.1.2.10) of the equations of motion form a complete set of orthogonal functions if q is restricted to the first Brillouin zone. Hence, the actual displacement of an atom (κl) can be represented by a linear combination of the : Since this displacement is an observable quantity, it must correspond to a real vector, not a complex one. Hence, the coefficients obey the relationand equation (2.1.2.39) reduces towhere If the displacement vectors are combined to form a 3Ndimensional vector in analogy to the formation of the eigenvector from the individual polarization vectors [equation (2.1.2.18)] we obtainThus, the atomic displacement is a linear combination of the eigenvectors of the dynamical matrix. The coefficients are called normal coordinates. They reflect the relative weight and amplitude of a particular vibrational mode () which is temperaturedependent and may be determined by statistical methods.
In terms of these normal coordinates, the Hamiltonian of the lattice (2.1.2.8) is reduced to a sum of independent harmonic oscillators. These are called phonons and may be regarded as quantum quasiparticles.( is the number of primitive cells within the crystal.)
Lattice vibrations that are characterized by both the frequencies and the normal coordinates are elementary excitations of the harmonic lattice. As long as anharmonic effects are neglected, there are no interactions between the individual phonons. The respective amplitudes depend on the excitation level and can be determined by quantum statistical methods. The energy levels of a lattice vibration () are those of a single harmonic oscillator: as illustrated in Fig. 2.1.2.6. The levels are equidistant and the respective occupation probabilities are given by Boltzmann statistics: In the quasiparticle description, this quantity is just the probability that at a temperature T there are n excited phonons of frequency . Moreover, in thermal equilibrium the average number of phonons is given by the Bose factor:and the corresponding contribution of these phonons to the lattice energy is The meansquare amplitude of the normal oscillator coordinate is obtained asAt high temperatures (), the phonon number, the corresponding energy and the amplitude approach the classical values of respectively. Note that occupation number, energy and amplitude merely depend on the frequency of the particular lattice vibration. The form of the corresponding eigenvector is irrelevant.
The total energy stored in the harmonic phonon system is given by the sum over all phonon states (): Related thermodynamic quantities like the internal energy or the heat capacity are determined by the frequency distribution of the lattice vibrations rather than by details of the phonon dispersion. Hence, it is useful to introduce the phonon density of states in such a way that is the number of phonons with frequencies between and . Using this quantity, the sum in (2.1.2.52) may be replaced by an integral expression:Here, is the energy at . The derivative with respect to temperature provides the lattice heat capacity at constant volume: As an example, Fig. 2.1.2.7 displays the phonon dispersion of GaAs as determined by inelastic neutron scattering along with the phonon density of states. Obviously, even in this relatively simple substance (DOS) exhibits a rather complicated multipeak structure. Integral properties like the heat capacity are, however, not very sensitive to details of . There are two well known approximations for that are able to reproduce some prominent features of :
In the Einstein model, it is assumed that all phonons exhibit the same frequency (Einstein oscillator) and is represented by a delta function: Consequently, the heat capacity turns out to bewhere we use the abbreviation which is the Einstein temperature.
At low temperatures, this model predicts an exponential temperature dependence of the heat capacity , which does not correspond to the experimental findings in most substances. Here, the Debye model provides a significant improvement.
In contrast to the Einstein model, which takes only one optic mode into account, the Debye model is restricted to acoustic modes that exhibit a linear dispersion close to the Γ point (see Section 2.1.2.4). Neglecting any deviation from linear behaviour, we get the simple result that the density of states is proportional to the square of the phonon frequency. The total number of phonon states is, however, given by , which is the number of all dynamical degrees of freedom of the whole system. Consequently, the frequency spectrum is assumed to be limited to frequencies below a particular value according to This limiting frequency is called the Debye frequency and is related to an appropriate average of (longitudinal and transverse) sound velocities and exhibits large values for hard materials. Fig. 2.1.2.8 compares schematically the true phonon density of states with the Debye approximation. The density of phonon states may thus be represented by and, correspondingly, the heat capacity is yielding a temperature dependence as shown in Fig. 2.1.2.9.

Schematic representation of the true phonon density of states (solid line) along with the Debye approximation (dotted line). Note that the areas under the two curves are identical. 
is the Debye temperature, which is defined as At low temperatures, the heat capacity is proportional to , in excellent agreement with most experiments:
It is not surprising that the Debye model provides a reasonable description of the lowtemperature heat capacity, since in this temperature regime well below the Debye temperature, optical phonons are hardly excited and the heat capacity is dominated by the lowfrequency acoustic modes which are modelled exactly. At higher temperatures it is, however, necessary to take into account the thermal excitation of (in general less dispersive) optic modes. This can be achieved either by introducing a temperature dependence of the Debye temperature or by mixing a Debye term like (2.1.2.60) and Einstein terms like (2.1.2.56).
As an example, we consider the case of GaAs, the density of states of which is shown in Fig. 2.1.2.7. Obviously there are two very pronounced peaks at high frequencies, which are due to nearly dispersionless optical phonon branches. These modes may therefore be regarded as Einstein oscillators. The remaining acoustic branches lead to the more continuous part of the spectrum at lower frequencies, which may be approximated by a Debye law.
So far, we have always assumed that the crystal volume is constant. As long as we are dealing with harmonic solids, the thermal excitation of phonons does not result in a mean displacement of any atom. Consequently, thermal expansion cannot be understood in the harmonic approximation. It is due to the fact that there are anharmonic contributions to the lattice energy, i.e. third and higherorder terms in the expansion with respect to atomic displacements [equation (2.1.2.4)]. Moreover, in an anharmonic lattice phonons are no longer independent elementary excitations. Rather, different lattice vibrations interact with each other leading to temperaturedependent frequency shifts, damping etc. Quantitatively, anharmonic effects may be analysed by means of perturbation theory, which is, however, beyond the scope of the present article. Details may be found, for example, in the monograph The Physics of Phonons (Reissland, 1973).
Some aspects of anharmonicity can, however, be discussed on the basis of the quasiharmonic model. This approach makes use of the fact that the atomic interactions vary with the interatomic spacing and, hence, with the volume or, more generally, with any type of lattice deformation. The phonon frequencies will therefore depend on the deformation as well. Using the deformed lattice as a new reference frame for lattice dynamical calculations, the corresponding frequencies may be obtained again on the basis of a harmonic model with modified effective force constants. The comparison of phonons of both the original and the (arbitrarily) deformed lattice finally yields the partial derivatives of the frequencies with respect to the components of the strain tensor.
If the deformation is exclusively due to a change of temperature, the phonon frequency shifts can thus be related to the coefficients of thermal expansion. In this approximation, any intrinsic temperature dependence of phonon frequencies due to phonon interactions is neglected. Note that just those effects are, however, of particular importance if displacive phase transitions that are associated with soft phonon modes are considered.
In the quasiharmonic approximation, we use the thermodynamic relation between the Helmholtz free energy A and the partition function Z: Z is given in terms of the energy levels of the independent harmonic oscillators: where Φ is the potential energy of the crystal if all atoms occupy their equilibrium positions. Hence, the following expression for A results: Elementary thermodynamics yields the pressure p as the partial derivative of A with respect to the volume at constant temperature: This relation may be generalized if not only volume changes are taken into account but also arbitrary deformations as described by the strain tensor : with the stress tensor . Note that the hydrostatic pressure p and the relative volume change are given by the traces of and , respectively: Taking the temperature derivative of (2.1.2.66), we obtain Using Euler's relations, the lefthand side may be written as where the tensor of the elastic stiffnesses and the tensor of thermal expansion have been introduced.
The free energy depends on the lattice deformations via the phonon frequencies. Hence, using (2.1.2.65), the righthand side of (2.1.2.69) is evaluated as under the assumption that the phonon frequencies do not depend explicitly on the temperature. [ is the Bose factor, (2.1.2.47)].
Let us denote the contribution of a single phonon () to the heat capacity at constant volume by Then the combination of (2.1.2.69), (2.1.2.70) and (2.1.2.71) yields the result with the generalizedmode Grüneisen parametersThe set of equations (2.1.2.75) (for ) provide relations between the variation of phonon frequencies with the lattice deformations on the one hand and the tensors of elastic stiffnesses and thermal expansion on the other hand.
For cubic crystals, the tensor of the thermal expansion is diagonal,where β represents the coefficient of volume expansion. If Voigt's notation is used for the elastic stiffnesses (cf. Section 1.3.3.2.2 ), equation (2.1.2.75) reduces toIf there is an isotropic deformation, the shift of phonon frequencies may be described by an averagedmode Grüneisen parameter: and (2.1.2.78) may be rewritten asIntroducing the mean Grüneisen parameter γ by the summation over all phonon states,we arrive at Remembering that in cubic crystals the expression represents the isothermal compressibility, we find the commonly used scalar form of equation (2.1.2.75): which relates the thermodynamic quantities thermal expansion, compressibility and heat capacity with the mean Grüneisen parameter. For most substances, γ exhibits values between 1 and 4 which are hardly temperature dependent. Hence, equation (2.1.2.84) may be regarded as an equation of state for solid systems.
Experimentally, it is almost impossible to determine the heat capacity at constant volume since the thermal expansion cannot be easily compensated. The more convenient quantity is , the heat capacity at constant pressure. There is a simple thermodynamic relation between the two quantities, and hence the following equation is obtained:
References
Reissland, J. A. (1973). The physics of phonons. London: Wiley. (ISBN 0–471–71585–9.)