Hex
1.0
Hydrogenelectron collision solver

HexECS computes partial Tmatrices for elastic, excitation and ionization electronhydrogen scattering. Together with the interface program HexDB it offers the possibility of generating many scatterign quantities like differential or integral cross sections.
Hex is written in C++11 to make use of comfort of the modern C++ extensions, so one may need a newer compiler. Tested compilers are:
Both worked with the same Makefile, just by setting the variable CPP to "g++" or "icpc". The program also uses following external packages (tested versions are given in parentheses):
The next libraries are optional:
All listed libraries are opensource and easily obtainable on the internet and in the repositories of some Linux distributions (openSUSE 12.3 at least).
Equations in this documentation use MathJax, which should work in every uptodate JavaScriptenabled web browser. Tested browsers are:
Unknown scattering state, eigenfunction of the full system hamiltonian, is split into two parts, asymptotic “incoming particle” state and the scattered part, which is a solution of the driven Schrödinger equation,
\[ \left(E  \hat{H}\right) \Psi_{\mathrm{sc}} = \hat{H}_{\mathrm{int}} \Psi_{\mathrm{inc}} \ . \]  (1) 
The state \( \Psi_{\mathrm{inc}} \) is a product of initial atomic state and a projectile plane wave,
\[ \Psi_{\mathrm{inc}}(\mathbf{r}_1, \mathbf{r}_2) = \frac{1}{k_i r_1 r_2} P_{n_i l_i}(r_1) Y_{l_i m_i}(\mathbf{\hat{r}}) \cdot (2\pi)^{3/2} \sum_{lm} \sqrt{4\pi} \mathrm{i}^l \hat{j}_l(k_i r_2) Y_{lm}(\mathbf{\hat{k}}_i) Y_{lm}^\ast(\mathbf{\hat{r}}_2) \ , \]  (2) 
and anti/symmetrized with respect to electron exchange, according to total spin,
\[ \Psi^S_{\mathrm{inc}}(\mathbf{r}_1, \mathbf{r}_2) = \Psi_{\mathrm{inc}}(\mathbf{r}_1, \mathbf{r}_2) + (1)^S \Psi_{\mathrm{inc}}(\mathbf{r}_2, \mathbf{r}_1) \ . \]  (3) 
Hamiltonian \( \hat{H} \) is a sum of free hamiltonian and interaction hamiltonian. Free hamiltonian is the hamilton operator for electron and atom separated far away,
\[ \hat{H}_{\mathrm{free}} = \frac{\nabla_1^2}{2}  \frac{\nabla_2^2}{2} \frac{1}{r_1} \ , \]  (4) 
whereas the interaction hamiltonian contains the interaction between projectile and the atomic constituents, that is
\[ \hat{H}_{\mathrm{int}} = \frac{1}{r_{12}}  \frac{1}{r_2} \]  (5a) 
for direct case and
\[ \hat{H}_{\mathrm{int}} = \frac{1}{r_{12}}  \frac{1}{r_1} \]  (5b) 
for exchange (anti/symmetrized) case. Having the solution of equation (1), we can extract the scattering amplitude, which is a matrix element of \( \hat{H}_{\mathrm{int}} \) between the solution \( \Psi = \Psi_{\mathrm{sc}} + \Psi_{\mathrm{inc}} \) and some asymptotic state \( \Psi_{\mathrm{out}} \) that is experimentally inquired. Conventionally, \( \Psi_{\mathrm{out}} \) is also a product of (now final) atomic state, possibly different from the original, with an outgoing plane wave of the projectile. Any exchange effect are already included in the original anti/symmetrization, so they need not be considered now. The projection
\[ f = 4\pi^2 \left<\Psi_{\mathrm{out}}\right\hat{H}_{\mathrm{int}}\left\Psi\right> \]  (6) 
can be also written as
\[ f = 4\pi^2 \left<\Psi_{\mathrm{out}}\rightE  \hat{H}_{\mathrm{free}}\left\Psi_{\mathrm{sc}}\right> \ , \]  (7) 
where the following was used: \( \hat{H}_{\mathrm{int}} = \hat{H}  \hat{H}_{\mathrm{free}} \), \( \Psi = \Psi_{\mathrm{sc}} + \Psi_{\mathrm{inc}} \) and \( \hat{H}\Psi = E\Psi \).
All computations are done in timeindependent way, and the exterior complex scaling is used instead of boundary condition fitting. Radial part of sought wavefunctions is expanded in a Bspline basis \( \left\{B_i\right\}_{i = 0}^{\mathrm{Nspline}1} \) of a given order,
\[ \Psi_{\mathrm{sc}}^{LM}(\mathbf{r}_1, \mathbf{r}_2) = \sum_{l_1 l_2} \psi_{l_1 l_2}^{LM}(\mathrm{r}_1, \mathrm{r}_2) \mathcal{Y}_{l_1 l_2}^{LM}(\mathbf{\hat{r}}_1, \mathbf{\hat{r}}_2) \ , \]  (8) 
\[ \psi^{LM}_{\mathrm{sc},l_1 l_2}(\mathrm{r}_1, \mathrm{r}_2) = \frac{1}{r_1 r_2} \sum_{ij} \psi_{l_1 l_2,ij}^{LM} B_i(r_1) B_j(r_2) \ , \]  (9) 
and when projecting the equation (1) on a bipolar spherical function \( \left<\mathcal{Y}_{l_1 l_2}^{LM}\right \) to get rid of angular dependence, and on a pair of Bsplines to get rid of radial dependence and keep only matrix elements, one arrives at a matrix equation for components of \( \psi_{l_1 l_2}^{LM}(\mathrm{r}_1,\mathrm{r}_2) \),
\[ \left[\left(ES_{ik}S_{jl}  \frac{1}{2} D_{ik}S_{jl}  \frac{1}{2} S_{ik} D_{jl}  \frac{1}{2}l_1(l_1+1) M_{ik}^{(2)}S_{jl}  \frac{1}{2}l_2(l_2+1) S_{ik}M_{jl}^{(2)} + M_{ik}^{(1)}S_{jl} + S_{ik}M_{jl}^{(1)}\right) \delta_{l_1}^{l_1'}\delta_{l_2}^{l_2'}  \sum_{\lambda} f_{l_1l_2l_1'l_2';L}^\lambda R_{ijkl}^\lambda \right] \psi_{l_1' l_2',kl}^{LMS} = \chi_{l_1l_2,kl}^{LMS} \ , \]  (10) 
or symbolically
\[ \left[ \mathsf{Id}_1 \otimes \mathsf{Id}_2 \otimes \left(E\mathsf{S}\otimes\mathsf{S}  \frac{1}{2} \mathsf{D}\otimes\mathsf{S}  \frac{1}{2} \mathsf{S}\otimes\mathsf{D}  \frac{1}{2} l_1 (l_1+1) \mathsf{M}^{(2)}\otimes\mathsf{S}  \frac{1}{2} l_2 (l_2+1) \mathsf{S}\otimes\mathsf{M}^{(2)} + \mathsf{M}^{(1)} \otimes \mathsf{S} + \mathsf{S} \otimes \mathsf{M}^{(1)}\right) + \sum_{\lambda} \mathsf{f}_L^\lambda \otimes \mathsf{R}^\lambda \right] \mathsf{\psi}^{LMS} = \mathsf{\chi}^{LMS} \ . \]  (10*) 
Symbol \( \otimes \) stands for Kronecker product (“flattened tensor product”) and matrices have following meanings:
\[ \left<B_i\right \left(\frac{\mathrm{d}^2}{\mathrm{d} x^2}\right) \leftB_k\right> = + \int_a^b \frac{\mathrm{d}B_i}{\mathrm{d}x} \frac{\mathrm{d}B_k}{\mathrm{d}x} \equiv + D_{ik} \ . \]
\[ S_{ik} = \left<B_i\right\left.\!B_k\right> \ . \]
\[ R_{ijkl}^\lambda = \int_a^b\int_a^b B_i(r_1) B_j(r_2) \frac{r_<^\lambda}{r_>^{\lambda + 1}} B_k(r_1) B_l(r_2) dr_1 dr_2 \ , \]
flattened so that “i” and “j” form one multiindex [ij] and the other indices the multiindex [kl].\[ \left<l_1 l_2 \right\! \frac{1}{r_{12}} \!\leftl_1' l_2'\right>_L = \sum_{\lambda} f_{l_1 l_2 l_1' l_2'; L}^\lambda \frac{r_<^\lambda}{r_>^{\lambda + 1}} \ . \]
flattened so that “l₁” and “l₂” form one multiindex [l₁l₂] and the other two indices the second.Finally, the symbol \( \chi^{LMS} \) stands for projection of the right hand side, which is
\[ \chi_{l_1 l_2, ij}^{LMS} = \frac{1}{k_f} \sum_\ell \mathrm{i}^\ell \sqrt{2\pi(2\ell+1)} \left\{ \left(\sum_\lambda f_{l_1 l_2 l_i \ell; L}^\lambda R_{ijkl}^\lambda  \delta_{l_1}^{l_i} \delta_{l_2}^{\ell} S_{ik}M_{jl}^{(1)}\right) \left[P_{n_il_i}(r_1)\right]_k \left[\hat{j}_\ell(k_i r_2)\right]_l + (1)^{S+\Pi} \left(\sum_\lambda f_{l_1 l_2 \ell l_i; L}^\lambda R_{ijkl}^\lambda  \delta_{l_1}^{\ell} \delta_{l_2}^{l_i} M_{ik}^{(1)} S_{jl}\right) \left[\hat{j}_\ell(k_i r_1)\right]_k \left[P_{n_il_i}(r_2)\right]_l \right\} \]  (11) 
or symbolically
\[ \chi^{LMS} = \frac{1}{k_f} \sum_\ell \mathrm{i}^\ell \sqrt{2\pi(2\ell+1)} C_{l_i m_i \ell 0}^{L M} \left\{ \left(\sum_\lambda \mathsf{f}_L^\lambda \otimes \mathsf{R}^\lambda  \mathsf{Id}_1 \otimes \mathsf{Id}_2 \otimes \mathsf{S} \otimes \mathsf{M}^{(1)}\right) \cdot \Delta_{l_1}^{l_i} \otimes \Delta_{l_2}^{\ell} \otimes \mathsf{P}_{n_i l_i} \otimes \mathsf{j}_{\ell,k_i} + (1)^{S+\Pi} \left(\sum_\lambda \mathsf{f}_L^\lambda \otimes \mathsf{R}^\lambda  \mathsf{Id}_1 \otimes \mathsf{Id}_2 \otimes \mathsf{M}^{(1)} \otimes \mathsf{S}\right) \cdot \Delta_{l_1}^{\ell} \otimes \Delta_{l_2}^{l_i} \otimes \mathsf{j}_{\ell,k_i} \otimes \mathsf{P}_{n_i l_i} \right\} \ . \]  (11*) 
Here, the onedimensional (!) vectors \( \Delta \) are zero vectors with only one element equal to one at position \( l_1 = l_i \) (etc.). P and j vectors are components of respective function in chosen Bspline basis.
Factor C in the expression (11*) is the ClebschGordan coefficient and the zero projection of \( \ell \)momentum reflect the deliberate choice of scattering axis along the projectile moemntum, so that the angular momentum projection is zero.
Exterior complex scaling of right hand side poses a serious problem for typical (not exponentially decreasing) potentials. One of the factors in the right hand side is the RiccatiBessel functions \( \hat{j} \), which exponentially diverges under ECS transformation. To avoid this, potential \( \hat{H}_{\mathrm{int}} \) is artificially truncated at (or before) the turning point \( R_0 \), which has several consequences for the numerical construction:
As was said above, the scattering amplitude is
\[ f = 4\pi^2 \left<\Psi_{\mathrm{out}}\rightE  \hat{H}_{\mathrm{free}}\left\Psi_{\mathrm{sc}}\right>_{R_0}\ , \]  (12) 
where the subscript \( R_0 \) means, that the radial integration is done only for radii less than \( R_0 \) (because the original matrix, \( \hat{H}_{\mathrm{int}} \), has to be truncated at such distance to avoid farregion divergence. The outgoing – detected – wavefunction \( \Psi_{\mathrm{out}} \) has form simillar to (2), with replaced initial to final quantum numbers. Substituting such expansion into the equation (12) and once again using zero boundary trait of chosen Bspline basis when doing per parts integration one easily arrives at the formula for cross section
\[ \sigma = \frac{4}{k_i k_f} \sum_{\ell L L'} C_{l_f m_f \ell 0}^{L m_f} C_{l_f m_f \ell 0}^{L' m_f} \left \psi_{l_f \ell, ij}^{LMS} W[P]_i S[j]_j + \psi_{l_f \ell, ij}^{L'MS} S[P]_i W[j]_j \right^2 \ , \]  (13) 
where \( W[P]_i \) and \( W[j]_j \) stand for wronskian (evaluated at \( R_0  \varepsilon \)) of the ith Bspline and the (final) hydrogenic or RiccattiBessel function, and \( S[P]_i \) and \( S[j]_i \) stand for (truncated) overlap integrals of the ith Bspline and the (final) hydrogenic or RiccattiBessel function.
The code runs along the following outline:
mpirun np 4 hexecs mpi