Debye vs. Einstein Solids

Previously, both in comp lab and in lecture, we studied the Einstein solid, a model of a solid that treats each constituent atom as an independent harmonic oscillator. We found that the heat capacity per particle in dimensionless units was given by \[\frac{C_V}{Nk}=\frac{3(\eps/kT)^2\e^{\eps/kT}}{(1-\e^{\eps/kT})^2}\,,\] where \(\eps = hf\) is the energy of the first excited state of the oscillator. Although the Einstein model gives qualitatively similar predictions to what is observed experimentally, the exponential drop in heat capacity for low temperatures does not agree quantitatively. The main problem is that each constituent particle in a solid does not oscillate independently from its neighbors.

The Debye model of a solid takes these collective oscillations into account, and more accurately reproduces the low-temperature behavior of solids with respect to their heat capacities. The main idea is to treat the occupation numbers in the oscillator modes with a Planck distribution, but with a high-frequency cutoff which is controlled by the lattice spacing between the particles. Although light waves can have arbitrarily short wavelengths, the shortest wavelength of sound waves in a solid is one half of the lattice spacing.

In the Debye approximation, we define a dimensionless variable \(x\) which gives an upper bound on which frequecies of modes we will sum over, and we find that this bound corresponds to the Debye temperature \(T_D\), defined implicitly by \[x_\max = \frac{T_D}{T} = \frac{1}{T} \frac{h c_s}{2k} \biggl(\frac{6 N}{\pi V}\biggr)^{1/3} \,.\] Here \(c_s\) is the speed of sound in the solid, analogous to the speed of light with our treatment of the Planck distribution. The factor of \(\left(\tfrac{6N}{\pi}\right)^{1/3}\) comes from Debye's approximation of a square volume in mode space by a sphere with an equal number of modes. This leads to a formula for the total average energy given by \[U=\frac{9NkT^4}{T_D^3}\int_0^{T_D/T} \frac{x^3}{\e^x-1}\dd x\,,\] and the heat capacity given by \[C_V= 9Nk\biggl(\frac{T}{T_D}\biggr)^3\int_0^{T_D/T}\frac{x^4\e^x}{(\e^x-1)^2}\dd x\,.\] If you carefully check these formulas, you'll notice that the factor of 2 we used for counting polarization modes of photons has been replaced with a factor of 3 since sound waves can also be longitudinally polarized.

To accurately compare the Einstein and Debye solids, we have to make sure that we're summing over the same number of modes. Remember, Debye used a spherical integral to approximate the mode sum, whereas the Einstein model simply counts the number of modes in a cube. Therefore, we need a conversion factor \((\pi/6)^{1/3}\) to be able to compare the models. Let's plot the dimensionless heat capacities per particle predicted by both models for temperatures below the Debye temperature.

c = (pi/6)^(1/3);
CE = @(y) (c./y).^2.*exp(c./y)./(exp(c./y)-1).^2;
CE = @(y) (c./y).^2.*exp(c./y)./(exp(c./y)-1).^2;
% This vectorizes the integral. 
CD = @(y) 3*y.^3.*integral( @(z) (1./y).*(z./y).^4.*exp(z./y)./(exp(z./y)-1).^2, ... 
        0, 1,'ArrayValued',1);

hold on;
ezplot(CD,[0,1]);
ezplot(CE,[0,1]);
hold off;
title('Comparison of Einstein and Debye heat capacities');
xlabel('$T/T_D$','interpreter','latex');
ylabel('$C_V/3Nk$','interpreter','latex');
legend('Debye','Einstein','Location','southeast');

You might think, "What's the big deal? The Einstein model is close enough!" This is an artifact of the way we are plotting the data. In fact, the prediciton of the Einstein model is radically different from that of the Debye model. There are two ways to see this: looking at the percent error between the two models, and using a log plot to see the scaling behavior for low temperature.

Unfortunately, I don't know how to avoid this for loop in the following code.

r=.01:.001:1;
debye=CD(r); 
einstein=CE(r);
perror = abs(einstein-debye)./debye;
plot(r,perror)
title('Fractional error between Einstein and Debye heat capacities');
xlabel('$T/T_D$','interpreter','latex'); 
ylabel('$|C_E - C_D|/C_D$','interpreter','latex');
loglog(r,debye,r,einstein);
title('Log-Log plot of Einstein and Debye heat capacities');
xlabel('$T/T_D$','interpreter','latex'); 
ylabel('$C_V / 3NkT$','interpreter','latex');
legend('Debye','Einstein','Location','southeast');

As these two plots illustrate, below about half of the Debye temperature, the Einstein model really makes quite strikingly bad predictions (though it does still get the qualitative story correct).

Spin Waves and Ferromagnetism

A ferromagnet is a material that at small enough temperatures will spontaneously magnetize even in the absense of an external magnetic field. Iron is the classic example of a ferromagnet. Microscopically, the reason that some materials are ferromagnetic is that they have interactions between their constituent spins which energetically favor them to align along the same direction as their neighbors.

At zero temperature, the ferromagnet is perfectly polarized along a particular direction in space. (That direction is chosen spontaneously at random by quantum and thermal fluctuations as the magnet is cooled below a certain critical temperature.) The total magnetization for a system with \(N\) magnetic dipoles is then given by \(M = 2\mu_B N\), where \(\mu_B\) is the Bohr magneton.

When \(T \gt 0\), there are excitations that take the form of spin waves. These are periodic precessions of the spins with a characteristic wavelength, as shown in the figure below taken from the text. Just like sound waves, these spin waves can only take energies in discrete quanta of energy, where the exact value is determined by the microscopic details of the ferromagnet such as the coupling strength and the lattice spacing. We can think of these discrete packets of energy as being particles (much like we did for phonons). These particles are called magnons.

One magnon changes the spin of the ferromagnet by \(h/2\pi\), hence the total magnetization decreases by about \(2 \mu_B\). In addition, each magnon carries with it an energy, and it has the same dispersion relation as a typical massive nonrelativistic particle, namely \[\eps = \frac{p^2}{2 m^*} \,,\] where \(m^*\) is an effective mass that again depends on the microscopic details of the material, and \[ p = \frac{hn}{2L} \] for a magnon with mode number \(n = (n_x,n_y,n_z)\) and \(L\) the size of the ferromagnet. This is in contrast to phonons, which behave in analogy with photons and have a linear energy-momentum dispersion relation. Magnons may also only be polarized in one possible direction (compared to photons which carry two polarizations or phonons which carry three).

The total number of magnons at a (low) temperature \(T\) should be given by the Planck distribution summed over all the magnon modes. We can do the same trick that we did for blackbody radiation and the Debye model and convert that sum over modes into an integral, and we get \[N_m = \frac{\pi}{2} \int_0^\infty \frac{n^2}{\e^{\eps/kT}-1}\dd n \,,\] where the upper limit is set to infinity because this picture only really holds at low temperatures. By changing variables to \(x = \eps/kT\), we can write this as \[N_m = 2\pi V \biggl(\frac{2m^*kT}{h^2}\biggr)^{3/2} \int_0^\infty\frac{\sqrt{x}}{\e^x-1} \dd x \,.\]

a = integral(@(x) x.^(1/2)./(exp(x)-1),0,Inf)
a =

    2.3152

If each magnon reduces the average magnetism by \(2\mu_B\), then the fractional decrease in magnetism at temperature \(T\) is simply \[ \frac{2 \mu_B N_m}{2 \mu_B N} = 2 \pi a \frac{V}{N} \biggl(\frac{2m^*kT}{h^2}\biggr)^{3/2} = \biggl(\frac{T}{T_0}\biggr)^{3/2}\,,\] where \[T_0 = \frac{h^2}{2 m^* k} \biggl(\frac{N}{2 \pi a V}\biggr)^{2/3}\,.\]

Magnons carry energy, hence they will contribute to the heat capacity of a material. If we follow the same mode-sum-to-integral method as before, we can compute the average total energy first to be \[U = 2\pi V k T \biggl(\frac{2m^*kT}{h^2}\biggr)^{3/2} \int_0^\infty\frac{x^{3/2}}{\e^x-1} \dd x \,.\]