No estoy convencido de que el valor medio del conjunto de un operador de aniquilación no desaparezca en el gran conjunto canónico. Aquí está mi argumento de por qué se desvanece:
Sea $H$ sea el Hamiltoniano y $\mathscr N = \sum_{i=0}^\infty a_i^\dagger a_i$ sea el operador numérico. Sea $$ |\mathbf n \rangle = |n_0, n_1, \dots\rangle $$ sea la base numérica de ocupación que satisface \begin{align} H|\mathbf n \rangle &= E_{\mathbf n}|\mathbf n\rangle,\qquad E_\mathbf n = \sum_{i=0}^\infty \epsilon_in_i \\ \mathscr N|\mathbf n \rangle &= N_{\mathbf n}|\mathbf n\rangle,\qquad N_\mathbf n = \sum_{i=0}^\infty n_i \end{align} y \begin{align} a_i|\mathbf n\rangle &= \sqrt{n_i}|n_0, \dots, n_{i-1}, n_i-1,n_{i+1}, \dots\rangle \\ a^\dagger_i|\mathbf n\rangle &= \sqrt{n_i+1}|n_0, \dots, n_{i-1}, n_i+1,n_{i+1}, \dots\rangle \end{align} La media del conjunto de cualquier operador de aniquilación $a_i$ es \begin{align} \frac{1}{Z}\mathrm{tr}(e^{-\beta(H-\mu \mathscr N)} a_i) &= \frac{1}{Z}\sum_{\mathbf n}\langle \mathbf n| e^{-\beta(H-\mu \mathscr N)}a_i|\mathbf n\rangle\\ &= \frac{1}{Z}\sum_{\mathbf n}e^{-\beta(E_\mathbf n-\mu N_\mathbf n)}\langle \mathbf n |n_0, \dots, n_{i-1}, n_i-1,n_{i+1}, \dots\rangle\\ &= 0 \end{align} donde la última igualdad se deduce de la ortogonalidad de los vectores base $|\mathbf n\rangle$ .
Avísame si encuentras algún error.