$\mathcal{I}\psi_{n,k}(r)=\psi_{n,-k}(-r)$ is this the true behaviour of Bloch waves under parity operators, or
is this one the correct behaviour
$\mathcal{I}\psi_{n,k}(r)=\psi_{n,k}(-r)=\psi_{n,-k}(r)$?
Parity defined as,
$ \mathcal{I}p\mathcal{I}^{-1}=-p$ and $\mathcal{I}x\mathcal{I}^{-1}=-x$, so I think under parity operator both $k$ and $r$ should change sign.
Best Answer
I think the answer is that you need the lattice to have inversion symmetry to begin with to get a simple answer, otherwise in principle you might get a sum of bloch waves of different momenta.
Following some Cornell notes here by F. Rana which I'll summarise in case the link dies, and add some extra bits.
So supposing you have a lattice which has symmetry $S$ then this reflects itself in the potential having the same symmetry, $V(Sr)=V(r)$. In the case of inversion we have $V(-r)=V(r)$.
Now the bloch functions solve the schrodinger equation,
$$\left(-\frac{\hbar}{2m}\nabla^2+V(r)\right)\psi_{n,k}(r)=E_n(k)\psi_{n,k}(r)$$
If we replace $r$ everywhere by $r'=S(r)=-r$ then we have that \begin{align} \left(-\frac{\hbar}{2m}\nabla^2_{r'}+V(r')\right)\psi_{n,k}(r')=E_n(k)\psi_{n,k}(r')\\ \left(-\frac{\hbar}{2m}\nabla^2_{r}+V(r)\right)\psi_{n,k}(r')=E_n(k)\psi_{n,k}(r')\\ \end{align}
Which shows that $\psi_{n,k}(r')$ has the same energy eigenvalue as $\psi_{n,k}(r)$. To check what happens with the momentum, we remember that
$$\psi_{n,k}(r+R)=e^{ikR}\psi_{n,k}(r) \quad \text{for } R \text{ a lattice vector} \quad*$$
So translating $\psi_{n,k}(r')$ we have that
$$T_R\psi_{n,k}(-r)=\psi_{n,k}(-(r+R))=\psi_{n,k}(-r-R)$$
but due to the symmetry $-R$ is also a lattice vector so we should have from $*$
$$\psi_{n,k}(-r-R)=e^{ik(-R)}\psi_{n,k}(-r)=e^{i(-k)R}\psi_{n,k}(-r)$$
which is an eigenstate of the translation operator $T_R$ with eigenvalue $e^{i(-k)R}$. So $S\psi_{n,k}$ has the same energy as $\psi_{n,k}$ and an opposite momentum value, which is the same as $\psi_{n,-k}$. The conclusion is then that
$$S\psi_{n,k}(r)=\psi_{n,k}(-r)=\psi_{n,-k}(r)$$
This also suggests a symmetry in the energy $$E(-k)=E(k)$$
The other to see that it might be one or the other is to write
$$\langle r | n,k\rangle=\psi_{n,k}(r)$$
On the one hand, using your $\mathcal{I}$ for inversion, we have that $$\psi'(r)=\langle r | \mathcal{I} |n,k\rangle=\langle -r |n,k\rangle=\psi_{n,k}(-r)$$
by acting $\mathcal{I}$ on the left and using $\mathcal{I}=\mathcal{I}^\dagger=\mathcal{I}^{-1}$. On the other hand we can check the crystal momentum of the state after inversion by acting with a translation operator $T_R$ by a lattice vector $R$, writing $T_R=e^{-i\hat{P}R/\hbar}$. We have that $T_R|r\rangle=|r+R\rangle$ or $$\psi_{n,k}(r+R)=\langle r | T_R^\dagger |n,k\rangle=e^{ikR}\psi_{n,k}(r)$$
So we know that $T^\dagger_R|n,k\rangle=e^{ikr}|n,k\rangle$. If we also use that $$T_R\mathcal{I}=\mathcal{I}T_{-R}$$we have that \begin{align} \psi'(r+R)=\langle r | T_R^\dagger\mathcal{I} |n,k\rangle&=\langle r | \mathcal{I}T_{-R}^\dagger |n,k\rangle\\ &=\langle r | \mathcal{I}e^{-ikR} |n,k\rangle\\ &=e^{-ikR}\langle r | \mathcal{I} |n,k\rangle\\ &=e^{i(-k)R}\psi'(r) \end{align}
So we should have a bloch function of crystal momentum $-k$, and we could call $\psi'$ instead $\psi_{n,-k}(r)$. Then we again assuming the same form of the Hamiltonian and symmetry of the potential we should have that $[H,\mathcal{I}]=0$, to get that the state has the same energy.
Another note would be that if you have a (orthogonal) symmetry $S$ of the crystal lattice, then writing the wave function in the other form of bloch's theorem,
$$\psi_{n,\vec{k}}(\vec{r})=e^{i\vec{k}\cdot\vec{r}}u_{n,\vec{k}}(\vec{r})$$
We can construct $\psi_{n,S(\vec{k})}$ from $\psi_{n,k}(S^{-1}\vec{r})=\psi_{n,k}(S^{T}\vec{r})$.
\begin{align} \psi_{n,\vec{k}}(S^{-1}\vec{r})&=e^{i\vec{k}\cdot(S^{-1}\vec{r})}u_{n,\vec{k}}(S^{-1}\vec{r})\\ &=e^{i\vec{k}\cdot(S^{T}\vec{r})}u_{n,\vec{k}}(S^{-1}\vec{r})\\ &=e^{i(S\vec{k})\cdot\vec{r}}u_{n,\vec{k}}(S^{-1}\vec{r})\\ \end{align} or if we call $u_{n,S(\vec{k})}(r)\equiv u_{n,\vec{k}}(S^{-1}\vec{r})$ $$\psi_{n,S(\vec{k})}(\vec{r})=e^{i(S\vec{k})\cdot\vec{r}}u_{n,S(\vec{k})}(r)$$