write_bt_H_Aq Subroutine

public subroutine write_bt_H_Aq(id, namelist_file, n_max, m_max)

Write Hydrogen-Impurity reaction rates to a HDF5 file

Arguments

Type IntentOptional AttributesName
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


Calls

proc~~write_bt_h_aq~~CallsGraph proc~write_bt_h_aq write_bt_H_Aq h5ltmake_dataset_int_f h5ltmake_dataset_int_f proc~write_bt_h_aq->h5ltmake_dataset_int_f h5ltset_attribute_string_f h5ltset_attribute_string_f proc~write_bt_h_aq->h5ltset_attribute_string_f h5ltmake_dataset_double_f h5ltmake_dataset_double_f proc~write_bt_h_aq->h5ltmake_dataset_double_f interface~h5ltmake_compressed_dataset_double_f h5ltmake_compressed_dataset_double_f proc~write_bt_h_aq->interface~h5ltmake_compressed_dataset_double_f h5gclose_f h5gclose_f proc~write_bt_h_aq->h5gclose_f proc~h5ltmake_compressed_dataset_double_f_7 h5ltmake_compressed_dataset_double_f_7 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_7 proc~h5ltmake_compressed_dataset_double_f_6 h5ltmake_compressed_dataset_double_f_6 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_6 proc~h5ltmake_compressed_dataset_double_f_5 h5ltmake_compressed_dataset_double_f_5 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_5 proc~h5ltmake_compressed_dataset_double_f_4 h5ltmake_compressed_dataset_double_f_4 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_4 proc~h5ltmake_compressed_dataset_double_f_3 h5ltmake_compressed_dataset_double_f_3 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_3 proc~h5ltmake_compressed_dataset_double_f_2 h5ltmake_compressed_dataset_double_f_2 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_2 proc~h5ltmake_compressed_dataset_double_f_1 h5ltmake_compressed_dataset_double_f_1 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_1 proc~h5ltmake_compressed_dataset_double_f_7->h5ltmake_dataset_double_f h5dwrite_f h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_7->h5dwrite_f h5dclose_f h5dclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5dclose_f h5sclose_f h5sclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5sclose_f h5pset_chunk_f h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_7->h5pset_chunk_f h5pcreate_f h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_7->h5pcreate_f proc~chunk_size chunk_size proc~h5ltmake_compressed_dataset_double_f_7->proc~chunk_size h5pset_deflate_f h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_7->h5pset_deflate_f h5pset_shuffle_f h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_7->h5pset_shuffle_f h5screate_simple_f h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_7->h5screate_simple_f h5pclose_f h5pclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5pclose_f h5dcreate_f h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_7->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_6->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_6->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_6->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_6->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_6->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_6->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_6->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_6->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_5->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_5->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_5->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_5->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_5->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_5->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_5->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_5->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_5->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_4->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_4->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_4->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_4->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_4->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_4->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_4->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_4->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_3->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_3->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_3->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_3->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_3->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_3->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_3->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_3->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_3->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_2->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_2->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_2->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_2->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_2->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_2->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_2->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_2->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_2->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_1->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_1->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_1->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_1->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_1->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_1->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_1->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_1->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_1->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5dcreate_f

Called by

proc~~write_bt_h_aq~~CalledByGraph proc~write_bt_h_aq write_bt_H_Aq program~generate_tables generate_tables program~generate_tables->proc~write_bt_h_aq

Contents

Source Code


Source Code

subroutine write_bt_H_Aq(id, namelist_file, n_max, m_max)
    !+ Write Hydrogen-Impurity reaction rates 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

    integer :: q
    real(Float64) :: mass

    real(Float64) :: emin
    real(Float64) :: emax
    integer :: nenergy

    real(Float64) :: tmin
    real(Float64) :: tmax
    integer :: ntemp

    real(Float64) :: eb
    real(Float64) :: dlogE
    real(Float64), dimension(:), allocatable :: ebarr

    real(Float64) :: ti
    real(Float64) :: dlogT
    real(Float64), dimension(:), allocatable :: tarr

    real(Float64), dimension(:,:,:,:), allocatable :: ioniz, cx
    real(Float64), dimension(:,:,:,:,:), allocatable :: excit

    integer(HID_T) :: gid
    integer(HSIZE_T), dimension(1) :: dim1
    integer(HSIZE_T), dimension(2) :: dim2
    integer(HSIZE_T), dimension(3) :: dim3
    integer(HSIZE_T), dimension(4) :: dim4
    integer(HSIZE_T), dimension(5) :: dim5

    integer :: ie, it, ia, n, m, error, cnt
    real(Float64) :: rate
    integer, parameter :: n_bt_amu = 2
    real(Float64), dimension(2,n_bt_amu) :: a

    character(len=10) :: aname
    character(len=5) :: asym
    logical :: exis

    NAMELIST /H_Aq_rates/ q, mass, nenergy, emin, emax, ntemp, tmin, tmax

    q = 6 ; mass = C_amu
    nenergy = 100; emin = 1.d-3 ; emax = 4.d2
    ntemp = 100; tmin = 1.d-3 ; tmax = 2.d1
    inquire(file=namelist_file,exist=exis)
    if(.not.exis) then
        write(*,'(a,a)') 'WRITE_BT_H_Aq: Input file does not exist: ',trim(namelist_file)
        write(*,'(a)') 'Continuing with default settings...'
    else
        open(13,file=namelist_file)
        read(13,NML=H_Aq_rates)
        close(13)
    endif

    allocate(ebarr(nenergy))
    allocate(tarr(ntemp))
    allocate(ioniz(n_max,nenergy,ntemp,n_bt_amu))
    allocate(cx(n_max,nenergy,ntemp,n_bt_amu))
    allocate(excit(n_max,m_max,nenergy,ntemp,n_bt_amu))


    select case (q)
      case (5)
          aname = "Boron"
          asym = "H_B5"
      case (6)
          aname = "Carbon"
          asym = "H_C6"
      case DEFAULT
          aname = "Impurity"
          asym = "H_Aq"
    end select

    ebarr = 0.d0
    ioniz = 0.d0
    cx = 0.d0
    excit = 0.d0
    a(:,1) = [H1_amu, mass]
    a(:,2) = [H2_amu, mass]

    dlogE = (log10(emax) - log10(emin))/(nenergy - 1)
    do ie=1, nenergy
        ebarr(ie) = 10.d0**(log10(emin) + (ie-1)*dlogE)
    enddo

    dlogT = (log10(tmax) - log10(tmin))/(ntemp - 1)
    do it=1, ntemp
        tarr(it) = 10.d0**(log10(tmin) + (it-1)*dlogT)
    enddo
    write(*,'(a)') "---- H-"//trim(adjustl(aname))//" reaction rates settings ----"
    write(*,'(T2,"q = ", i2)') q
    write(*,'(T2,"mass = ",f7.2, " amu")') mass
    write(*,'(T2,"Emin = ",e9.2, " keV")') emin
    write(*,'(T2,"Emax = ",e9.2, " keV")') emax
    write(*,'(T2,"Nenergy = ", i4)') nenergy
    write(*,'(T2,"Tmin = ",e9.2, " keV")') tmin
    write(*,'(T2,"Tmax = ",e9.2, " keV")') tmax
    write(*,'(T2,"Ntemp = ", i4)') ntemp
    write(*,*) ''

    cnt = 0
    !$OMP PARALLEL DO private(ie, it, ia, n, m, eb, ti, rate)
    do ie=1, nenergy
        eb = ebarr(ie)
        do it=1, ntemp
            ti = tarr(it)
            do ia=1, n_bt_amu
                do n=1, n_max
                    do m=1, m_max
                        if(m.gt.n) then
                            call bt_maxwellian(Aq_excit_n_m, q, ti, eb, &
                                               a(2,ia), a(1,ia), n, m, rate)
                            excit(n,m,ie,it,ia) = rate

                            call bt_maxwellian(Aq_excit_n_m, q, ti, eb, &
                                               a(2,ia), a(1,ia), n, m, &
                                               rate, deexcit=.True.)
                            excit(m,n,ie,it,ia) = rate
                        endif
                    enddo
                    call bt_maxwellian(Aq_cx_n, q, ti, eb, &
                                       a(2,ia), a(1,ia), n, rate)
                    cx(n,ie,it,ia) = rate

                    call bt_maxwellian(Aq_ioniz_n, q, ti, eb, &
                                       a(2,ia), a(1,ia), n, rate)
                    ioniz(n,ie,it,ia) = rate
                enddo
            enddo
            cnt = cnt + 1
            WRITE(*,'(f7.2,"%",a,$)') 100*cnt/real(nenergy*ntemp),char(13)
        enddo
    enddo
    !$OMP END PARALLEL DO

    call h5gcreate_f(id, trim(adjustl(asym)), gid, error)

    dim1 = [1]
    dim4 = [n_max, nenergy, ntemp, n_bt_amu]
    dim5 = [n_max, m_max, nenergy, ntemp, n_bt_amu]
    call h5ltmake_dataset_int_f(gid, "nenergy", 0, dim1, [nenergy], error)
    call h5ltmake_dataset_int_f(gid, "ntemp", 0, dim1, [ntemp], error)
    call h5ltmake_dataset_int_f(gid, "n_bt_amu", 0, dim1, [n_bt_amu], 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)
    call h5ltmake_dataset_double_f(gid, "dlogT", 0, dim1, [dlogT], error)
    call h5ltmake_dataset_double_f(gid, "tmin", 0, dim1, [tmin], error)
    call h5ltmake_dataset_double_f(gid, "tmax", 0, dim1, [tmax], error)

    dim2 = [2,n_bt_amu]
    call h5ltmake_compressed_dataset_double_f(gid, "bt_amu", 2, dim2, a, error)
    dim1 = [nenergy]
    call h5ltmake_compressed_dataset_double_f(gid, "energy", 1, dim1, ebarr, error)
    dim1 = [ntemp]
    call h5ltmake_compressed_dataset_double_f(gid, "temperature", 1, dim1, tarr, error)
    call h5ltmake_compressed_dataset_double_f(gid, "cx", 4, dim4, cx, error)
    call h5ltmake_compressed_dataset_double_f(gid, "ionization", 4, dim4, ioniz, error)
    call h5ltmake_compressed_dataset_double_f(gid, "excitation", 5, dim5, excit, error)

    call h5ltset_attribute_string_f(id, trim(adjustl(asym)), "description", &
         "Beam-Target reaction rates for Hydrogen(beam)-"//trim(adjustl(aname))// &
         "(target) interactions", error)
    call h5ltset_attribute_string_f(gid, "nenergy", "description", &
         "Number of energy values", error)
    call h5ltset_attribute_string_f(gid, "ntemp", "description", &
         "Number of target temperature values", error)
    call h5ltset_attribute_string_f(gid, "n_bt_amu", "description", &
         "Number of beam-target amu combinations", 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, "bt_amu", "description", &
         "Combinations of beam-target amu's e.g. b_amu, t_amu = bt_amu[:,i]", error)
    call h5ltset_attribute_string_f(gid, "energy", "description", &
         "Energy values", error)
    call h5ltset_attribute_string_f(gid, "energy", "units", "keV", 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)", error)
    call h5ltset_attribute_string_f(gid, "emin","description", &
         "Minimum energy", error)
    call h5ltset_attribute_string_f(gid, "emin", "units", "keV", error)
    call h5ltset_attribute_string_f(gid, "emax","description", &
         "Maximum energy", error)
    call h5ltset_attribute_string_f(gid, "emax", "units", "keV", error)

    call h5ltset_attribute_string_f(gid, "temperature", "description", &
         "Target temperature values", error)
    call h5ltset_attribute_string_f(gid, "temperature", "units", "keV", error)
    call h5ltset_attribute_string_f(gid, "dlogT", "description", &
         "Temperature spacing in log-10", error)
    call h5ltset_attribute_string_f(gid, "dlogT", "units", "log10(keV)", error)
    call h5ltset_attribute_string_f(gid, "tmin","description", &
         "Minimum temperature", error)
    call h5ltset_attribute_string_f(gid, "tmin", "units", "keV", error)
    call h5ltset_attribute_string_f(gid, "tmax","description", &
         "Maximum temperature", error)
    call h5ltset_attribute_string_f(gid, "tmax", "units", "keV", error)

    call h5ltset_attribute_string_f(gid, "cx", "description", &
         "n-resolved charge exchange reaction rates: cx(n,energy,temp,bt_amu)", error)
    call h5ltset_attribute_string_f(gid, "cx", "units", "cm^3/s", error)
    call h5ltset_attribute_string_f(gid, "cx", "reaction", &
         "A(q+) + H(n) -> A((q-1)+) + H(+)", error)

    call h5ltset_attribute_string_f(gid, "excitation", "description", &
         "n/m resolved (de-)excitation reaction rates: excitation(n,m,energy,temp,bt_amu)", error)
    call h5ltset_attribute_string_f(gid, "excitation", "units", "cm^3/s", error)
    call h5ltset_attribute_string_f(gid, "excitation", "reaction", &
         "A(q+) + H(n) -> A(q+) + H(m); m > n excitation, m < n de-excitation", error)

    call h5ltset_attribute_string_f(gid, "ionization", "description", &
         "n resolved ionization reaction rates: ionization(n,energy,temp,bt_amu)", error)
    call h5ltset_attribute_string_f(gid, "ionization", "units", "cm^3/s", error)
    call h5ltset_attribute_string_f(gid, "ionization", "reaction", &
         "A(q+) + H(n) -> A(q+) + H(+) + e-", error)

    call h5gclose_f(gid, error)

    deallocate(ebarr, tarr, excit, ioniz)

end subroutine write_bt_H_Aq