While studying many-body theory (in particular, field theory at positive temperatures) I came across the following exercise: given the grand-canonical Hamiltonian for a bosonic system \begin{equation} \hat{K}=\sum_a[(\epsilon_a-\mu)\hat{c}^\dagger_a\hat{c}_a+\gamma^*_a\hat{c}_a+\gamma_a\hat{c}^\dagger_a] \end{equation} where $\hat{c}_a$ and $\hat{c}^\dagger_a$ are the destruction and creation operators associated to a complete orthonormal system indexed by $a$, I have to calculate the expectation value \begin{equation} \langle \hat{c}_a\rangle_K= \frac{1}{Z}\operatorname{tr}(e^{-\beta \hat{K}}\hat{c}_a). \end{equation}
Denoting by $\hat{c}_a(\tau)_K=e^{\hat{K}\tau/\hbar}\hat{c}_ae^{-\hat{K}\tau/\hbar}$ the evolution in "imaginary time" of the operator, that is the solution of the Heisenberg equation \begin{equation} \hbar\frac{\mathrm{d}}{\mathrm{d}\tau}\hat{c}_a(\tau)_K=[\hat{K},\hat{c}_a(\tau)_K], \end{equation} which, if my computation is correct, equals \begin{equation} \hat{c}_a(\tau)_K= e^{-(\epsilon_a-\mu)\tau/\hbar}\hat{c}_a-\gamma_a\frac{\tau}{\hbar}\hat{1}, \end{equation} I can write \begin{equation} \langle\hat{c}_a\rangle_K= \frac{1}{Z}\operatorname{tr}(e^{-\beta \hat{K}}\hat{c}_a)= \frac{1}{Z}\operatorname{tr}(\hat{c}_ae^{-\beta \hat{K}})= \frac{1}{Z}\operatorname{tr}(e^{-\beta \hat{K}}e^{\beta \hat{K}}\hat{c}_ae^{-\beta \hat{K}})= \frac{1}{Z}\operatorname{tr}(e^{-\beta \hat{K}}\hat{c}_a(\hbar\beta)_K) \end{equation} and at the same time \begin{equation} \langle\hat{c}_a\rangle_K= \frac{1}{Z}\operatorname{tr}(e^{-\beta \hat{K}}\hat{c}_a)= \frac{1}{Z}\operatorname{tr}(e^{-\beta \hat{K}}\hat{c}_ae^{\beta \hat{K}}e^{-\beta \hat{K}})= \frac{1}{Z}\operatorname{tr}(\hat{c}_a(-\hbar\beta)_Ke^{-\beta \hat{K}})= \frac{1}{Z}\operatorname{tr}(e^{-\beta \hat{K}}\hat{c}_a(-\hbar\beta)_K). \end{equation} This is absurd, since \begin{equation} \hat{c}_a(\hbar\beta)_K-\hat{c}_a(-\hbar\beta)_K= -2\sinh\bigl(\beta(\epsilon_a-\mu)\bigr)\hat{c}_a-2\beta\gamma_a\hat{1} \end{equation} therefore \begin{equation} \operatorname{tr}\bigl(e^{-\beta \hat{K}}[\hat{c}_a(\hbar\beta)_K-\hat{c}_a(-\hbar\beta)_K]\bigr)\ne 0 \end{equation} while it should be zero, following from the equations above. I must have made some mistake, but I can't see where... I need some external advice!
