The thermal state corresponding to a photon mode $a$ can be written as
$\rho_{a}(\beta) = \frac{e^{-\beta a^\dagger a}}{Tr(e^{-\beta a^\dagger a})}$. Since it is easy (I suppose) to talk in terms of the photon number state, I am interested in following:
1) How can one write a similar state for two modes $a_1$ and $a_2$?, That is, what would be $\rho_{a_1a_2}(\beta)?$
2) How can these states (one mode and multimode) be represented in terms of photon number states $|n>$?
1 Answer
A general thermal state has the form $$ \rho = \frac{e^{-\beta H} }{Z} $$ where $H$ is the Hamiltonian and $Z = \mathrm{Tr}\, e^{-\beta H}$ is the partition function.
For a single photon mode $$ H = \omega a^\dagger a $$ which, for $\omega = 1$ gives the expression in your question
For a system with 2 modes the Hamiltonian is simply the sum of the Hamiltonians for each mode so $$ H =\omega_1 a_1^\dagger a_1 + \omega_2 a_2^\dagger a_2 $$ and so $$ \rho = \frac{e^{-\beta(\omega_1 a_1^\dagger a_1 + \omega_2 a_2^\dagger a_2)}}{Z} $$
To express this in terms of the photon number basis we place the resolution of the identity, $\sum_{n_1, n_2}|n_1,n_2\rangle\langle n_1,n_2|$ before and after the density matrix \begin{align} \rho &= \frac{1}{Z}\sum_{n_1, n_2, n_1',n_2'}|n_1,n_2\rangle\langle n_1,n_2| e^{-\beta(\omega_1 a_1^\dagger a_1 + \omega_2 a_2^\dagger a_2)}|n_1',n_2'\rangle\langle n_1',n_2'| \\ &= \frac{1}{Z}\sum_{n_1, n_2}|n_1,n_2\rangle e^{-\beta(\omega_1n_1+\omega_2n_2)}\langle n_1,n_2| \end{align}