Calculates Maxwellian reaction rate for a n
(\rightarrow)m
transition due to a beam with atomic mass ab
and energy eb
firing into a target with atomic mass am
, temperature T
, and charge q
which has a cross section given by the function fqnm
Type | Intent | Optional | Attributes | Name | |||||||||||||||||||||||||||||||||||||
---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|---|
public function fqnm(a, b, c, d)Cross section function Arguments
Return Value real(kind=8) |
|||||||||||||||||||||||||||||||||||||||||
integer, | intent(in) | :: | q | Target charge |
|||||||||||||||||||||||||||||||||||||
real(kind=Float64), | intent(in) | :: | T | Target temperature [keV] |
|||||||||||||||||||||||||||||||||||||
real(kind=Float64), | intent(in) | :: | eb | Beam energy [keV] |
|||||||||||||||||||||||||||||||||||||
real(kind=Float64), | intent(in) | :: | am | Target atomic mass [amu] |
|||||||||||||||||||||||||||||||||||||
real(kind=Float64), | intent(in) | :: | ab | Beam atomic mass [amu] |
|||||||||||||||||||||||||||||||||||||
integer, | intent(in) | :: | n | Initial atomic energy level/state |
|||||||||||||||||||||||||||||||||||||
integer, | intent(in) | :: | m | Final atomic energy level/state |
|||||||||||||||||||||||||||||||||||||
real(kind=Float64), | intent(out) | :: | rate | Reaction Rate [] |
|||||||||||||||||||||||||||||||||||||
logical, | intent(in), | optional | :: | deexcit | Calculate de-excitation reaction rate |
subroutine bt_maxwellian_q_n_m(fqnm, q, T, eb, am, ab, n, m, rate, deexcit)
!+ Calculates Maxwellian reaction rate for a `n`\(\rightarrow)`m` transition due to a beam with atomic mass `ab` and energy `eb`
!+firing into a target with atomic mass `am`, temperature `T`, and charge `q` which has a cross section given by the function `fqnm`
interface
function fqnm(a, b, c, d)
!+Cross section function
real(8) :: fqnm !sigma
real(8), intent(in) :: a !eb
integer, intent(in) :: b !q
integer, intent(in) :: c !n
integer, intent(in) :: d !m
end function fqnm
end interface
integer, intent(in) :: q
!+Target charge
real(Float64), intent(in) :: T
!+Target temperature [keV]
real(Float64), intent(in) :: eb
!+Beam energy [keV]
real(Float64), intent(in) :: am
!+Target atomic mass [amu]
real(Float64), intent(in) :: ab
!+Beam atomic mass [amu]
integer, intent(in) :: n
!+Initial atomic energy level/state
integer, intent(in) :: m
!+Final atomic energy level/state
real(Float64), intent(out) :: rate
!+Reaction Rate [\(cm^3/s\)]
logical, intent(in), optional :: deexcit
!+Calculate de-excitation reaction rate
logical :: dxc
integer, parameter :: n_vr = 30
real(Float64) :: vr_max, dvr
real(Float64), dimension(n_vr) :: vr
real(Float64), dimension(n_vr) :: fr
integer, parameter :: n_vz = 60
real(Float64) :: vz_max,dvz
real(Float64), dimension(n_vz) :: vz
real(Float64), dimension(n_vz) :: fz
real(Float64) :: T_per_amu, eb_per_amu, ared, sig, sig_eff
real(Float64) :: zb, u2_to_erel, u2, erel, dE, factor, En, Em, v_therm
integer :: i, j
if(present(deexcit)) then
dxc = deexcit
else
dxc = .False.
endif
vr_max = 4.d0
dvr = vr_max/(n_vr - 1.d0)
do i=1,n_vr
vr(i) = (i-1)*dvr
enddo
vz_max = 4.d0
dvz = 2.0*vz_max/(n_vz - 1.d0)
do i=1,n_vz
vz(i) = (i-1)*dvz - vz_max
enddo
En = (13.6d-3)*(1.0 - (1.d0/n)**2.0)
Em = (13.6d-3)*(1.0 - (1.d0/m)**2.0)
dE = Em - En
T_per_amu = max(T, 1.d-6)/am
eb_per_amu = eb/ab
ared = am*ab/(am + ab)
v_therm = 1.384d6 * sqrt(T_per_amu*1.d3)
zb = sqrt(eb_per_amu/T_per_amu)
u2_to_erel = ared*T_per_amu
if(ared.lt.0.5) ared = 1.0
fz = 0.d0
fr = 0.d0
do i=1,n_vz
do j=1,n_vr
u2 = (zb - vz(i))**2.0 + vr(j)**2.0
erel = u2_to_erel*u2
if(dxc) then
factor = (erel + dE)/erel
erel = erel + dE
else
factor = 1.0
endif
if(erel.ge.dE) then
sig = fqnm(erel/ared, q, n, m)
else
sig = 0.d0
endif
fr(j) = factor*sig*sqrt(u2)*exp(-(vz(i)**2.0 + vr(j)**2.0))*vr(j)
enddo
fz(i) = simpsons_rule(fr, dvr)
enddo
sig_eff = (2.0/sqrt(PI))*simpsons_rule(fz, dvz)
rate = sig_eff*v_therm
if(dxc) rate = rate*(real(n)/real(m))**2.0
end subroutine bt_maxwellian_q_n_m