Mahan's derivation of energy current of a free-particle system

Question

On p.25 of the 3rd Ed. analogously to the polarisation operator $\textbf{P}$ for particle currents $$\textbf{P}=\int\textbf{r}\rho(\textbf{r})d^3r$$ he defines an operator $$\textbf{R}_E=\frac{1}{2}\int\textbf{r}\mathcal{H}(\textbf{r})+\mathcal{H}(\textbf{r})\textbf{r}d^3r$$ He defines the energy current, $\textbf{j}_E=\partial_t\textbf{R}_E$, then claims that for a free-particle system, the energy current is

$$ \textbf{j}_E=\sum_{\textbf{p}\sigma}\textbf{v}_{\textbf{p}}\epsilon_{\textbf{p}}c_{\textbf{p}\sigma}^\dagger c_{\textbf{p}\sigma} $$

My question is how does he get to this? I think since $\partial_t$ doesn't make sense for a free particle system, he's using the commutator relation to get the time evolution, but if you plug $\psi=a_ke^{ikr}$ into the definition of $\textbf{R}_E$, using $H=\nabla^2$ the first term gives $$\frac{\hbar^2}{2m}\psi^\dagger \textbf{r}\nabla^2\psi=\frac{\hbar^2}{2m}\psi^\dagger \textbf{r}\nabla^2e^{ikr}=-\frac{\hbar^2k^2}{2m}\textbf{r}\psi^\dagger\psi$$ for the second term, we use the product rule (and reusing an old short-hand he's used previously) to find $$\frac{\hbar^2}{2m}\psi^\dagger\nabla^2(\textbf{r}\psi)=\frac{\hbar^2}{2m}\left[\psi^\dagger \textbf{r}\nabla^2\psi+2\psi^\dagger(\nabla \textbf{r}\cdot\nabla\psi)\right]=\frac{\hbar^2}{2m}\left[-\textbf{r}\psi^\dagger\psi+2\psi^\dagger\nabla\psi\right]$$ where the only difference is an additional $\frac{i\hbar^2k}{2m}\psi^\dagger\psi$ term. However, once we carry out the integral terms with $\textbf{r}$ dominate and the integral diverges. Using the alternative definition of $$\frac{\partial}{\partial t}\textbf{R}_E=\frac{1}{2}\int\textbf{r}\frac{\partial}{\partial t}\mathcal{H}(\textbf{r})+\frac{\partial}{\partial t}\mathcal{H}(\textbf{r})\textbf{r}d^3r=\textbf{j}_E$$ doesn't seem to change anything either. Even if we got to a stage where we had something to put into the commutator $\partial_t\textbf{R}_E=i[H,\textbf{R}_E]$ the number operators would commute, meaning the current would be zero. How does he get the energy current for a free particle system?

Answer 1

TLDR: $\textbf{r}\mathcal{H}$ is $\textbf{r}\psi^* H\psi$.

In more detail, this is really several questions in a trench coat. First of all,

Q: What does the energy current actually represent?

If we think about the particle current, we define the sources/sinks in a said volume as the change in mass, i.e. the continuity equation $\nabla\cdot\textbf{j}+\dot{\rho}(t)=0$. The question now becomes, how to implement this for energy? We could take the time derivative of the Hamiltonian $\dot{H}$, but operators represents measurements. If we want to find a continuity equation for energy, we would need the actual measured value of energy. As such, it makes sense to start with the expected value of energy $$\nabla\cdot\textbf{j}+\partial_t\int\psi^*H\psi d^3r=0$$ Now let us look at what this would give for a travelling wave wave-function. This corresponds to the free particle Hamiltonian $$H=\frac{-\hbar^2}{2m}\nabla^2$$ As such, by Leibnitz' rule we can move the $\partial_t$ into the integral, which gives $$\dot{\psi}^*H\psi+\psi^*H\dot{\psi}$$ from the Schroedinger equation $$i\hbar\dot{\psi}=H\psi$$ we can substitute this in to get $$\left(\frac{1}{-i\hbar}H\psi^*\right)(H\psi)+\psi^*H\left(\frac{1}{i\hbar}H\psi\right)$$ Now, we can use the eigenfunction property of $\psi,\psi^*$ to rewrite $H\psi=E\psi$ and $\psi^*H=E\psi^*$. Now, this brings us the next question

Q: Why are we treating the Hamiltonian differently based on whether it came from the Schroedinger equation or not?

The answer to this is because we are trying to find the energy current, which we know is related by $\nabla\cdot\textbf{j}$. As such, we'd expect to be able to express our RHS as an overall divergence. Since $\nabla^2u$ for a scalar function $u$ is defined as $\nabla\cdot(\nabla u)$, any higher order $\nabla$ will likely not be needed. To simplify the equation, we've then chosen to get rid of two orders of $\nabla$ by using the eigenfunction property. This then leaves us the $\nabla^2$ in the Hamiltonian, which we can guess can be split using $$\begin{split} \nabla\cdot(\psi^*\nabla\psi)&=\psi^*\nabla^2\psi+(\nabla\psi^*)\cdot(\nabla\psi)\\ \nabla\cdot(\psi\nabla\psi*)&=\psi\nabla^2\psi^*+(\nabla\psi^*)\cdot(\nabla\psi) \end{split}$$ which gives $$\nabla\cdot{\textbf{j}}=\frac{E}{i\hbar}\cdot\frac{\hbar^2}{2m}[\nabla\cdot(\psi^*\nabla\psi)-\nabla\cdot(\psi\nabla\psi^*)]$$ Since $\hat{p}=-i\hbar\nabla$, we can rewrite this as $$\frac{E}{2m}[\psi^*\hat{p}\psi+\psi\hat{p}\psi^*]$$ which is simply the real part of the expected value of the momentum. Thus, this is simply $$\epsilon_k\textbf{v}_k$$ from which we can follow (How do you get from the polarisation operator to particle current?) to reproduce the energy current from the $\int\textbf{r}\mathcal{H}d^3r$ form. This leads to our last question

Q: Does it make sense to look at the energy current for a plane-wave, given that it exists over all space?

Yes! The wavefunction gives a probability distribution, but the particle collapses to a single point after measurement. We use the wavefunction because that is where information such as particle momentum or energy is described. Although the wavefunction is spread out over all space, the particle itself has a defined momentum, so when we take the integral over a surface it gives the probability that the particle carries the energy out of that surface.