Deriving the Pauli-Schrödinger equation from the Dirac equation Since the Schrödinger Pauli equation describes a non-relativistic spin ½ particle. This equation must be an approximation of the Dirac equation in an electromagnetic field. I was trying to derive this but I got stuck at a point. The free particle Dirac can be reduced to the equations
\begin{align}
σ^{i}(p_{i}+eA_{i}) u_B & = (E-m+eA_0)u_A. \\
\sigma^{i}(p_{i}+eA_{i})u_A & = (E+m+A_{0})u_B
\end{align}
I multiplied both sides of the first equation by $(E+m+eA_0)$ to get the Schrödinger Pauli equation. I was not able to eliminate $u_B$ completely from the equation. Can someone help me with the derivation?
 A: Dirca's equation has is the following form:
$$i\hbar\gamma^{\mu}\partial_{\mu}\psi - mc\psi = 0$$
where $\mu = 0,1,2,3$ and $\gamma^{\mu}$ are the $4\times4$ matrices(in Dirac's representation):
\begin{align}
\gamma^0 = 
\begin{bmatrix}
\mathbb{1} & 0 \\
0 & -\mathbb{1}
\end{bmatrix}
\gamma^k = 
\begin{bmatrix}
0 & \sigma_k \\
-\sigma_k & 0
\end{bmatrix}
\end{align}
Where $\sigma_k$ are the Pauli Matrices ($k = 1,2,3$):
\begin{align}
\sigma_1 =
\begin{bmatrix}
0 & 1 \\
1 & 0
\end{bmatrix}
\sigma_2 =
\begin{bmatrix}
0 & -i \\
i & 0
\end{bmatrix}
\sigma_3 =
\begin{bmatrix}
1 & 0 \\
0 & -1
\end{bmatrix}
\end{align}.
Using the momentum e energy in quantum mechanics differential operators; $\vec{p} = -i\hbar\vec{\nabla}$, $E = ih\partial_t = i\hbar c\partial_0$ the Dirac's equation become:
$$(E\gamma^0 -c\gamma^{k} p_k - mc^2)\psi = 0$$
Breaking the four dimensional spinnor $\psi$ in two components the equation in a matrix for become:
\begin{align}
\begin{bmatrix}
(E - mc^2)\mathbb{1}& -\sigma_k p_k c\\
\sigma_k p_k c & -(E + mc^2)\mathbb{1}
\end{bmatrix}
\begin{bmatrix}
u\\
v
\end{bmatrix}
=0
\end{align}
To coupling the EM field the energy and the momentum change 
\begin{align}
E \rightarrow E - e \Phi \\
\vec{p} \rightarrow \vec{p} - e \vec{A}
\end{align}
So the EM weak coupling Dirac's equation in matricial form become:
\begin{align}
\begin{bmatrix}
(E - e \Phi - mc^2)\mathbb{1}& -\sigma_k (p_k - e A_k) c\\
\sigma_k (p_k - e A_k) c & -(E - e \Phi + mc^2)\mathbb{1}
\end{bmatrix}
\begin{bmatrix}
u\\
v
\end{bmatrix}
=0
\end{align}
Now writing it in function of $u$ and $v$:
\begin{align}
(E - e \Phi - mc^2) u -(\sigma_k (p_k - e A_k) c)v =0 \\
\sigma_k (p_k - e A_k)c u -(E - e \Phi + mc^2)v =0
\end{align}
In the second equation the $u$ and $v$ relation become:
\begin{align}
\frac{\sigma_k (p_k - e A_k)c}{E - e \Phi + mc^2}u=v
\end{align}
For the non relativistic approach $e \Phi << mc^2$ and $E = mc^2$
So the $u$ and $v$ relation became:
\begin{align}
\frac{\sigma_k (p_k - e A_k)c}{2mc^2}u=v
\end{align}
Plugging this non relativistic relation in the first equation:
\begin{align}
(E - e \Phi - mc^2) u -\frac{(\sigma_k (p_k - e A_k)c)^2}{2mc^2}u =0 
\end{align}
Renaming $E - mc^2= E_{NR}$ and reorganizing:
\begin{align}
\left( e \Phi  +\frac{(\sigma_k (p_k - e A_k)c)^2}{2mc^2}\right)u =E_{NR}u 
\end{align}
Focus in the $(\sigma_k (p_k - e A_k))^2$ term we have:
\begin{align}
(\sigma_k (p_k - e A_k))^2 &= (\sigma_i (p_i - e A_i)(\sigma_j (p_j - e A_j)\\
&= \sigma_i \sigma_j \Pi_i \Pi_j
\end{align}
Where $\Pi_i = p_i - e A_i$
From Pauli anti-commutation and commutation relations:
$$\sigma_i \sigma_j = \delta_{ij} + i \varepsilon_{ijk}\sigma_k $$ 
With that the quadratic term become:
\begin{align}
(\sigma_k (p_k - e A_k))^2 &= \Pi_i \Pi_i + i \varepsilon_{ijk}\sigma_k \Pi_i \Pi_j
\end{align}
The last term in $(\sigma_k (p_k - e A_k))^2 u$ is:
\begin{align}
i \varepsilon_{ijk}\sigma_k \Pi_i \Pi_j u &= i \varepsilon_{ijk}\sigma_k\left[(-i\hbar\partial_i - e A_i)(-i\hbar \partial_j - e A_j)\right] \\
&= i \varepsilon_{ijk}\sigma_k\left[-\hbar^2 \partial_i \partial_j + e^2 A_i A_j + i\hbar e(\partial_i A_j + A_i \partial_j)\right]u
\end{align}
The first two terms are symmetrical so they become:
\begin{align}
\varepsilon_{ijk} \partial_i \partial_j &= \frac{1}{2}\varepsilon_{ijk}(\partial_i \partial_j + \partial_j \partial_i) \\
&= \frac{1}{2}(\varepsilon_{ijk}\partial_i \partial_j +\varepsilon_{ijk}\partial_j \partial_i) \\
&= \frac{1}{2}(\varepsilon_{ijk}\partial_i \partial_j -\varepsilon_{jik}\partial_j \partial_i)\\
& = 0
\end{align}
The same thing are value for the $A_i A_j$ so the term become:
\begin{align}
i \varepsilon_{ijk}\sigma_k \Pi_i \Pi_j u &= -\hbar e \varepsilon_{ijk}\sigma_k\left[ \partial_i( A_j u) + A_i \partial_j (u)\right] \\
&= -\hbar e \varepsilon_{ijk}\sigma_k\left[ \partial_i( A_j)u+ A_j \partial_i (u) + A_i \partial_j (u)\right]
\end{align}
Where were applied the product rule in the first derivative. Now looking in the last two terms we can show that they vanish each other:
\begin{align}
\varepsilon_{ijk}\left[ A_j \partial_i + A_i \partial_j \right] &= \varepsilon_{ijk}A_j \partial_i + \varepsilon_{ijk}A_i \partial_j \\
&= \varepsilon_{ijk}A_j \partial_i - \varepsilon_{jik}A_i \partial_j \\
&= 0
\end{align}
The only term that does not vanish is the $k$-th component os the curl of $\vec{A}$
\begin{align}
i \varepsilon_{ijk}\sigma_k \Pi_i \Pi_j &= -\hbar e \varepsilon_{ijk}\sigma_k  \partial_i( A_j) \\
&=  -\hbar e \sigma_k \varepsilon_{kij}\partial_i A_j \\
&= -\hbar e \vec{\sigma}\cdot \vec{B}
\end{align}
With that we can write:
\begin{align}
(\sigma_k (p_k - e A_k))^2 &= \frac{(p_k - e A_k)^2}{2m} - \frac{\hbar e}{2m}\vec{\sigma}\cdot \vec{B}
\end{align}
Now finally putting all together the non relativistic Dirac's equation is:
\begin{align}
\left( e \Phi + \frac{(p_k - e A_k)^2}{2m} - \frac{\hbar e}{2m}\vec{\sigma}\cdot \vec{B}\right)u =E_{NR}u 
\end{align}
A: Let's start with the Dirac equation as you've written
(but I choose to use lower indexes only, since they denote 3-vectors here, not 4-vectors):
$$\begin{align}
\sigma_i(p_i+eA_i)u_B & = (E-m+eA_0)u_A \\
\sigma_i(p_i+eA_i)u_A & = (E+m+eA_0)u_B
\end{align} \tag{1}$$
As @Sunyam suggested in his comment
you can solve the second equation of (1) for $u_B$
$$u_B=\frac{1}{E+m+eA_0}\sigma_j(p_j+eA_j)u_A$$
and then you plug that into the first equation of (1):
$$\frac{1}{E+m+eA_0}\sigma_i(p_i+eA_i)\sigma_j(p_j+eA_j)u_A = (E-m+eA_0)u_A \tag{2}$$
Since you want to derive the Pauli-Schrödinger equation
you restrict to the non-relativistic situation.
Therefore the total energy $E$ is positive and very close to the rest energy $m$,
in other words
$$E=m+E_{NR} \quad \text{with } E_{NR} \ll m \tag{3}$$
where $E_{NR}$ denotes the non-relativistic energy.
Inserting (3) into (2) you get
$$\frac{1}{2m+E_{NR}+eA_0}\sigma_i(p_i+eA_i)\sigma_j(p_j+eA_j)u_A = (E_{NR}+eA_0)u_A \tag{4}$$
This equation is still exact.
But in the non-relativistic situation, the non-relativistic energy $E_{NR}$
and the electric energy $eA_0$ both are much smaller than the rest energy $m$.
Therefore you can neglect these in the denominator on the left-hand-side of (4) and get
$$\frac{1}{2m}\sigma_i(p_i+eA_i)\sigma_j(p_j+eA_j)u_A \approx (E_{NR}+eA_0)u_A. \tag{5}$$
Now you're nearly finished.
Using $\sigma_i\sigma_j=\delta_{ij}+i\epsilon_{ijk}\sigma_k$,
$\vec{p}=-i\vec{\nabla}$, and doing some simple but tedious math
you can rewrite equation (5) to
$$\frac{1}{2m}\left((p_i+eA_i)(p_i+eA_i)+e\sigma_i\epsilon_{ijk}\nabla_j A_k\right)u_A \approx (E_{NR}+eA_0)u_A. \tag{6}$$
Here you recognize $\epsilon_{ijk}\nabla_j A_k$ as the magnetic field $B_i$
and have the Pauli-Schrödinger equation
$$\frac{1}{2m}\left((p_i+eA_i)(p_i+eA_i)+e\sigma_i B_i\right)u_A \approx (E_{NR}+eA_0)u_A. \tag{7}$$
