$N$ bosons in two energy levels What is the partition function for $N$ bosons in a two state system with $E_1=0$ and $E_2= E$?  I know that bosons don't obey the Pauli exclusion principle and they are indistinguishable but I am having doubts when it comes to the fact they are indistinguishable so I have to exclude cases where $x$ amount of particles are in state $E_1$ and $y$ amount particles in $E_2$ is the same thing that $y$ amount of particles are in $E_1$ and $x$ amount are in $E_2$. How do take that in account when writing the expression?
 A: For a system of $N$ indistinguishable bosons and two single-particle states (e.g. lattice sites), we have $N+1$ basis states that we can enumerate as follows:
$\{|N-m,m\rangle\}_m$, where $m=0,\ldots,N$. In the basis state $|N-m,m\rangle$ there are $N-m$ bosons in the single-particle state $1$ and $m$ bosons in the single-particle state $2$. The Hamiltonian
$$\hat{H} = E_1\, \hat{n}_1 + E_2\, \hat{n}_2 \quad,$$
describes a system of non-interacting bosons. Here, the operator $\hat{n}_i$ counts the number of bosons in the single-particle state $i=1,2$. Note that in particular the basis introduced above is an eigenbasis of this Hamiltonian. We define the partition function associated with $\hat{H}$ as
$$ Z\equiv \mathrm{Tr}\,e^{-\beta\hat{H}} \quad .$$
Using the above considerations, we obtain
$$ Z= e^{-\beta N E_1}\,\sum\limits_{m=0}^{N}\left(e^{-\beta (E_2-E_1)}\right)^m \quad .\tag{1}\label{1}$$
We note that the $n$-th partial sum of the geometric series for $q\neq1$ is given by:
$$ s_n \equiv \sum\limits_{k=0}^n q^k= \frac{1-q^{n+1}}{1-q} \quad .$$
Substituting $n\rightarrow N$, $q\rightarrow e^{-\beta (E_2-E_1)}$ and $k\rightarrow m$ eventually yields
$$Z = e^{-\beta N E_1} \, \left(\frac{1-e^{-\beta (E_2-E_1)(N+1)}}{1-e^{-\beta (E_2-E_1)}}\right) \quad. $$
Edit: Regarding your comment: How do we obtain equation $(1)$?
First, note that the energy of a basis state $|N-m,m\rangle$ is given by
$$ \epsilon_m\equiv \langle N-m,m|\hat{H}|N-m,m\rangle \quad.$$
This yields $\epsilon_m = E_1\, (N-m) + E_2 \,m$. These are the eigenvalues of the Hamiltonian.
Second, the partition function is defined as the trace of the exponentiated Hamiltonian. Roughly speaking, this means that we must sum over all expectation values of this exponentiated Hamiltonian, i.e. we have to calculate:
$$Z = \sum\limits_{m=0}^N \langle N-m,m|e^{-\beta\hat{H}}|N-m,m\rangle \quad.  $$ This can be rewritten to (this should be covered in any textbook of quantum mechanics or statistical physics)
$$ Z = \sum\limits_{m=0}^N e^{-\beta \epsilon_m} \quad .$$
In other words, the partition function is given by the sum of the exponentiated eigenvalues of the Hamiltonian.
As a final step, replace the expression of $\epsilon_m$ with the value we have derived earlier. Then use simple rules of the exponential function to obtain equation $(1)$.
A remark: Actually, the trace of an operator is independent of the basis we choose. But the basis chosen is of course the most convenient one.
