The Fermi-Dirac distribution function.

[Jezuz]

Hi.
Does anyone know if it is possible to start from the thermal density matrix

\hat \rho_T = \frac{e^{-\hat H_0/kT}}{\mathrm{Tr}e^{-\hat H_0/kT}}

and from that derive that the single particle density matrix can be written as

\rho(p ; p') = \delta_{p,p'} f(\epsilon_p)

just by using the antisymmetry of fermions?

I have tried to start with the many particle density matrix

\rho^{(N)}(p_1, p_2 ,... p_N ; p&#039;_1, p&#039;_2 ,... p&#039;_N) = <br /> \left&lt; p_1 \wedge p_2 ... p_N \right| \hat \rho_T <br /> \left| p&#039;_1 \wedge ... p&#039;_N \right&gt;

and traced out all but one particle momenta, but this gives me quite difficult sums that are dependent on each other.

Here \hat H_0 = \sum_{i=1}^N \hat p^2_i /2m where N is the number of fermions. \epsilon_p = p^2/2m and f is the Fermi-Dirac distribution function. The states \left| p_1 \wedge p_2 ,... p_N\right&gt; are totally antisymmetric states.

 

[olgranpappy]

Hi.
Does anyone know if it is possible to start from the thermal density matrix

\hat \rho_T = \frac{e^{-\hat H_0/kT}}{\mathrm{Tr}e^{-\hat H_0/kT}}

and from that derive that the single particle density matrix can be written as

\rho(p ; p&#039;) = \delta_{p,p&#039;} f(\epsilon_p)

just by using the antisymmetry of fermions?

Try and do it for a single fermionic oscillator first:
<br /> Z=1+e^{-\beta \omega}<br />
<br /> \rho_{0,0}=1/Z<br />
<br /> \rho_{1,1}=\frac{1}{e^{\beta\omega}+1}<br />
<br /> \rho_{1,0}=\rho_{0,1}=0<br />

 

[Jezuz]

I tried that but when the momentum basis is used the partition function is
Z = \sum_p \exp(-p^2/2mkT).
I don't know how to evaluate this sum to get the right result. Of course I can assume that the particle is in a large volume V and convert the sum to an integral with a volume factor V/(2\pi \hbar)^3. This yields a gaussian integral and solving this yields a result quite different from the required result.

The problem becomes even harder when the full N-body problem is considered.

 

[olgranpappy]

I tried that but when the momentum basis is used the partition function is
Z = \sum_p \exp(-p^2/2mkT).
I don't know how to evaluate this sum to get the right result.

Well, these are two different problems. The example I gave was for a single fermionic oscillator. The example you give above is for a single free particle. You won't end up with a fermi-dirac distribution.

The exact problem that you are interested in is worked out in detail in a very cute and original way in the text "Solid State Physics" by Ashcroft and Mermin. Pages 40 thru 42.

 

[Jezuz]

Ok. It's been a while since I was working on this, but now I am giving it another go. I will try to be more specific about what I am asking exactly. The problem is:

Consider a system of N particles in thermal equilibrium at temperature T. The density operator is then given by:
\hat \rho_T = \frac{e^{-\hat H_0/kT}}{\mathrm{Tr}e^{-\hat H_0/kT}}.

Here \hat H_0 = \hat{\mathbf p}^2/2m.

Derive that the single particle density matrix is given by

\rho(\mathbf{p} ; \mathbf{p&#039;}) = \delta_{\mathbf{p},\mathbf{p}&#039;} f(\epsilon_p).

Ok. We all have seen the derivation of this using the Grand Canoncial Ensemble, but here I want to do it in a rather straight forward way instead.

I want to do it using the representation of the particles in momentum space and where the states are totally antisymmetric (to account for the fact that they're fermions). This means that the state of the particles are given by

\left| \mathbf p_1 \wedge \mathbf p_2 \wedge ... \mathbf p_N \right&gt; = <br /> \frac{1}{N!}\sum_P (-1)^{\sigma_P} \left| \mathbf p_1 \right&gt; \otimes \left|\mathbf p_2 \right&gt; ... \left| \mathbf p_N \right&gt;,

where the sum is over all possible permutations P of the indices and \sigma_P is the number of transpositions of the permutation.

From this you can prove easily prove that

\left&lt; \mathbf p_1 \wedge \mathbf p_2 ... \mathbf p_N \right| \left. \mathbf p&#039;_1 \wedge \mathbf p&#039;_2 ... \mathbf p&#039;_N \right&gt;

is \pm 1 if \{ \mathbf p_i \} is a permutation of \{ \mathbf p&#039;_i \} and 0 othervice.

The density matrix for the many particle state is given by:

\rho(\mathbf p_1, \mathbf p_2 ... \mathbf p_N ; \mathbf p_1 , \mathbf p_2 ... \mathbf p_N) = \left&lt; \mathbf p_1 \wedge \mathbf p_2 ... \mathbf p_N \right| <br /> \hat \rho \left| \mathbf p&#039;_1 \wedge \mathbf p&#039;_2 ... \mathbf p&#039;_N \right&gt;

where \hat \rho is the density operator. The trace of the density matrix then is given by

\mathrm{Tr} \hat \rho = \sum_{\mathbf p_1, \mathbf p_2 , ... , \mathbf p_N} <br /> \left&lt; \mathbf p_1 \wedge \mathbf p_2 ... \mathbf p_N \right| \hat \rho \left| \mathbf p_1 \wedge \mathbf p_2 ... \mathbf p_N \right&gt; = N! \sum_{|\mathbf p_1| &lt; |\mathbf p_2 |&lt; ... &lt; |\mathbf p_N| } \left&lt; \mathbf p_1 \wedge \mathbf p_2 ... \mathbf p_N \right| \hat \rho \left| \mathbf p_1 \wedge \mathbf p_2 ... \mathbf p_N \right&gt;,

where the last equality follows quite easily from the properties of the antisymmetric states.

Ok, this is what I have to use. From this I would like to calculate the reduced density matrix defined by

\rho(\mathbf p; \mathbf p&#039;) = N \sum_{\mathbf p_2, \mathbf p_3,...\mathbf p_N} <br /> \rho(\mathbf p , \mathbf p_2, \mathbf p_3, ... , \mathbf p_N; \mathbf p&#039; , \mathbf p_2, \mathbf p_3 ,... \mathbf p_N)

which is just tracing out N-1 particles (does not matter which particle I keep fixed since the state are antisymmetric).

When calculating \hat \rho_T the nominator is not a problem. It will easily be calculated to \delta_{\mathbf p, \mathbf p&#039;} e^{-p^2/2mkT}.

For the denominator, however, I run into problems. This is what I have done.

\textrm{Tr} e^{-H_0/kT} = \sum_{\mathbf p_1, \mathbf p_2, ... \mathbf p_N} \left&lt; \mathbf p_1 \wedge ... \mathbf p_N\right| e^{-\sum_{i=1}^N p_i^2/2mkT} \left| \mathbf p_1 \wedge \mathbf p_2 ... \mathbf p_N \right&gt; = <br /> N! \sum_{|\mathbf p_1|&lt; |\mathbf p_2| &lt; ... |\mathbf p_N| } e^{-\sum_i p_i^2/2mkT} \left&lt; \mathbf p_1 \wedge ... \mathbf p_N \right| \left. \mathbf p_1 \wedge ... \mathbf p_N \right&gt; .

Now I can use that the inner product between the states is zero and further transform the sums into integrals:

\textrm{Tr} \hat e^{-\hat H_0/kT} = N! \left(\frac{4\piV}{(2\pi \hbar)^3} \right)^N \int_0^\infty p_1^2 dp_1 \int_{p_1}^\infty p_2^2 dp_2 ... \int_{p_{N-1}}^\infty p^2_N dp_N e^{-\sum_i p_i^2/2mkT} = \left(\frac{4\pi V}{(2\pi \hbar)^3} \right)^N \left( \int_0^\infty p^2 dp e^{-p^2/2mkT} \right)^N<br />.

Here I have used that \lim_{V\rightarrow \infty} \frac{1}{V} \sum_{\mathbf p} = <br /> \int \frac{d \mathbf p}{(2\pi\hbar)^2} and then integrated out all the angular dependence to give the factor 4 \pi. Further on I have extended all the lower integration limits to zero and divided by N! to account for the extra contribution. This integral can now be solved to get

\textrm{Tr} \hat e^{-\hat H_0/kT} = \left(\frac{4\pi V}{(2 \pi \hbar)^3} \right)^N \left( \frac{\pi^{1/2} (2mkT)^{3/2} }{4}\right)^N.

This, however, does not look anything like the required result. Can someone please help and point out where I make a mistake?

 

[Jezuz]

Ok. I think I have found my error. When I go from the sums to integrals I actually no longer account for the exclusion principle and the result I get in the end is the classical Maxwell-Boltzmann distribution.

Apparently the trace in the denominator cannot be calculated exactly. Instead one has to use the fact that the denominator is the partition function and then use the relation between the partition function, the free energy and the chemical potential.