\beta 
\alpha - \beta \alpha \alpha$ and $\chi_2 = 2 \alpha \alpha \beta  - \beta \alpha \alpha - \alpha \beta \alpha$. Here the notations $\alpha$ and $\beta$ are the one-electron spin-up and spin-down functions, respectively (their definition can be found, e.g., in [@Dir]). The notations $A$ and $B$ in Eq.(\[psi\]) mean that the two sets of non-linear parameters associated with radial functions $\psi$ and $\phi$ can be optimized independently. In the general case, each of the radial basis functions explicitly depends upon all six interparticle (relative) coordinates $r_{12}, r_{13}, r_{23}, r_{14}, r_{24}, r_{34}$. It is clear that in actual bound state calculations only one spin function, e.g., the $\chi_1$ function, can be used. Note also that when the corresponding radial function has been constructed, then one can use an artifice called the ‘doubling’ of the wave function (see, e.g., [@Fro2011]). This trick is based on the use of the same set of non-linear parameters in the two radial functions in Eq.(\[psi\]). Obviously, this is not an optimal choice of the non-linear parameters in the wave function, but in many cases this simple handling leads to a substantial improvement of the answer.

The semi-exponential variational expansion of the radial function $\psi_{L=0}(A; 