Thermodynamics+of+Harmonic+Bosons

toc =Introduction=
 * Class Session 07 : Thermodynamics of Harmonic Bosons **

Bosons are indistinguishable particles, which implies that the wave function describing an ensemble of bosons is //symmetric// with respect to any permutation of the individual positions of the bosons. Consequently, the density matrix of //N// non-interacting bosons does not simply factorise into a product of //N// single particle density matrices.

In this Class Session we examine some properties of the partition function of non-interacting bosonic particles.

=The partition function of //N// non-interacting bosons=

Canonical partition function of bosonic particles
math Z_N= \frac 1{N!}\sum_P P\!Z_N^{\text{dist}} math > where the sum is over the //N//! permutations //P// of (1,...,//N//), and //P////Z N // dist is the the partition function of //N// distinguishable particles where the final positions has been permuted by //P//: math P\!Z_N^{\text{dist}} = \int d^N x\ \rho^{\text{dist}}(\{x_1,\ldots,x_N\},\{x_{P(1)},\ldots,x_{P(N)}\},\beta) math > In this expression, //ρ// dist is the density matrix of //N// distinguishable particles at inverse temperature is //β// = 1 / //T//. math Y_N = \sum_P \text{weight}(P) math > Any permutation //P// of (1,...,//N//) can be decomposed as a unique product of commuting cycle permutations. By partitionning the sum Σ // P // over the length //k// of the cycle to which the last element of the permutation //P// belongs, show the following recursion relation: math Y_N = \sum_{k=1}^N z_k \binom{N-1}{k-1} (k-1)!\: Y_{N-k} \qquad \text{with}\quad Y_0=1 math > where //z k // is the single particle partition function at inverse temperature //kβ//. math Z_N= \frac 1N \sum_{k=1}^N z_k Z_{N-k} \qquad\qquad \qquad\qquad \qquad\qquad \qquad(1) math > with //Z// 0 = 1. math \pi_k = \frac 1N \frac {z_k}{Z_{N}}{Z_{N-k}} \qquad\qquad \qquad\qquad \qquad\qquad \qquad(2) math
 * Write the density matrix of //N// distinguishable particles described by their positions (//x// 1 ,...,//x N //) in terms of the single particle density matrix.
 * Explain why the same relation does not hold for bosonic particles.
 * Find a relation between the density matrix of //N// bosonic particles and the density matrix of //N// distinguishable particles.
 * Let's now consider the canonical partition function //Z N // of //N// bosonic (//i.e.// indistinguishable) particles described by positions (//x// 1 ,...,//x N //). Show that //Z N // writes as:
 * We note that the partition function writes //Z N // writes as //Z N //= //Y N // / //N//! with //Y N // a weight which writes as
 * Show now [Borrmann and Franke 1993] the following recursion relation for the partition //Z N //:
 * In other words, the probability //π k // that the last particle belongs to a cycle of length //k// in the permutation //P// is:

The case of bosons in an harmonic trap
We chose a harmonic trap the such that the energy levels in each of the three spatial directions of the trap are //E n // = //n//. The (three-dimensional) single particle partition //z k // writes: math z_k= \Big(\frac{1}{1-e^{-k\beta}}\Big)^3 math We define a reduced temperature //T// ٭ as (see also SMAC p.202): math T_\star = \frac{T}{N^{1/3}} math

Algorithm
code format="python" from pylab import * from math import sqrt, sinh, tanh, exp from random import uniform as ran, gauss def z(beta,k): sum = 1 / (1 - exp(-k*beta))**3 return sum def canonic_recursion(beta,N): Z = [1.] for M in range(1, N+1): Z.append( sum(Z[k] * z(beta,M-k) for k in range(M)) / M ) return Z def pi_list_make(Z,M): pi_list = [0] + [ z(beta,k)*Z[M-k]/Z[M]/M for k in range(1,M+1) ] pi_sum = [0] for k in range(1, M+1): pi_sum.append(pi_sum[k-1] + pi_list[k]) return pi_sum def tower_sample(data,Upsilon): #naive tower sampling, cf. SMAC Sect. 1.2.3  for k in range(len(data)): if Upsilon < data[k]: break return k def levy_harmonic_path(Del_tau,N): # beta = N * Del_tau xN = gauss(0.,1./sqrt(2*tanh(beta/2.))) x = [xN] for k in range(1,N): Upsilon_1 = 1./tanh(Del_tau) + 1./tanh((N-k)*Del_tau) Upsilon_2 = x[k-1]/sinh(Del_tau) + xN/sinh((N-k)*Del_tau) x_mean = Upsilon_2 / Upsilon_1 sigma = 1. / sqrt(Upsilon_1) x.append(gauss(x_mean,sigma)) return x

N = 512 T_star = .1 T = T_star * N**(1./3.) beta = 1. / T Z = canonic_recursion(beta,N) M = N x_config = [] y_config = [] while M > 0: pi_sum = pi_list_make(Z,M) Upsilon = ran(0,1) k = tower_sample(pi_sum,Upsilon) x_config += levy_harmonic_path(beta,k) y_config += levy_harmonic_path(beta,k) M -= k figure(1) axis('scaled') xlim(-5.,5.) ylim(-5.,5.) plot(x_config,y_config,'ro') xlabel('$x$') ylabel('$y$') title('3d bosons in trap (2d projection):'+' $N$= '+str(N)+' T* =' + str(T_star)) show code
 * Familiarize yourself with the above algorithm, and run it for values of temperatures of //T// ٭ . What is the advantage of defining //T// ٭ ?
 * Modify the program such that, together with the scatter plot of particles, it also draws the longest cycle in the configuration. Plot snapshots at various values of //T// ٭.
 * Modify the program so that it computes the cycle-length distribution at the values of //T// ٭ you studied before. Compare the values of a Monte Carlo sampling, and a direct plot of the distribution using equation (1) above. Show that the cycle distribution has a flat part, and comment.
 * Text books of theoretical physics tell us that the quantum nature of bosonic fluids becomes noticeable when the inter-particle distance is of the order of the thermal de Broglie wavelength. We would rather say that the quantum nature of bosonic fluids becomes noticeable when long cycles appear. Why are the two statements closely related?

Modified Algorithm
The following program enables to compare the theoretical and numerical distribution of cycles: code format="python" from pylab import * from math import sqrt, sinh, tanh, exp from random import uniform as ran, gauss def z(beta,k): sum = 1 / (1 - exp(-k*beta))**3 return sum def canonic_recursion(beta,N): Z = [1.] for M in range(1, N+1): Z.append( sum(Z[k] * z(beta,M-k) for k in range(M)) / M ) return Z def pi_list_make(Z,M): pi_list = [0] + [ z(beta,k)*Z[M-k]/Z[M]/M for k in range(1,M+1) ] pi_sum = [0] for k in range(1, M+1): pi_sum.append(pi_sum[k-1] + pi_list[k]) return pi_sum def tower_sample(data,Upsilon): #naive tower sampling, cf. SMAC Sect. 1.2.3  for k in range(len(data)): if Upsilon < data[k]: break return k def levy_harmonic_path(Del_tau,N): # beta = N * Del_tau xN = gauss(0.,1./sqrt(2*tanh(beta/2.))) x = [xN] for k in range(1,N): Upsilon_1 = 1./tanh(Del_tau) + 1./tanh((N-k)*Del_tau) Upsilon_2 = x[k-1]/sinh(Del_tau) + xN/sinh((N-k)*Del_tau) x_mean = Upsilon_2 / Upsilon_1 sigma = 1. / sqrt(Upsilon_1) x.append(gauss(x_mean,sigma)) return x

N = 512 T_star = .1 T = T_star * N**(1./3.) beta = 1. / T Z = canonic_recursion(beta,N) M = N x_config = [] y_config = [] while M > 0: pi_sum = pi_list_make(Z,M) Upsilon = ran(0,1) k = tower_sample(pi_sum,Upsilon) x_config += levy_harmonic_path(beta,k) y_config += levy_harmonic_path(beta,k) M -= k figure(1) axis('scaled') xlim(-5.,5.) ylim(-5.,5.) plot(x_config,y_config,'ro') xlabel('$x$') ylabel('$y$') title('3d bosons in trap (2d projection):'+' $N$= '+str(N)+' T* =' + str(T_star)) show code



url}?f=print| [Print this page]

=References=
 * Borrmann P., Franke G. (1993) Recursion formulas for quantum statistical partition functions, Journal of Chemical Physics 98, 2484–2485
 * Landsberg P. T. (1961) Thermodynamics with quantum statistical illustrations, Interscience Publishers
 * SMAC part 4.2.3 - 4.2.6