Calculates an array of cross sections for a proton-Hydrogen impact excitation transitions from
the state to the state at energy eb
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
real(kind=Float64), | intent(in) | :: | eb | Relative collision energy [keV/amu] |
||
integer, | intent(in) | :: | m_max | Number of m states to calculate |
Array of cross sections where the m'th index refer to the transition from to m []
function p_excit_2_janev(eb, m_max) result(sigma)
!+Calculates an array of cross sections for a proton-Hydrogen impact excitation transitions from
!+the \(n=2\) state to the \(m=1..{m_max}\) state at energy `eb`
!+
!+###Equation
!+$$ H^+ + H(2) \rightarrow H^+ + H(m=1..m_{max}), m \gt 2$$
!+
!+###References
!+* Eq. 32 and Table 6 in Ref. 2 for \(m \le 5\) [[atomic_tables(module)]]
!+* Eq. 33 and Table 6 in Ref. 2 for \(m = 6-10\) [[atomic_tables(module)]]
!+* Eq. 34 and Table 6 in Ref. 2 for \(m \gt 10\) [[atomic_tables(module)]]
!+
real(Float64), intent(in) :: eb
!+ Relative collision energy [keV/amu]
integer, intent(in) :: m_max
!+ Number of m states to calculate
real(Float64), dimension(m_max) :: sigma
!+ Array of cross sections where the m'th index refer to the transition
!+ from \(n=2\) to m [\(cm^2\)]
integer :: m
sigma = 0.d0
do m=1,m_max
select case (m)
case (1)
sigma(1) = 0.d0
case (2)
sigma(2) = 0.d0
case (3)
sigma(3) = p_excit_2_3_janev(eb)
case (4)
sigma(4) = p_excit_2_4_janev(eb)
case (5)
sigma(5) = p_excit_2_5_janev(eb)
case (6)
sigma(6) = p_excit_2_6_janev(eb)
case (7)
sigma(7) = p_excit_2_7_janev(eb)
case (8)
sigma(8) = p_excit_2_8_janev(eb)
case (9)
sigma(9) = p_excit_2_9_janev(eb)
case (10)
sigma(10) = p_excit_2_10_janev(eb)
case DEFAULT
sigma(m) = sigma(10)*(10.0/real(m))**3.0
end select
enddo
end function p_excit_2_janev