Derivation of the "Bethe sum rule" I am trying to work out the steps of the proof of the expression: $$\sum_n (\mathcal{E_n}-\mathcal{E_s})|\langle n|e^{i\mathbf{q}\cdot\mathbf{r}}|s \rangle|^2 = \frac{\hbar^2q^2}{2m}$$ from Eq. (5.48) in the book Principles of the Theory of Solids by Ziman. In the book it is mentioned that this can be shown by expanding out: $$[[\mathcal{H},e^{i\mathbf{q}\cdot\mathbf{r}}],e^{-i\mathbf{q}\cdot\mathbf{r}}]$$ where $$\mathcal{H} = -\frac{\hbar^2}{2m}\nabla^2+\mathcal{V}(\mathbf{r})$$ with $\mathcal{V}(\mathbf{r})$ being periodic. What I did is a simple expansion: $$[[\mathcal{H},e^{i\mathbf{q}\cdot\mathbf{r}}],e^{-i\mathbf{q}\cdot\mathbf{r}}] = 2\mathcal{H}-e^{i\mathbf{q}\cdot\mathbf{r}}\mathcal{H}e^{-i\mathbf{q}\cdot\mathbf{r}}-e^{-i\mathbf{q}\cdot\mathbf{r}}\mathcal{H}e^{i\mathbf{q}\cdot\mathbf{r}}$$ Then I took the inner product with the eigenstates of the Hamiltonian $\mathcal{H}|s\rangle = E_s |s\rangle$ to get $$\langle s|[[\mathcal{H},e^{i\mathbf{q}\cdot\mathbf{r}}],e^{-i\mathbf{q}\cdot\mathbf{r}}]|s\rangle = 2E_s - \langle s|e^{i\mathbf{q}\cdot\mathbf{r}}\mathcal{H}e^{-i\mathbf{q}\cdot\mathbf{r}}|s\rangle - \langle s|e^{-i\mathbf{q}\cdot\mathbf{r}}\mathcal{H}e^{i\mathbf{q}\cdot\mathbf{r}}|s\rangle$$ Now, what's obstructing my calculation is the fact that I cannot justify the last two terms in the above expression being equal. I really need them to be equal to show the top identity (also known as the Bethe sum rule). The main obstacle is the fact that $e^{i\mathbf{q}\cdot\mathbf{r}}$ is non-Hermitian.
I have found this identity in many books and journal articles. But I cannot find a satisfactory proof anywhere. One such example is this article:

Sanwu Wang. Generalization of the Thomas-Reiche-Kuhn and the Bethe sum rules. Phys. Rev. A 60 no. 1, pp. 262–266 (1999).

Th result is stated in Eq. (3) and proof is given in section III. A. At the end of the proof they set $F = e^{i.\mathbf{q}.\mathbf{r}}$ to recover Eq. (3). The fact that they did not assume any form for $F$ means it must hold for any function. There is, however, one step that I cannot justify. In Eq. (9) how can they write: $$\langle 0|F(x)|l \rangle = \langle l|F(x)|0 \rangle$$ If I can even justify the above equality (for my case) then I'm set.
 A: Start with the expression
$$\sum_n (\mathcal{E_n}-\mathcal{E_s})|\langle n|e^{i\mathbf{q}\cdot\mathbf{r}}|s \rangle|^2 =
\sum_n (\mathcal{E_n}-\mathcal{E_s})\langle s|e^{-i\mathbf{q}\cdot\mathbf{r}}|n \rangle\langle n|e^{i\mathbf{q}\cdot\mathbf{r}}|s \rangle. $$ 
The first trick is to realize that $(\mathcal{E_n}-\mathcal{E_s})e^{i\mathbf{q}\cdot\mathbf{r}}$ equals $\pm[\mathcal{H},e^{i\mathbf{q}\cdot\mathbf{r}}]$ when it is inside either of the two brackets:
$$\sum_n (\mathcal{E_n}-\mathcal{E_s})|\langle n|e^{i\mathbf{q}\cdot\mathbf{r}}|s \rangle|^2 
=
\sum_n \langle s|e^{-i\mathbf{q}\cdot\mathbf{r}}|n \rangle\langle n|[\mathcal{H},e^{i\mathbf{q}\cdot\mathbf{r}}]|s \rangle
=
\langle s|e^{-i\mathbf{q}\cdot\mathbf{r}}[\mathcal{H},e^{i\mathbf{q}\cdot\mathbf{r}}]|s \rangle, $$ 
by summing $|n \rangle\langle n|$ to $1$. Analogously, you can do this on the first factor to get
$$\sum_n (\mathcal{E_n}-\mathcal{E_s})|\langle n|e^{i\mathbf{q}\cdot\mathbf{r}}|s \rangle|^2 
=-
\langle s|[\mathcal{H},e^{-i\mathbf{q}\cdot\mathbf{r}}]e^{i\mathbf{q}\cdot\mathbf{r}}|s \rangle. $$

You now need to calculate the commutator. Since $\mathcal{V}(\mathbf{r})$ commutes with $e^{\pm i\mathbf{q}\cdot\mathbf{r}}$, you only need to worry about the kinetic term. Thus
$$
[\mathcal{H},e^{\pm i\mathbf{q}\cdot\mathbf{r}}]=\frac{1}{2m}[\mathbf{p}^2,e^{\pm i\mathbf{q}\cdot\mathbf{r}}]
=\frac{1}{2m}\mathbf{p}\cdot[\mathbf{p},e^{\pm i\mathbf{q}\cdot\mathbf{r}}]
+\frac{1}{2m}[\mathbf{p},e^{\pm i\mathbf{q}\cdot\mathbf{r}}]\cdot\mathbf{p},
$$
so
$$
[\mathcal{H},e^{\pm i\mathbf{q}\cdot\mathbf{r}}]=
\pm\frac{1}{2m}\mathbf{p}\cdot(\hbar \mathbf{q}e^{\pm i\mathbf{q}\cdot\mathbf{r}})
\pm\frac{1}{2m}(\hbar \mathbf{q}e^{\pm i\mathbf{q}\cdot\mathbf{r}})\cdot\mathbf{p}
=\pm\frac{\hbar}{2m}\mathbf{q}\cdot\left(e^{\pm i\mathbf{q}\cdot\mathbf{r}}\mathbf{p}+\mathbf{p}e^{\pm i\mathbf{q}\cdot\mathbf{r}}\right).
$$
To phrase this just right, you need to have $e^{+i\mathbf{q}\cdot\mathbf{r}}$ to the left or $e^{- i\mathbf{q}\cdot\mathbf{r}}$ to the right. Using the commutator $[\mathbf{p},e^{\pm i\mathbf{q}\cdot\mathbf{r}}]=\pm \hbar \mathbf{q}e^{\pm i\mathbf{q}\cdot\mathbf{r}}$, as above, you get
$$
[\mathcal{H},e^{+ i\mathbf{q}\cdot\mathbf{r}}]=
+\frac{\hbar}{2m}\mathbf{q}\cdot\left(e^{+ i\mathbf{q}\cdot\mathbf{r}}\mathbf{p}+\mathbf{p}e^{+ i\mathbf{q}\cdot\mathbf{r}}\right)
=
\frac{\hbar}{2m}e^{+ i\mathbf{q}\cdot\mathbf{r}}\mathbf{q}\cdot\left(2\mathbf{p}+\hbar\mathbf{q}\right)
$$
and
$$
[\mathcal{H},e^{- i\mathbf{q}\cdot\mathbf{r}}]=
-\frac{\hbar}{2m}\mathbf{q}\cdot\left(e^{- i\mathbf{q}\cdot\mathbf{r}}\mathbf{p}+\mathbf{p}e^{- i\mathbf{q}\cdot\mathbf{r}}\right)
=-
\frac{\hbar}{2m}\mathbf{q}\cdot\left(2\mathbf{p}+\hbar\mathbf{q}\right)e^{- i\mathbf{q}\cdot\mathbf{r}}.
$$

Putting either of these into the corresponding formula above, you get
$$\sum_n (\mathcal{E_n}-\mathcal{E_s})|\langle n|e^{i\mathbf{q}\cdot\mathbf{r}}|s \rangle|^2 
=
\frac{\hbar}{2m}\langle s|\mathbf{q}\cdot\left(2\mathbf{p}+\hbar\mathbf{q}\right)|s \rangle
=
\frac{\hbar^2}{2m}\mathbf{q}^2+\frac{\hbar}{m}\langle s|\mathbf{q}\cdot\mathbf{p}|s \rangle. $$ 
You can then always force the mean value of $\mathbf{p}$ to vanish.
A: I) Let us work in the Schrödinger picture. Assume that the Hamiltonian operator 
$$\tag{1} \hat{H}~=~T(\hat{\bf p})+ V(\hat{\bf r}), \qquad T(\hat{\bf p})~:=~ \frac{\hat{\bf p}^2}{2m}, $$ 
does not depend explicit on time. In particular, the point spectrum of energy eigenvalues $E_n$ does not depend on time. Let $\psi_n({\bf r})$ be the corresponding energy eigen-functions in position space, which satisfies the TISE 
$$\tag{2} \hat{H}\psi_n({\bf r})~=~E_n \psi_n({\bf r}).$$ 
We may assume that $\psi_n({\bf r})\in{\mathbb{R}}$ is a real function, see also this Phys.SE post.
Let us consider a fixed time $t_0$. Let us pick a basis of time-dependent solutions
$$\tag{3} \langle n, t |{\bf r} \rangle~=~  \Psi_n({\bf r},t)~=~ \psi_n({\bf r})e^{-\frac{i}{\hbar}E_n (t-t_0)} $$
that becomes real at the time $t=t_0$. It follows that for any (possibly complex) function $F=F({\bf r})$ of position ${\bf r}$, that the matrix element
$$\tag{4} \langle n, t_0 |F(\hat{\bf r}) | m, t_0 \rangle ~=~ \int \! d^3r ~ \psi^{*}_n({\bf r}) ~F({\bf r})~\psi_m({\bf r})~=~ n \leftrightarrow m  $$
is symmetric in $n\leftrightarrow m$ at $t=t_0$. The time dependence of the matrix element is just a phase factor. Hence the
absolute value  
$$\tag{5} |\langle n, t |F(\hat{\bf r}) | m, t \rangle|~=~ n \leftrightarrow m  $$
of the matrix element is symmetric in $n\leftrightarrow m$ at any time $t$.
II) Next
$$2T(\hbar{\bf k})~\stackrel{(1)}{=}~T(\hat{\bf p}+\hbar{\bf k})+T(\hat{\bf p}-\hbar{\bf k})-2T(\hat{\bf p}) $$
$$~=~e^{-i{\bf k}\cdot\hat{\bf r}}T(\hat{\bf p}) e^{i{\bf k}\cdot\hat{\bf r}} +e^{i{\bf k}\cdot\hat{\bf r}}T(\hat{\bf p}) e^{-i{\bf k}\cdot\hat{\bf r}} -2T(\hat{\bf p}) $$
$$\tag{6}  ~=~e^{-i{\bf k}\cdot\hat{\bf r}}\hat{H} e^{i{\bf k}\cdot\hat{\bf r}} +e^{i{\bf k}\cdot\hat{\bf r}}\hat{H} e^{-i{\bf k}\cdot\hat{\bf r}}-2\hat{H} .$$
In the last step we used that $[\hat{\bf r},V(\hat{\bf r})]=0$.
Hence we get the Bethe sum rule
$$2T(\hbar{\bf k}) ~=~ 2\langle m, t |T(\hbar{\bf k}) | m, t \rangle $$
$$~\stackrel{(6)}{=}~\langle m, t |e^{-i{\bf k}\cdot\hat{\bf r}}\hat{H} e^{i{\bf k}\cdot\hat{\bf r}}| m, t \rangle +\langle m, t | e^{i{\bf k}\cdot\hat{\bf r}}\hat{H} e^{-i{\bf k}\cdot\hat{\bf r}}| m, t \rangle-2\langle m, t |\hat{H}| m, t \rangle $$
$$~=~\sum_{n} (E_n-E_m) |\langle n, t |e^{i{\bf k}\cdot\hat{\bf r}}| m, t \rangle |^2 + \sum_{n} (E_n-E_m) |\langle m, t |e^{i{\bf k}\cdot\hat{\bf r}}| n, t \rangle |^2$$
$$\tag{7} ~\stackrel{(5)}{=}~2\sum_{n} (E_n-E_m) |\langle n, t |e^{i{\bf k}\cdot\hat{\bf r}}| m, t \rangle |^2 .$$
A: Using $[p,f(r)] = - i \hbar \nabla f(r)$ we find
\begin{align}
[[p^2,e^{i q \cdot r}],e^{-i q \cdot r}] 
&=  2\hbar q \cdot[ p,e^{-i q \cdot r}] 
\\ &= - 2 \hbar^2 q^2 \; .
\end{align}
On the other hand,
$$[[H,e^{i q \cdot r}],e^{-i q \cdot r}]= 2H -e^{-i q \cdot r}He^{iq \cdot r} -e^{-i q \cdot r}He^{iq \cdot r} $$
Also note that,
$$1 = \left<s|e^{-iq \cdot r}e^{+iq \cdot r} |s\right>= \sum_n  \left<s \right|e^{-iq \cdot r}\left|n\right>\left<n \right|e^{+iq \cdot r} \left|s \right> =  \sum_n \left| \left<n \right|e^{iq \cdot r}\left|s\right> \right|^2.$$
Now take expectation values of either side to get. Let $E_s = \left<s|H|s\right>$, then
\begin{align}
-2\frac{\hbar^2 q^2}{2m} 
&= 2E_s - \left<s \right| \left(   e^{-i q \cdot r}He^{iq \cdot r} + h.c.  \right) \left|s \right> \; , \\
 &= 2E_s -\left( \sum_n \left<s \right|  e^{-i q \cdot r}E_n \left|n \right> \left<n \right|e^{iq \cdot r} \left|s \right> + c.c.  \right) \; , \\
&=2E_s -\left( \sum_n E_n \left| \left<n \right|e^{iq \cdot r} \left|s \right> \right|^2 + c.c.  \right) \; , \\
&=2E_s -2 \sum_n E_n \left| \left<n \right|e^{iq \cdot r} \left|s \right> \right|^2 \; , \\
\end{align}
Using the identity and canceling the factor of $2$ we find
$$\frac{\hbar^2 q^2}{2m} = \sum_{n} \left( E_n - E_s\right) \left|\left<n \right|e^{i q \cdot r} \left| s \right> \right|^2 \; .$$
