Compared to Maxwell-Boltzmann particles, Fermi-Dirac (FD) particles behave in very mysterious ways. This is due to their half integer spin which permits only two particles, spin up and spin down, in any single (spatial) wave state. Here we will limit our animation to 2 dimensions because three dimensions inovlve too many particles and its states are hard to display. For an idea of how states in 3 dimensions would look see Simple Fermi Surface. The physics of solids are mainly due to the electrons in solids which have half integer spin which make Fermi-Dirac particles as mentioned above. So I have decided to start with the simplest case where the Fermi-Dirac particles are confined to a container and therefore have simple quantum states. I've found even this to be a challenging task to graphically and physically show the changes of state occupation due to effective temperarture above absolute zero but I think I've succeeded. Here I do the numerical interaction calculations that result in the momentum state changes due to the effective temperature. The resulting state changes then are displayed as a density of momentum states. If we know this density, then we can compute macroscopic quantities like the pressure, heat capacity, and thermal and electrical conductivity. This dispenses with need for the term called the entropy which was conceived when all materials were thought of as fluids before the particle and quantum natures of materials were realized.
The first requirement is that the wave function be zero at the `x=0` and `x=L` boundaries so a 1D solution must be `psi(x)=sin kx` for any allowed `k`. Therefore `k` would be `k=pin/L=ndeltak` where n is an integer. The second requirement is that all wave functions be orthogonal on the interval `x=0->L`. Orthogonal means that the interference product oftwo states integrated over `x=0->L` must be zero. That results in the following expression
`int_0^Lsin npix/L sin mpix/L dx=0\ \ m ne n\ \ \ \(8)`
where m and n are integers. Expanding the integrand into the cosines of the sums and difference of the arguments we have:`sin npix/L Lsin mpix/L=1/2 {cos[(m-n)pix/L]-cos[(m+n)pix/L]}\ \ \ \ \(9)`
which obviously integrates to zero on the interval `x=0->L`. Totally analogous results apply for the `y` and `z` coordinates.Suppose the container is square with side `L`. The basic unit of wave vector is then
`delta k=pi/L`
Then the `vec k` states of lowest magnitude form in rings of radius `k_n=ndeltak` in `k` space. If there are `m` states in the `nth` `k_r` ring then the state locations in `k` space can be written
`vec k_("ring")=(ndeltak)sum_(i=1)^m [cos(2pii/m)hatx+sin(2pii/m)haty]`
The states in `k_x` and `k_y` space are shown in Canvas 1. The quadrilaterals are of the area of the state (`deltak^2`) The red and black dots are symbols for spin up and spin down FD particles. In this display, states in occupied rings each contain a pair of FD particles so the occupation of the states shown is 100%. However, if there are not enough FD particles to fill a given ring, then its occupation may not be 100% and some states may contain just one FD particle. Generally the number of spin up and spin down FD particles in a given container will be equal. Otherwise the container will have a net magnetic moment. Simple solids do not usually have a magnetic moment.Starting here in this document we will assume that the FD particle mass, m, is unity and that the Planck constant hbar, ℏ, is also unity so since momentum is`vec p=mℏ veck` it becomes `vec p=vec k` and the energy `E=(vec p*vec p)^2/(2m)` becomes `E=(vec k*vec k)^2/2`;
When two particles, p1 and p2, interact their momenta and energies generally change. If we let the initial k values be
`vec k_1` and `vec k_2` and the final k values be `vec k'_1` and `vec k'_2` the conservation laws are:
Momenta: `vec k'_1+vec k'_2=vec k_1+vec k_2`
Energy: `1/2(vec k'_1*vec k'_1+vec k'_2*vec k'_2)=
1/2(vec k_1*vec k_1+vec k_2*vec k_2)`
`Momenta^2: vec k'_1*vec k'_1+vec k'_2*vec k'_2+2vec k'_1*vec k'_2=vec k_1*vec k_1+vec k_2*vec k_2+2vec k_1*vec k_2`
Subtracting twice the Energy equation from the `Momenta^2` equation we have`(Momenta^2-2*"Energy:") \ \ \ 2vec k'_1*vec k'_2=2vec k_1*vec k_2`
The nature of the interaction that we will use here is that it occurs as a repelling force between the centers of any two particles and the direction will have a unit vector`hat bb"u"_12=hat bb"x"cos psi+hat bb"y"sinpsi`
where `psi` will be chosen as a random value between 0 and `2pi`. After a lot of algebra we come to the result:
`vec k'_1=vec k_1-deltaKhat bb"u"_12`
`vec k'_2=vec k_2+deltaKhat bb"u"_12`
`deltaK=hat bb"u"_12*(vec k_1-vec k_2)`
The `vec k` values of the two particles have changed so we must find the new allowed FD spatial states for them. We must therefore find un-occupied states of the same spin for both of them. If the temperature is absolute zero then all rings out to `k_(Fermi)` are fully occupied so we need to provide an additional ring that has a gap between itself and the `k_(Fermi)` ring so there become some un-occupied states in the gap as shown in Canvas 1. The learner can easily verify that no state change occurs when the ring gap is zet to zero, When random interactions start, one of the resulting 'vec k' magnitudes will be smaller than its parent and the other will be larger. In order to make a state change, the smaller one needs to find an un-occupied state with the smaller 'vec k' magnitude. That state will not usually be in the ring gap but sometimes it will. The wider the gap, the more likely it will find a state so the time to quasistatic equilibrium becomes shorter as can also be vetified by the learner. The resulting momenta from the interaction will not generally fit into any state precisely. However, it is adequate to search for the closest un-occupied ring index and azimuth index, `phi_m`, for the resulting momentum. Having found the pair of un-occupied states, the occupation of the parent states is set to a spin of 0 (meaning that the parent state is no longer occupied) and the spins of the new states is set to the spins of the parent states. Both the average energy, Eavg, and the total spin remain constant.
Now we have the task of finding the analog of temperature when we have only the Fermi-Dirac particle kinetic energy. It is obvious that temperature is a linear measure of the increase of FD particles above their kinetic energy at absolute zero. First we need the kinetic energy of all the particles in a single ring. The number of particles in occupied ring `n` is:
`N_n="Integer"(2pin)`
and all of the particles have the kinetic energy `E_n=1/2(ndeltak)^2` If all rings out to rin `n` are fully occupied, the total energy of particles out to ring n is`E_("sum")=sum_1^n(2pi i^3delta k^2)`
and for large n that sum can be given by its integral:`E_("sum")=2pin^4/4delta k^2`
The total number of particles out to ring `n` is obtained by its integral`N_("sum")=2pin^2/2`
so the average energy of particles out to ring `n_(Fermi)` is`E_(avg)=E_("sum")/N_("sum")=(n_(Fermi)^2delta k^2)/2=E_(Fermi)/2`
Their kinetic energy at absolute zero is a simple function of the Fermi energy `E_(Fermi)` which is related to the Fermi momentum `k_(Fermi)=n_(Fermi)deltak` as`E_(Fermi)=1/2 k_(Fermi)^2`
When we add the extra occupied gapped ring with gap `deltan_(gap)`to the rings at `n_F` we increase the total energy to:`E'_("sum")=2pin_(Fermi)^4/4delta k^2+2pi(n_(Fermi)+deltan_(gap))^3delta k^2`
and we increase the total number of occupied states to`N'_("sum")=2pin_(Fermi)^2/2+2pi((n_(Fermi)+deltan_(gap))`
The new average state energy then becomes:`E'(avg)=(E'_("sum"))/(N'_("sum"))`
Then the Boltzmann temperature energy is proportional to`k_BT=E'(avg)-E(avg)=E'(avg)-E_(Fermi)/2\ \ \ \ (1)`
where `k_B` is the Boltzmann constant and `T` is absolute temperature.The expression of the Fermi-Dirac temperature distribution is
`FD(E,T)=1/(exp((E-mu)/(k_BT))+1)`
where E is particle energy, T is Kelvin temperature and `mu` is the chemical potential.For FD particles the chemical potential is just the maximum energy of occupied states at absolute zero which is `E_(Fermi)`.
Then we can use the expression for `k_BT` to complete our expression for `FD`:`FD(E,E'(avg),E_(Fermi))=1/(exp(E-E_(Fermi))/(E'(avg)-E_(Fermi)/2)+1) \ \ \ \ (2)`
To get the actual density of states Vs `|k|` of the rings, equation 2 must be multiplied by the number of states in each ring.`(dN(k))/(dk)="Integer"((2pik)/(deltak))FD(E,E'(avg),E_(Fermi))\ \ \ \ \(3)`
Equation 3 is plotted in blue on Canvas 2. Its plot maximum is set to the height of Canvas 2. After pressing Start and letting the program run for a significant time, the learner should note how well the red statistical distribution agrees with the blue plot of equation 3. The statistical plot maximum is also set to the height of Canvas 2 so as to be comparable to equation 3. That is what is expected.