The solutions (\[sol\_phi\]) are chosen to be the eigenfunctions of the $ z $-component of the fermion particle angular momentum operator $ J_{z} =
-i \frac{{\partial}}{{\partial}\varphi} + \frac{1}{2} \Sigma_{3} $, which corresponds to the axial symmetry existing in our problem: $$J_{z} \psi = \left( l - {{\textstyle}\frac{1}{2}} \right) \psi, \qquad
    l = n - s.$$

In the standard (or Dirac) representation of the $ \gamma $-matrices the coefficients $ \{
c_{\mathfrak a} \} $, which meet the normalization requirement for the wave functions $$\int r dr \, d\varphi \, \phi^{\dag} \phi^{ } = 1, $$ can be written as follows: $$\label{sol_coeffs}
    \setlength\arraycolsep{1pt}
    \begin{pmatrix}
        c_{1} \\
        c_{2} \\
        c_{3} \\
        c_{4} \\
    \end{pmatrix}
    = \frac{1}{2 \sqrt{2}}
    \begin{pmatrix}
        A ( P \alpha + \epsilon \zeta Q \beta )        \\
        -\zeta B ( P \alpha - \epsilon \zeta Q \beta ) \\
        A ( P \beta - \epsilon \zeta Q \alpha )        \\
        \zeta B ( P \beta + \epsilon \zeta Q \alpha )  \\
    \end{pmatrix},$$ where $$\label{sol_ABPQ}
    \begin{aligned}
        A &= \sqrt{1 + \frac{\tilde{m}}{\Pi}}, & \qquad &
        P &= \sqrt{1 + \frac{\tilde{p}}{E}}, \\
        B &= \sqrt{1 - \frac{\tilde{m}}{\Pi}}, & \qquad &
        Q &= \sqrt{1 - \frac{\tilde{p}}{E}}, \\
    \end{aligned}$$ and $$\label{sol_alpha_beta}
    \alpha = \cos{{\textstyle}\frac{\delta}{2}} 