{'fraction_non_alphanumeric': 0.050472270592580866, 'fraction_numerical': 0.016595430320299722, 'mean_word_length': 4.277499303815094, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
General solution vs spin invariant eigenstates of the Dirac equation with the Coulomb potential

Solutions of the Dirac equation for an electron in the Coulomb potential are obtained using operator invariants of the equation, namely the Dirac, Johnson-Lippmann and recently found new invariant. It is demonstrated that these operators are the spin invariants. The generalized invariant is constructed and the exact general solution of the Dirac equation are found. In particular, the explicit expressions of the bispinors corresponding to the three complete sets of the invariants, their eigenvalues and quantum numbers are calculated. It is shown that the general solution of one center Coulomb Dirac equation contains free parameters. Changing one or more of these parameters, one can transform one solution of the Dirac equation into any other. It is shown for the first time that these invariants determine electron spatial probability amplitude and spin polarization in each quantum state. Electron probability densities and spin polarizations are explicitly calculated in the general form for several electron states in the hydrogen-like energy spectrum. Spatial distributions of these characteristics are shown to depend essentially on the invariant set, demonstrating, in spite of the accidental degeneracy of energy levels, physical difference of the states corresponding to different spin invariants.


Introduction
The well-known and comprehensively studied Dirac equation (DE) is the main equation of the relativistic theory of particles in external electromagnetic field. In Eq. (1) ϕ(r) and A(r) are scalar and vector potentials of the field. Basic role of the DE was revealed soon after its discovery in 1928 [1] (see [2,3,4]). In the same year C. G. Darwin wrote "... Dirac has brilliantly removed the defects before existing in the mechanics of the electron, and has shown how the phenomena usually called the " spinning electron " fit into place in the complete theory" [5]. The study of particle motion in the Coulomb potential within the DE [5,6] resulted in prediction of the fine-structure splitting of hydrogen-like spectrum due to spin-orbit interaction. Solutions of the DE were investigated in numerous papers [5,6,7,8,9,10,11] whose results are summarized in textbooks (see, e.g., [2,3,4,12,13,14,15]).
The DE, as was demonstrated by Darwin [5] and Gordon [6], admits exact stationary solutions, which means that Dirac problem with the Coulomb potential describes the integrable system. It can be characterized by the set of independent commuting operators H = {Ĥ 1 ,Ĥ 2 , . . . ,Ĥ n }, where n corresponds to the number of the degrees of freedom, and where the set H includes the Hamiltonian and operators of the integrals of motion, i.e., invariants. In quantum case the operators of observables in such sets have a joint system of eigenstates, and their eigenvalues determine quantum numbers of each state. It is well known that in the relativistic Coulomb problem the Hamiltonian and the angular momentum do not generate a complete set of commuting observables, and an additional operator is needed to complete the set. In particular, the integrability of the non-relativistic Schrödinger equation is secured by a set of three independent commuting operators including the Hamiltonian, so that the quantum states acquire three quantum numbers. The state space of a particle with spin 1/2 is represented by the tensor product of the spatial, including coordinate and momentum components, and the two-dimensional spin states. In other words, integrability of the DE requires four such operators in the full set. This gives an additional quantum number which takes two values. Its meaning in the non-relativistic limit is the spin quantum number σ.
For the central-symmetric field integrals of motion of the DE include the total angular momentumĴ and operatorK D used by Dirac (Dirac invariant, DI). The angular momentum is given by the sum of the orbitalLÎ and intrinsicΣ angular momenta,Ĵ = LÎ + ( /2)Σ, andÎ is 4 × 4 unit matrix. This means that the set H D = {Ĥ D ,Ĵ 2 ,Ĵ z ,K D } withĤ D being Dirac Hamiltonian, provides complete integrability of the DE. This permits to reduce the system of four partial differential equations for spinor amplitudes to the system of the first order ordinary differential equations for radial and angular functions [16].
There is also the so called Johnson-Lippmann invariant (JLI),Â JL , constructed by Johnson and Lippmann (JLI) in 1950 [17]. In classical and quantum physics it is established that presence of integrals of motion additional to the number of the degrees of freedom, makes system degenerate. For instance, the accidental degeneracy takes place for a particle in the central field V ∼ 1/r.
The degeneracy is substantially lowered for the relativistic Kepler problem. Nevertheless, the removal of the dynamic symmetry is not complete and most of the energy levels remain to be degenerate. The degeneracy is accounted for by means of the Johnson-Lippmann constant of motion [17] whose operator does not commute with the DI and results in the accidental degeneracy. It allows to present equations for radial functions in various representations [8,9,16], as well as to obtain a purely algebraic solution of the DE using the methods of supersymmetric quantum mechanics (see [18,19,20,21,22]).
The presence of two noncommutative invariants, resulting in energy degeneracy, allows to describe stationary states by different systems of Hamiltonian eigenfunctions corresponding to different invariants and, therefore, with different physical meaning. For example, in Ref. [23] the observation of the nodal structure of Stark hydrogen levels has been reported. Theoretical explanation of the obtained data was based on the solution of the Schrödinger equation in parabolic coordinates, i.e. on the joint system of eigenfunctions of the set H par = {Ĥ,Â z ,L z } whereÂ z is the operator of z-component of Laplace-Runge-Lenz vector. Direct comparison of the experimental and calculated quantities has revealed good agreement [23] and validated that hydrogen discrete states correspond to the set H par , not to the set H sph = {Ĥ,L 2 ,L z } which gives solution in spherical coordinates. In view of this one can hope that in future precise experiments, depending on the chosen experimental conditions, hydrogen atom quantum states corresponding to each set of invariants, can be measured.
In the present paper we aim to calculate the solutions of the one center Coulomb Dirac equation using the complete set of invariants, which includes the Hamiltonian, square of the total angular momentum, its projection and one of the spin invariants, namely Dirac, Johnson-Lippmann or recently found new invariant [24]. We aim to prove that spatial distributions of the electron charge probability and spin expectations depend essentially on the invariant set, demonstrating, in spite of the accidental degeneracy of energy levels, physical difference of the states corresponding to different spin invariants.
The paper, which is the continuation of our previous paper [24], is organized as follows. First, in Section 2, we discuss solution of the stationary DE with the Coulomb potential in a way different from the standard textbooks and Ref. [24] in order to demonstrate that the DE eigen bispinors are not completely determined without the four invariants. Then, in Section 3, we analyze recently found new invariant of the DE [24], called BEL invariant (BELI), and study the eigenvalue problem for the sets of the independent commuting operators {Ĥ,Ĵ 2 ,Ĵ z ,Î inv } whereÎ inv is DI, JLI or BELI. Explicit expressions of the eigen bispinors for each set are calculated and a full set of quantum numbers of the stationary states is defined. In Section 4 we construct the generalized invariant and calculate its eigen bispinors. In Section 5 the main attention is devoted to electron probability densities and spin polarizations for some hydrogen-like states. Finally, in Section 6, the summary of the obtained results and their discussion are presented.

Invariants of the Dirac Equation and its solution
In time-independent Coulomb field of a positive charge Ze electron states Ψ(r, t) = Ψ(r) exp (−iEt/ ) are determined by the DÊ in which E is the energy, operatorĤ is the Dirac Hamiltonian, is a 4-component amplitude of the spinor field (4-spinor, or bispinor), c is speed of light , p = −i ∇ is momentum operator,β and componentsα j (j = x, y, z) of the vector-matrix α = j e jαj together with the unit matrixÎ are 4 × 4 Hermitian Dirac matrices (DM). In fact, Eq. (2) is a system of four equations for the components ψ ν (ν = 1, 2, 3, 4), which acquire the obvious form after choosing DM. Below we use their standard form, represent 4 × 4-dimensional matrices via 2 × 2 matrices and write the Hamiltonian in the formĤ whereÎ 2 is a unit 2 × 2 matrix,σ j (j = x, y, z) are Pauli matrices,σ = j e jσj , and e j are the corresponding orths.
The block form of DM allows to represent the bispinor as where ψ (u/d) are its upper/lower spinors, respectively, with the components ψ (ν) and ψ (ν+1) (ν = 1 for the upper bispinor, and ν = 3 for the lower one). In representation (3) and (4) the DE (2) is transformed to To obtain these equations in polar coordinates, one can use transformation of the kinetic energy operatorσp using a unit matrix and the identitŷ This transformation is equivalent to transformations of DE using DM algebra [12,8].
Functions f ν (ϑ) satisfy differential equations which can be reduced to one equation of the second order: So, spinor components are described by equation known in theory of spherical functions. It admits solutions that satisfy finiteness and unambiguity conditions, provided equality is valid, where l = 0, 1, 2, . . . and m ν takes integer numbers in the interval −l ≤ m ν ≤ l. The solution is given by the associated Legendre polynomial P mν l (cos ϑ). Relation (12) gives two eigenvalues λ, one of which is positive, and another one -negative due to l positiveness. Therefore, spinor components, with the accuracy of a constant multiplier, are Substituting these expressions into Eqs. (10) and requiring their identical validity, one gets the relation for the constants C 1 = (λ + m 1 ) C 2 which depends on the eigenvalue λ ± and gives two solutions for spinors corresponding to the positive and negative eigenvalue, respectively. Constants C 2± are defined by the normalization condition. Note, normalization for spinors determines the incoming constants up to the phase multiplier. Therefore, the spinors are spherical ones, Their components are characterized by integer numbers m ν which are connected via relations (11). Two numbers with the given difference can be represented by m 2 + m 1 = 2M whence M which takes half-integer values, distinguishes spherical spinors and is included into their subscripts (15) Using its exact forms (15)- (16), one can define coefficient in the relationσ r χ λ ∼ χ −λ and come to the following equalities [4] σ r χ l,M,+ = iχ l+1,M,− ,σ r χ l,M,− = −iχ l−1,M,+ .
Spherical harmonics determine angular dependence of spinors ψ (u/d) . Solution of Eqs. (9) can be represented as where the coefficients at χ lν ,M,± in the upper F (±) and lower G (±) spinors can depend on r, only, and their upper indeces (±) indicate positive or negative eigen value ofΛ of the corresponding spherical spinor. Substitution of (19) into Eqs. (9) with account of (17) and equalities (18) demonstrates that solution requires fulfilment of the relationŝ So, positive integer numbers l 1 and l 2 in (19) should be connected by the equality and Eq. (9) gives the following equalities: Since the spinors χ l,M,± are independent, the relations (21) are valid at zero values of the "coefficients", which gives four equations for functions F (±) (r) and G (±) (r). The pair of numbers l 1 and l 2 has the given difference (20), and, thus, one can represent them as l 1 = j − 1/2 and l 2 = j + 1/2, where j ≡ (l 1 + l 2 )/2 takes positive half-integer values. So, the positive numbers l 1 + 1 and l 2 are l 1 + 1 = j + 1/2 and l 2 = j + 1/2, and, hence, Eqs. (21) for radial functions involve the number which takes non-zero positive integer values.
System of equations for radial functions separates into two independent pairs, onefor functions F (+) (r) and G (−) (r), and the second one -for functions F (−) (r) and G (+) (r). Similar equations for such functions are well known and their solutions are expounded in many textbooks (e.g., [3,4,12]).
It is convenient to introduce dimensionless coordinate ξ = (mc/ ) r for functions rF (±) and rG (∓) , where /mc is electron Compton wavelength. Taking into account the boundary condition at infinity, solutions for the bound states can be searched for in the form where are dimensionless damping κ and energy ε, ε = E/mc 2 . Substitution of expressions (23) into equations for radial functions following from Eqs. (21), leads to equations one of which gives functions u 1 and v 1 , and the second oneu 2 and v 2 : where the upper sign corresponds to functions with ν = 1, lower -with ν = 2, and α = e 2 / c is Sommerfeld fine structure constant. This system can be reduced to the second order equations for each function: from where it follows that functions u and v should satisfy hypergeometric differential equation. Finitness of the solution requires that hypergeometric series should terminate at some value of n which imposes two conditions: γ 2 − κ 2 j + Z 2 α 2 = 0 and γ − Zαε/κ = n, where n is zero or positive integer. Thus, equations for u and v become the differential equation for generalized Laguerre polynomial L 2γ n (ρ) for variable 2κξ = ρ. From the first condition we get where the square root has positive sign only, as it follows from the condition of the convergence of the normalization. The second one is a well-known equation for the relativistic hydrogen-like energy spectrum. So the searched functions can be expressed via generalized Laguerre polynomials: where the radial quantum number n r determines the degree of the Laguerre polynomials.
Equation (27) gives hydrogen-like energy spectrum for bound states where κ nr,j = Zα/N nr,j and the notation is used. From the functional relations for Laguerre polynomials the relations follow: Thus, the solution of DE (21) is given by spinors (19) with radial functions where Here the notation for dimensionless radial variable is used, in which r B = 2 /me 2 is Bohr radius. Constants A and B (31) can be found from the general normalization Introducing the common normalization constant 2Z/r B N nr,j C nr,j , where one can represent A and B as where β 1 and β 2 comply the relation Therefore, the solution of the DE (2) is given by the bispinor (4), Ψ = Ψ nr,j,M (r) with its upper and lower spinors defined by expressions (19) with spherical harmonics (15)- (16) and functions (30), whose coefficients are (34). It is described by numbers M , j and n r : Here spinors χ j∓1/2,M,± are given in Eqs.(15)- (16) and radial functions are with the constant defined in Eq. (33). Bispinor (36) provides the general solution for the bound states of the DE (28). This solution is however not completely determined due to arbitrary parameters β 1 and β 2 with condition (35), which leaves the ambiguity in its choice.

Eigenvalues of the invariants
As it was mentioned above, stationary states are joint eigenstate vectors of a complete set of independent commuting operators which includes the Hamiltonian and invariant operators. For one center Coulomb DE the invariants are given by the square of the total angular momentum,Ĵ =LÎ + 2Σ , its z-component, and DI,K D , or JLI,Â JL . The set of quantum numbers which characterize stationary states, reflects the eigenvalues of all commuting operators in the set {Ĥ,Ĵ 2 ,Ĵ z ,Î inv }, whereÎ inv is DI, JLI or BELI. The form of the latter is introduced below. So, the set of quantum numbers is determined by the following eigenvalue equationŝ Here notation ε = E/mc 2 is used, 2 j (j + 1) and m j are the eigenvalues of the opera-torsĴ 2 andĴ z respectively, defined algebraically from the commutation relations for the components of the operatorĴ, and inv is the eigenvalue of the operatorÎ inv . The first equation in the set (38) determines the energy spectrum (28) and gives the radial quantum number n r . Numbers j and M , whose sense is clear and generally known, in our calculations are not connected with the invariants. They appear in the solution due to two pairs of auxiliary numbers m 1 , m 2 and l 1 , l 2 with the given differences (11) and (20).
Taking into account the block-diagonal form of the matrices of the operatorsĴ 2 and J zĴ with equal blocks, we see from the second and third equations in (38) for the bispinor (36) that the upper ψ (u) and lower ψ (d) spinors in the eigenvalue problemsĴ z Ψ = M Ψ andĴ 2 Ψ = 2 j(j + 1)Ψ satisfy the same equations. The identityσ ·Lσ ·L =L 2 − σ ·L is valid, hence, the Legendre operator can be represented in the spinor form [9] and we getL 2 + σL + 3 Angular dependence of spinors is included in spinors χ j∓1/2,M,± , only, as it follows from Eqs.
where κ j is defined in (22). The physical meaning of numbers j and M ≡ m j follows from the equationŝ J 2 Ψ nr,j,M = 2 j (j + 1) Ψ nr,j,M ,Ĵ z Ψ nr,j,M = M Ψ nr,j,M which define them as the definite values of the total angular momentum and its projection on the polar axis in a given stationary state. The latter, of course, is a common knowledge, but here it was obtained in an original way.
(i) In relativistic problem, in addition to operators (39) and (40), there is DI [2] It commutes with the invariantsĴ 2 andĴ z , has a block-diagonal form and acts on angular variables only. Therefore, the fourth equation from the set (38),K D Ψ nr,j,M = D Ψ nr,j,M , can be written down as two simple spinor equationŝ  (β 1 = 0, β 2 = 1), the eigenvalue D = − κ j is negative. Thus, the last eigenvalue equation in Eqs. (38) fixes free parameters in the general expression (36) and results in two orthogonal normalized eigenbispinors with completely determined components. One can assign to the bispinors the number σ = ± which indicates the sign of the eigenvalueK D , and write them down in the form Bispinors (43) are the Darwin-Gordon solutions [5,6] and represent the joint eigenbispinor system of the set {Ĥ D ,Ĵ z ,Ĵ 2 ,K D } for which two independent but fixed pairs of the numbers β The assigned bispinor number σ = ± is the fundamental characteristic of stationary states and is the fourth quantum number. It was also presented in a recent paper [16] within symmetry algebra, in particular, using Yano-Killing operators.
(ii) In the Coulomb potential the JLI [17] in conventional representation of Dirac matrices has the form are the Laplace-Runge-Lenz vector operators. It is easy to see that JLI commutes witĥ J z andĴ 2 and does not commute withK D . Therefore, the eigenstates of the Dirac Hamiltonian can be described by eigenbispinors of the set H JL = {Ĥ,Ĵ z ,Ĵ 2 ,Â JL } which differ from eigenbispinors (43). In this case the fourth eigenvalue equation in Eqs. (38) isÂ JL Ψ nr,j,m j = JL Ψ nr,j,m j which can be written down as the system of two spinor equations. By calculating vector operators (46) explicitly and writing them in spherical coordinates, one can obtain that where operatorp r and matricesΛ andσ r were defined in (7), (8), and (6), respectively. This leads to the equations nr,j,m j .
Substituting spinors (19) and taking into account functional relations for Laguerre polynomials, we conclude that the bispinor (36) becomes eigen one for JLI if free constants A and B in functions (30) satisfy the relations The condition for their non-trivial solution gives two eigenvalues of the JLI JL = ±a nr,j , a nr,j = Zme 2 n r (n r + 2γ j ) and, respectively, two solutions for constants A and B, and, threrefore (see (34)), two solutions for coefficients (β 1 , β 2 ) first of which corresponds to the states with positive eigenvalue of the JLI, JL = a nr,j , and the second one -with negative, JL = −a nr,j . As a result, for the set H JL = {Ĥ,Ĵ z ,Ĵ 2 ,Â JL } we obtain the pair of orthonormal bispinors with the numbers σ = ± assigned to them according to the sign in the eigenvalue equationÂ JL Ψ nr,j,m j ,σ = σa nr,j Ψ nr,j,m j ,σ , where σ = ± characterizes the stationary states. These states are different from the states for the DE with DI, although the spectra of both sets coincide.
(iii) In addition to the invariants (42) and (45), we constructed in [24] the operator of a new invariant, BELI,Î which does commute withĤ, but does not commute with either operator it is constructed of. Hence, the set H BEL = {Ĥ,Ĵ z ,Ĵ 2 ,Î BEL } defines one more system of eigenbispinors, different from the previous two cases. The explicit expression for BELI can be obtained by direct calculation, which giveŝ whereΛ has been defined in Eq. (8), and the notation for the operator is used (cp. Eq. (46)). In polar coordinates the scalar productσ ·B ± in the invariant (50) has the formσ The last eigenvalue equation in Eqs. (38),Î BEL Ψ nr,j,m j = BEL Ψ nr,j,m j , leads to spinor equations nr,j,m j , nr,j,m j .
Substituting (19) and using explicit expressions for radial functions (30), we see that the bispinor (36) is the eigen bispinor of the BELI provided free constants A and B satisfy the relations (n r + 2γ j ) A + N nr ,j Respectively, the condition for non-trivial values of these constants determines two eigenvalues BEL = ± κ j a nr,j , where κ j and a nr,j are defined by Eqs. (22) and (48), respectively. Then from Eq. (53), one finds two solutions for constants which lead to two orthonormal bispinors, corresponding to two independent solutions for numbers (β 1 , β 2 ) (cp. (44) and (49)) one of which, β (+) ν , corresponds to the states with positive eigenvalue of BELI, BEL = κ j a nr,j , and the other one, β (−) ν , -with negative, BEL = − κ j a nr,j . The eigenstates, as before, are characterized by four quantum numbers (n r , j, m j , σ = ±).

General solution of the Dirac equation
Each set {Ĥ,Ĵ 2 ,Ĵ z ,Î inv }, whereÎ inv isK D ,Â JL orÎ BEL , results in its own system of eigen bispinors with certain values of the parameters β ν in the general solution (36). On the other hand, such a solution can be sought as a single eigen bispinor of the set with generalized invariantÎ which is given by arbitrary linear combination with coefficients c D , c JL and c BEL determining the weight of each invariant contribution. EqualityÎ gen Ψ nr,j,m j = gen Ψ nr,j,m j leads to the following equations for constants A and B: ( The condition for their non-trivial solution defines the eigenvalues gen = ± nr,j , nr,j = κ 2 j c 2 D + c 2 JL + κ 2 j c 2 BEL n r (n r + 2γ j ) N 2 and Eqs.(56) admit two solutions which correspond to positive and negative eigenvalues of the generalised invariant (55). Introducing angles θ and φ via relations where the first pair β (60) All sets for four independent operators include the invariantĴ 2 and differ by spin invariants only. Therefore, the DE can be exactly solved by variable separation in spherical coordinates. Obviously, the eigen bispinor of any of the invariants can be expanded over the eigen bispinors of another one. Due to anticommutation relations between spin invariants, an interconverting transformation includes only two bispinors corresponding to different signs of σ [20], and expression (60) is expansion of eigen bispinors ofÎ gen over the eigen bispinors of the DI.

Classification of quantum states
It is worth to note that in spite of several solutions, the eigen spectrum of the Dirac Hamiltonian with Coulomb interaction does not depend on the choice of the set {Ĥ D ,Ĵ z ,Ĵ 2 ,Î inv }. Indeed, electron states in this potential are defined by the full set of quantum numbers {n r , j, m j , σ}, where σ = ±. Instead of this, it is more convenient to use the principal quantum number [3] n = n r + κ j = 1, 2, ... and to characterize stationary states by the set of quantum numbers {n, j, m j , σ} where the principal quantum number n takes values n = 1, 2, . . ., quantum number of the total angular momentum j takes semi-integer values in the range 1/2 ≤ j ≤ n − 1/2 (1 ≤ κ j ≤ n), magnetic quantum number m j takes positive and negative semi-integer values in the range −j ≤ m j ≤ j and σ, in general, takes two values.
In this case the energy levels of the hydrogen spectrum, calculated from Eq. (28) with the accuracy of the lowest order with respect to the fine structure constant, take the form where W R = e 4 m/2 2 is the Rydberg energy and number κ j was defined in Eq. (22). Sub-levels with the maximal possible for the given principal quantum number, j = j max = n − 1/2, are very special in the fine structure. At these levels when κ jmax = n, the polynomials of the order n r = n − κ jmax = 0 become constant and according to definitions (37), polynomials P which do not depend on the choice of the eigen bispinors system and σ takes one value only. The corresponding energy levels of the fine structure are 2n times degenerate in accordance with the magnetic quantum number m j = ±1/2, . . . , ±(n − 1)/2. The rest n−1 sub-levels (at n ≥ 2, j < (2n−1)/2) of the fine structure of each multiplet are described by spinors (36) with polynomials of the order n r = n − κ j ≥ 1, in which m j takes 2(n − n r ) values, and σ can have both signs. The sub-level of the fine structure is degenerate with respect to the magnetic quantum number m j and spin number σ and, therefore, is 4(n − n r ) times degenerate. This means that 2n + 4 n−1 nr=1 (n − n r ) = 2n 2 states correspond to the same principal quantum number n, as it should be.
Using the found bispinors Ψ n,j,m j ,σ (r), it is possible to calculate the distribution of probability amplitude for the atomic orbitals w n,j,m j ,σ = Ψ † n,j,m j ,σ (r)Ψ n,j,m j ,σ (r) =| ψ n,j,m j ,σ | 2 , and their electron spin orientation n,j,m j ,σ .
In view of the smallness of the parameter Zα 1, contribution of the lower spinor of the order |Zα| 2 1 in these expressions is negligibly small and can be ignored. By the same reason the upper spinor can be expanded with respect to the small parameter, and terms ∼ Z 2 α 2 can be omitted also. Therefore, in calculations of the mean-values the lower spinor can be omitted. The upper spinor has the form where numbers β n,j (ρ n ), with the dimensionless radial variable ρ n = 2Zr/nr B (see (32)). Here the normalization constant and the polynomials in the limit Zα → 0 are defined in Eq. (33) and (37), respectively. In this approximation the function R (±) n,j coincides with the Schrödinger radial function except one difference: at j = 1/2 the value γ j is slightly less than 1 which leads to a small negative power of r. Of course, all integrals connected with the spectrum, are convergent. In fact, this corresponds to the non-relativistic approximation when an electron is described by the Pauli spinor (63). But it is necessary to keep in mind that excitations of the Dirac spinor field (particles and antiparticles) are described by the full 4-component field amplitude (bispinor) Ψ n,j,m j ,σ .
Electron states are characterised by the full set of quantum numbers {n, j, m j , σ}. In states (43) with the definite value of the DI with constants β n,j χ l,m j ,± involves spinors with integers l = 0, 1, 2, . . .. In view of the smallness of the lower spinor ∼ |Zα| 2 , the conventional spectroscopic classification of states by s, p, d, . . . based on number l, is used. But in our theory l defines the order of spherical harmonics in the main spinor and is not connected with the angular momentum operatorL 2 because angular momentum is no longer an integral of motion.
Worth mentioning is a peculiarity of the states on these levels which complete population of the n-th shell of the hydrogen-like spectrum. Hartree D.R. [27] has shown that summing of probabilities w n,j=n−1/2,m j by m j from m j = +1/2 up to m j = +j gives spherically symmetric charge distribution. In particular, one can see from Eq. (67) that w 2,3/2,±1/2 + w 2,3/2,±3/2 = R 2 20 (r) /2π. On the other sub-levels of the multiplet fine structure electron states are described by quantum numbers m j and σ which take two values. In these states the probability density and spin orientation significantly depend on the choice of the spin invariant in the set {Ĥ,Ĵ 2 ,Ĵ z ,Î inv }. For the states with definite value of the DI, the four-fold degenerate level E 2,1/2 includes 2S 1/2 states (n = 2, j = 1/2, m j = ±1/2, σ = +) which are described by spinors and 2P  cos ϑ e iϕ sin ϑ , In Ref. [28] the probability amplitude for Darwin solution is given in graphical presentation and it is shown that radial distributions, described by Shrödinger and Dirac functions, are almost identical. This means that the approximate expression via the Pauli spinors (63) reproduces |Ψ n,j,m j ,σ | 2 with high accuracy. For sub-level E 2,1/2 we have where e r and e ϑ are unit vectors in spherical coordinate system. It is seen, as it had been pointed out by H. E. White [28], "not only are all S states spherically symmetrical, as on the Schrödinger theory, but also one electron in a 2P 1/2 state". Nonetheless, the spin invariant itself is not determined a priori and an electron can, at least from the point of view of its energy, be in any of the discussed above states. Therefore, electron states should be described by the general solution (36) or (60) (or by Pauli spinors (63)) in the non-relativistic limit).
This means that on sub-levels with j < n−1/2 of the fine structure multiplet, the shape of electron cloud and the spin orientation in it take different forms. The deformation of electron cloud and changes in spin orientation are mutually connected and are controlled by variation of the values of β (σ) 1,2 , i.e., by parameters θ and φ in (59, that characterize spin degree of freedom in the spinor (63). In all cases the states with positive and negative values of the magnetic number m j have the same probability density distribution, but opposite spin orientations. For example, the states on sub-level E 2,1/2 are defined by Pauli spinors where the values of coefficients β In the case of JLI the symmetry with respect to reflection from the xy-plane perpendicular to the polar z−axis and, respectively, spherical symmetry are broken, charge density shifts along the polar axis in one or opposite direction depending on σ (see Fig.  5). Spin orientations s 2,1/2,m j ,σ (64) for each invariant are given by different expressions from which it follows that spin spatial dependencies differ significantly, as it is also clearly seen from the figures. Moreover, it follows from the above obtained expressions that for the states n = 2, j = 1/2, m = ±1/2 corresponding to the DI and JLI, spins lie in the plane (z, x), while in the case of BELI spins lie in this plane only at ϑ = 0, π, at other angles spin orientations have a component in ±e ϕ direction. This is shown in Fig. 6 for the state with m = 1/2, σ = + by arrows with symbols • and ×, respectively. In the state m = −1/2, σ = + spins satisfy the relation s BEL 2,1/2,−1/2,+ = −s BEL 2,1/2,1/2,+ , and, thus, they have opposite orientation to the one shown in Fig. 6 (c). In the state m = 1/2, σ = − spins have orientation in the plane (z, x) like the one, shown in Fig. 6 (c), but have opposite e ϕ component out of the plane (z, x). In the other state spin orientation satisfies the relation s BEL 2,1/2,−1/2,− = −s BEL 2,1/2,1/2,− .

Final remarks
From the above studies it follows that electron stationary states at the levels of hydrogen-like energy spectrum are determined by the full set of quantum numbers 1 {n, j, m j , σ}. It has been shown that the eigenvalue equationÎ gen Ψ nr,j,m j = gen Ψ nr,j,m j fixes free parameters in the general solution (60) and defines the two orthogonal eigen bispinors corresponding to the two signs σ = ± of the eigenvalue gen . As long as the free parameters in the general solution of the DE are not fixed, the spin vector is not defined. The free parameters and spin vector are determined by the choice of the particular invariantK D ,Â JL orÎ BEL . As a result, the spinors in the eigen bispinor become unambiguously determined and in this sense these invariants can be called "spin invariants" and σ -"spin quantum number" according to conventional meaning of this number.
Three quantum numbers n, j, m j correspond to spatial degree of freedom and σ to the spin one. In the non-relativistic limit electron is described by two-component Pauli spinor and σ is usually connected with the spin direction. Such statement is valid for systems with the translational symmetry when the integrals of motion are given by components of the momentum (quasi-momentum). In the central field when a constant of motion iŝ J z , the spin direction is determined by the sign of the magnetic quantum number.
Due to the hidden symmetry the energies (28) are degenerate with respect to magnetic m j and spin σ quantum numbers. External perturbations can break some symmetries and remove the degeneracy. For instance, an external magnetic field removes the degeneracy with respect to m j . Zeeman effect causes splitting of the corresponding states, in particular the 2P 3/2 state splits into four sub-levels, whereas each state 2S 1/2 and 2P 1/2 splits into two sub-levels. Some external factors (e.g., the interaction with quantum fluctuations of the vacuum electromagnetic field) remove the degeneracy with respect to σ and cause hyperfine splitting (the Lamb shift) of states with different σ: splitting between the states 2S 1/2 (σ = −) and 2P 1/2 (σ = +). The states of fine-structure multiplets with the maximal value j at the given principal quantum number are characterized by one value of the quantum number σ only and, therefore, stay unsplitted.
Finally, note, the expectation values can be calculated within the non-relativistic approximation with replacement of Dirac bispinors by Pauli spinors (63). Such calculations give the results which with high accuracy coincide with the exact results, and nevertheless, the important thing is that even in the non-relativistic approximation the electron is described by Dirac bispinor. Even though the spinors of this bispinor have different weights, neither of them can be ignored, and the complete description of an electron states requires both spinors. The DI has a diagonal block matrix structure, and the non-relativistic corrections are present in the lower spinor only. Its account results in the appearance in the non-relativistic approximation of the relativistic correction corresponding to the Thomas term, which in the case of the central field is known as spin-orbit interaction proportional to the product σ · L. If the matrix of a spin invariant is not diagonal, as it takes place for spin invariantsÎ JL andÎ BEL , relativistic corrections are present in both spinors, which results in another form of the spin-orbit interaction in the non-relativistic approximation. For example, as it can be easily shown, for the JLI the spin-orbit interaction operator acquires different form, namely, it is proportional to the product σ · A + (see (46)). Another example of non-conventional form of this interaction is quasi-two-dimensional electron gas [25].
Present study demonstrates spin state variability in the hydrogen-like spectrum. Since energies of the different spin states are degenerate, a small external factor, e.g., the one breaking the inversion symmetry, can be enough to cause electron transition from one spin orientation to another one. In this respect we recall that new technologies such as spintronics, spin chemistry etc., have attained a large amount of new data regarding the role of spin and its manifestation in various processes. Spin based effects take place in spincaloritronics devices, heat-to-spin conversion and energy harvesting systems, etc., spin-orbit torque devices for electric magnetization switching, as well as in the design of quantum materials [30,31,32], spin relativistic chemistry [33,34]. Precision measurements of spin related fine structure and Lamb shift in the anti-hydrogen spectrum have  been performed [29]. They can be used as tests of the charge-parity-time symmetry, as method of the determination of other fundamental quantities of such systems, etc.
Another field of spin sensitive processes is biological systems, in which electrons play essential role not only in the storage and transport of energy and information, but also in biorecognition connected with spin selectivity and chiral-induced spin selectivity. For instance, electron transfer through diamagnetic material (macromolecules, DNA) results in specific spin states of electrons [36,37,38,39] which is attributed to spin-orbit interaction. This underlines importance of further investigation of spin effects in different systems.    Fig.3 (c), showing spin orientations in the area of the highest electron probability for the state n = 2, j = 1/2, m = 1/2, σ = + for JLI. (b) -The insert of Fig.4 (c), showing spin orientations in the area of the highest electron probability for the state n = 2, j = 1/2, m = −1/2, σ = + for JLI. Figure 6: (a) -Spatial distribution of electron probability described by the DE with the BELI for the states n = 2, j = 1/2, m = ±1/2, σ = ± in the plane (z, ρ). (b) -A pattern of two 3D distributions of probabilities equal to the value 0.0015, in the Cartesian coordinate system. (c)-Charge 'isodensity' lines in the plane (z, x); arrows show spin orientation in the area of the highest charge density in the state m = 1/2, σ = +. Remarkably, with angle ϕ changing, spins acquire projection out of (z, x) plane in the up or down direction, indicated by symbols • and ×, respectively.