Write Hydrogen-Hydrogen interaction cross sections to a HDF5 file
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer(kind=HID_T), | intent(inout) | :: | id | HDF5 file or group object id |
||
character(len=*), | intent(in) | :: | namelist_file | Namelist file that contains settings |
||
integer, | intent(in) | :: | n_max | Number of initial atomic energy states to calculate |
||
integer, | intent(in) | :: | m_max | Number of final atomic energy states to calculate |
subroutine write_bb_H_H(id, namelist_file, n_max, m_max)
!+ Write Hydrogen-Hydrogen interaction cross sections to a HDF5 file
integer(HID_T), intent(inout) :: id
!+ HDF5 file or group object id
character(len=*), intent(in) :: namelist_file
!+ Namelist file that contains settings
integer, intent(in) :: n_max
!+ Number of initial atomic energy states to calculate
integer, intent(in) :: m_max
!+ Number of final atomic energy states to calculate
logical :: calculate
real(Float64) :: emin
real(Float64) :: emax
integer :: nenergy
real(Float64) :: eb
real(Float64) :: dlogE
real(Float64), dimension(:), allocatable :: ebarr
real(Float64), dimension(:,:), allocatable :: ioniz
real(Float64), dimension(:,:,:), allocatable :: cx, excit
integer(HID_T) :: gid
integer(HSIZE_T), dimension(1) :: dim1
integer(HSIZE_T), dimension(2) :: dim2
integer(HSIZE_T), dimension(3) :: dim3
integer :: i, cnt, error
logical :: exis
NAMELIST /H_H_cross/ calculate, nenergy, emin, emax
calculate = .True.; nenergy = 200; emin = 1.d-3 ; emax = 8.d2
inquire(file=namelist_file,exist=exis)
if(.not.exis) then
write(*,'(a,a)') 'WRITE_BB_H_H: Input file does not exist: ',trim(namelist_file)
write(*,'(a)') 'Continuing with default settings...'
else
open(13,file=namelist_file)
read(13,NML=H_H_cross)
close(13)
endif
if(.not.calculate) return
allocate(ebarr(nenergy))
allocate(ioniz(n_max,nenergy))
allocate(cx(n_max,m_max,nenergy))
allocate(excit(n_max,m_max,nenergy))
ebarr = 0.d0
ioniz = 0.d0
cx = 0.d0
excit = 0.d0
if(verbose) then
write(*,'(a)') "---- H-H cross sections settings ----"
write(*,'(T2,"Emin = ",e9.2, " keV")') emin
write(*,'(T2,"Emax = ",e9.2, " keV")') emax
write(*,'(T2,"Nenergy = ", i4)') nenergy
write(*,*) ''
endif
cnt = 0
dlogE = (log10(emax) - log10(emin))/(nenergy - 1)
!$OMP PARALLEL DO private(i, eb)
do i=istart, nenergy, istep
eb = 10.d0**(log10(emin) + (i-1)*dlogE)
ebarr(i) = eb
cx(:,:,i) = p_cx(eb, n_max, m_max)
excit(:,:,i) = p_excit(eb, n_max, m_max)
ioniz(:,i) = p_ioniz(eb, n_max)
cnt = cnt + 1
if(verbose) WRITE(*,'(f7.2,"%",a,$)') 100*cnt*istep/real(nenergy),char(13)
enddo
!$OMP END PARALLEL DO
#ifdef _MPI
call parallel_sum(ebarr)
call parallel_sum(cx)
call parallel_sum(excit)
call parallel_sum(ioniz)
#endif
if(verbose) then
call h5gcreate_f(id, "H_H", gid, error)
dim1 = [1]
dim2 = [n_max, nenergy]
dim3 = [n_max, m_max, nenergy]
call h5ltmake_dataset_int_f(gid, "nenergy", 0, dim1, [nenergy], error)
call h5ltmake_dataset_int_f(gid, "n_max", 0, dim1, [n_max], error)
call h5ltmake_dataset_int_f(gid, "m_max", 0, dim1, [m_max], error)
call h5ltmake_dataset_double_f(gid, "dlogE", 0, dim1, [dlogE], error)
call h5ltmake_dataset_double_f(gid, "emin", 0, dim1, [emin], error)
call h5ltmake_dataset_double_f(gid, "emax", 0, dim1, [emax], error)
dim1 = [nenergy]
call h5ltmake_compressed_dataset_double_f(gid, "energy", 1, dim1, ebarr, error)
call h5ltmake_compressed_dataset_double_f(gid, "cx", 3, dim3, cx, error)
call h5ltmake_compressed_dataset_double_f(gid, "ionization", 2, dim2, ioniz, error)
call h5ltmake_compressed_dataset_double_f(gid, "excitation", 3, dim3, excit, error)
call h5ltset_attribute_string_f(id, "H_H", "description", &
"Cross sections for Hydrogen-Hydrogen interactions", error)
call h5ltset_attribute_string_f(gid, "nenergy", "description", &
"Number of nucleon energy values", error)
call h5ltset_attribute_string_f(gid, "n_max", "description", &
"Number of initial energy levels", error)
call h5ltset_attribute_string_f(gid, "m_max", "description", &
"Number of final energy levels", error)
call h5ltset_attribute_string_f(gid, "energy", "description", &
"Nucleon energy values", error)
call h5ltset_attribute_string_f(gid, "energy", "units", "keV/amu", error)
call h5ltset_attribute_string_f(gid, "dlogE", "description", &
"Energy spacing in log-10", error)
call h5ltset_attribute_string_f(gid, "dlogE", "units", "log10(keV/amu)", error)
call h5ltset_attribute_string_f(gid, "emin","description", &
"Minimum energy", error)
call h5ltset_attribute_string_f(gid, "emin", "units", "keV/amu", error)
call h5ltset_attribute_string_f(gid, "emax","description", &
"Maximum energy", error)
call h5ltset_attribute_string_f(gid, "emax", "units", "keV/amu", error)
call h5ltset_attribute_string_f(gid, "cx", "description", &
"n/m resolved charge exchange cross sections: cx(n,m,energy)", error)
call h5ltset_attribute_string_f(gid, "cx", "units", "cm^2", error)
call h5ltset_attribute_string_f(gid, "cx", "reaction", &
"H(+) + H(n) -> H(m) + H(+)", error)
call h5ltset_attribute_string_f(gid, "excitation", "description", &
"n/m resolved excitation cross sections: excitation(n,m,energy)", error)
call h5ltset_attribute_string_f(gid, "excitation", "units", "cm^2", error)
call h5ltset_attribute_string_f(gid, "excitation", "reaction", &
"H(+) + H(n) -> H(+) + H(m), m > n", error)
call h5ltset_attribute_string_f(gid, "ionization", "description", &
"n resolved ionization cross sections: ionization(n,energy)", error)
call h5ltset_attribute_string_f(gid, "ionization", "units", "cm^2", error)
call h5ltset_attribute_string_f(gid, "ionization", "reaction", &
"H(+) + H(n) -> H(+) + H(+) + e-", error)
call h5gclose_f(gid, error)
endif
deallocate(ebarr, cx, excit, ioniz)
end subroutine write_bb_H_H