Reads in a nuclear reaction rates table from file and puts it into a NuclearRates type
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer(kind=HID_T), | intent(in) | :: | fid | HDF5 file ID |
||
character(len=*), | intent(in) | :: | grp | HDF5 group to read from |
||
type(NuclearRates), | intent(inout) | :: | rates | Atomic reaction rates |
subroutine read_nuclear_rates(fid, grp, rates)
!+ Reads in a nuclear reaction rates table from file
!+ and puts it into a [[NuclearRates]] type
integer(HID_T), intent(in) :: fid
!+ HDF5 file ID
character(len=*), intent(in) :: grp
!+ HDF5 group to read from
type(NuclearRates), intent(inout) :: rates
!+ Atomic reaction rates
integer(HSIZE_T), dimension(1) :: dim1
integer(HSIZE_T), dimension(3) :: dim3
logical :: path_valid, err
integer :: i, j, error
real(Float64) :: emin, emax, tmin, tmax, rmin
err = .False.
call h5ltpath_valid_f(fid, grp, .True., path_valid, error)
if(.not.path_valid) then
if(inputs%verbose.ge.0) then
write(*,'(a,a)') 'READ_NUCLEAR_RATES: Unknown nuclear interaction: ', trim(grp)
write(*,'(a)') 'Continuing without neutron calculation'
endif
inputs%calc_neutron=0
return
endif
dim1 = [2]
call h5ltread_dataset_double_f(fid, grp//"/bt_amu", rates%bt_amu, dim1, error)
if(abs(inputs%ab-rates%bt_amu(1)).gt.0.2) then
if(inputs%verbose.ge.0) then
write(*,'(a,f6.3,a,f6.3,a)') 'READ_NUCLEAR_RATES: Unexpected beam species mass. Expected ',&
rates%bt_amu(1),' amu got ', inputs%ab, ' amu'
endif
err = .True.
endif
if(abs(inputs%ai-rates%bt_amu(2)).gt.0.2) then
if(inputs%verbose.ge.0) then
write(*,'(a,f6.3,a,f6.3,a)') 'READ_NUCLEAR_RATES: Unexpected thermal species mass. Expected ',&
rates%bt_amu(2),' amu got ', inputs%ai, ' amu'
endif
err = .True.
endif
if(err) then
if(inputs%verbose.ge.0) then
write(*,'(a)') 'Continuing without neutron calculation'
endif
inputs%calc_neutron=0
return
endif
call h5ltread_dataset_int_scalar_f(fid, grp//"/nbranch", rates%nbranch, error)
call h5ltread_dataset_int_scalar_f(fid, grp//"/nenergy", rates%nenergy, error)
call h5ltread_dataset_double_scalar_f(fid, grp//"/emin", emin, error)
call h5ltread_dataset_double_scalar_f(fid, grp//"/emax", emax, error)
call h5ltread_dataset_double_scalar_f(fid, grp//"/dlogE", rates%dlogE, error)
call h5ltread_dataset_int_scalar_f(fid, grp//"/ntemp", rates%ntemp, error)
call h5ltread_dataset_double_scalar_f(fid, grp//"/tmin", tmin, error)
call h5ltread_dataset_double_scalar_f(fid, grp//"/tmax", tmax, error)
call h5ltread_dataset_double_scalar_f(fid, grp//"/dlogT", rates%dlogT, error)
rates%logemin = log10(emin)
rates%logemax = log10(emax)
rates%logtmin = log10(tmin)
rates%logtmax = log10(tmax)
allocate(rates%log_rate(rates%nenergy, &
rates%ntemp, &
rates%nbranch))
dim3 = [rates%nenergy, rates%ntemp, rates%nbranch]
call h5ltread_dataset_double_f(fid, grp//"/fusion", rates%log_rate, dim3, error)
rmin = minval(rates%log_rate, rates%log_rate.gt.0.d0)
where (rates%log_rate.le.0.d0)
rates%log_rate = 0.9*rmin
end where
rates%minlog_rate = log10(rmin)
rates%log_rate = log10(rates%log_rate)
end subroutine read_nuclear_rates