It is a standard procedure that you can find in any good quantum mechanical book. I suggest Messiah, since it's the best I know ! Everyone should check the Messiah before asking on this web-site ;-)

I call $\Phi$ the 4-spinor, it is defined as

$$\Phi\equiv\left(\begin{array}{c} u\\ v \end{array}\right)\;;\;\Phi^{\dagger}\equiv\left(\begin{array}{cc} u^{\dagger} & v^{\dagger}\end{array}\right)\;;\; u\equiv\left(\begin{array}{c} u_{\uparrow}\\ u_{\downarrow} \end{array}\right)$$

in my favorite representation.

The Schrödinger equation then reads

$$\mathbb{H}\Phi=\mathbf{i}\hbar\dfrac{\partial\Phi}{\partial t}$$

with $\mathbb{H}$ identified with your matrix. When the problem is time-independent, the usual substitution $\Phi=e^{-\mathbf{i}\varepsilon t/\hbar}\Psi$ gives back your BdG equation. Then you calculate

$$-\mathbf{i}\hbar\dfrac{\partial\Phi^{\dagger}}{\partial t}=\Phi^{\dagger}\mathbb{H}^{\dagger}$$

for the adjoint spinor. You can verify that $\mathbb{H}^{\dagger}$ is the Hermitian conjugate of the $\mathbb{H}$ indeed (transpose of the complex conjugate).

Then you define the probability of particle as $\rho_{p}=\left|\Phi\right|=\Phi^{\dagger}\Phi$ and you calculate

$$\dfrac{\partial\rho_{p}}{\partial t}=\dfrac{\partial\Phi^{\dagger}}{\partial t}\Phi+\Phi^{\dagger}\dfrac{\partial\Phi}{\partial t}$$

and you inject the previous relations. By virtue of the conservation of the number of particle, you must obtain something like $\partial_{t}\rho_{p}+\nabla\cdot\mathbf{j}_{p}=0$, so you next have to identify the current-probability $\mathbf{j}_{p}$ as some divergence in your expression for $\partial_{t}\rho_{p}$. For a s-wave superconductor and when $H_{0}=p^{2}/2m=-\hbar^{2}\nabla^{2}/2m$, I found

$$\mathbf{j}_{p}=\dfrac{\hbar}{m}\left[\Im\left\{ u^{\dagger}\nabla u-v^{\dagger}\nabla v\right\} \right]$$

which sounds correct. This result does not depend on the spin-interaction. Nevertheless, you may choose a $H_{0}=v_{0}\left(\boldsymbol{\sigma\cdot p}\right)$ and you will obtain something different of course (in that case the identification of the divergence of the current is straightforward, so I let you finding the details in the Messiah ... section about the Dirac equation ... )

The trick is that charge are not conserved. Indeed, if you reproduce the same calculation for the density of charge $\rho_{e}=\Phi^{\dagger}\tau_{3}\Phi=u^{\dagger}u-v^{\dagger}v$ you should find something like

$$\dfrac{\partial\rho_{e}}{\partial t}+\nabla\cdot\mathbf{j}_{e}=\dfrac{2e}{\hbar}\left[\Delta^{\ast}\left(v^{\dagger}\sigma_{2}u\right)+\Delta\left(u^{\dagger}\sigma_{2}v\right)\right]$$

with $\mathbf{j}_{e}=\dfrac{e\hbar}{m}\Im\left\{ u\nabla u^{\dagger}+v\nabla v^{\dagger}\right\} $, and the result does not depend on the spin interaction neither. The first paper to discuss this problem was the seminal paper by Tinkham, Blonder and Klapwick (I can send it to you if you have some difficulties to obtain it, follow the links in my profiles and you may find my email address out, but first check around you: library, colleagues, ...)

Only for the spin current might you find that a non-diagonal $H_{0}$ will induce fancy effects. The spin-current might be defined from the spin-density $\rho_{s}^{j}=\Phi^{\dagger}\sigma_{j}\Phi=u^{\dagger}\sigma_{j}u+v^{\dagger}\sigma_{j}v$ and you calculate the time-derivative as above. The current is not conserved one more time. If you really insist I may find the time to write out a generic expression for the spin-current for a generic Hamiltonian.

Please ask for more details if the above explanations are not sufficient.

It is a standard procedure that you can find in any good quantum mechanical book. I suggest Messiah, since it's the best I know ! Everyone should check the Messiah before asking on this web-site ;-)

I call $\Phi$ the 4-spinor, it is defined as

$$\Phi\equiv\left(\begin{array}{c} u\\ v \end{array}\right)\;;\;\Phi^{\dagger}\equiv\left(\begin{array}{cc} u^{\dagger} & v^{\dagger}\end{array}\right)\;;\; u\equiv\left(\begin{array}{c} u_{\uparrow}\\ u_{\downarrow} \end{array}\right)$$

in my favorite representation.

The Schrödinger equation then reads

$$\mathbb{H}\Phi=\mathbf{i}\hbar\dfrac{\partial\Phi}{\partial t}$$

with $\mathbb{H}$ identified with your matrix. When the problem is time-independent, the usual substitution $\Phi=e^{-\mathbf{i}\varepsilon t/\hbar}\Psi$ gives back your BdG equation. Then you calculate

$$-\mathbf{i}\hbar\dfrac{\partial\Phi^{\dagger}}{\partial t}=\Phi^{\dagger}\mathbb{H}^{\dagger}$$

for the adjoint spinor. You can verify that $\mathbb{H}^{\dagger}$ is the Hermitian conjugate of the $\mathbb{H}$ indeed (transpose of the complex conjugate).

Then you define the probability of particle as $\rho_{p}=\left|\Phi\right|=\Phi^{\dagger}\Phi$ and you calculate

$$\dfrac{\partial\rho_{p}}{\partial t}=\dfrac{\partial\Phi^{\dagger}}{\partial t}\Phi+\Phi^{\dagger}\dfrac{\partial\Phi}{\partial t}$$

and you inject the previous relations. By virtue of the conservation of the number of particle, you must obtain something like $\partial_{t}\rho_{p}+\nabla\cdot\mathbf{j}_{p}=0$, so you next have to identify the current-probability $\mathbf{j}_{p}$ as some divergence in your expression for $\partial_{t}\rho_{p}$. For a s-wave superconductor and when $H_{0}=p^{2}/2m=-\hbar^{2}\nabla^{2}/2m$, I found

$$\mathbf{j}_{p}=\dfrac{\hbar}{m}\left[\Im\left\{ u^{\dagger}\nabla u-v^{\dagger}\nabla v\right\} \right]$$

which sounds correct. This result does not depend on the spin-interaction. Nevertheless, you may choose a $H_{0}=v_{0}\left(\boldsymbol{\sigma\cdot p}\right)$ and you will obtain something different of course (in that case the identification of the divergence of the current is straightforward, so I let you finding the details in the Messiah ... section about the Dirac equation ... )

The trick is that charge are not conserved. Indeed, if you reproduce the same calculation for the density of charge $\rho_{e}=\Phi^{\dagger}\tau_{3}\Phi=u^{\dagger}u-v^{\dagger}v$ you should find something like

$$\dfrac{\partial\rho_{e}}{\partial t}+\nabla\cdot\mathbf{j}_{e}=\dfrac{2e}{\hbar}\left[\Delta^{\ast}\left(v^{\dagger}\sigma_{2}u\right)+\Delta\left(u^{\dagger}\sigma_{2}v\right)\right]$$

with $\mathbf{j}_{e}=\dfrac{e\hbar}{m}\Im\left\{ u\nabla u^{\dagger}+v\nabla v^{\dagger}\right\} $, and the result does not depend on the spin interaction neither. The first paper to discuss this problem was the seminal paper by Tinkham, Blonder and Klapwick (I can send it to you if you have some difficulties to obtain it, follow the links in my profiles and you may find my email address out, but first check around you: library, colleagues, ...)

Only for the spin current might you find that a non-diagonal $H_{0}$ will induce fancy effects. The spin-current might be defined from the spin-density $\rho_{s}^{j}=\Phi^{\dagger}\sigma_{j}\Phi=u^{\dagger}\sigma_{j}u+v^{\dagger}\sigma_{j}v$ and you calculate the time-derivative as above. The current is not conserved one more time. If you really insist I may find the time to write out a generic expression for the spin-current for a generic Hamiltonian.

Please ask for more details if the above explanations are not sufficient.

Important Edit: Contrary to what is said above, the charge current is obviously conserved in a superconductor, provided the gap parameter is treated self-consistently. See

• Scattering Theory of Transport for Mesoscopic Superconductors, by S. Datta, P. F. Bagwell, and M. P. Anantram (1996).

• L.-F. Chang, S. Chaudhuri, and P. Bagwell, Phys. Rev. B 54, 9399 (1996).

for instance.