-
Notifications
You must be signed in to change notification settings - Fork 5
/
radiation_sizes.f90
91 lines (60 loc) · 3.02 KB
/
radiation_sizes.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
! This module sets the sizes of radiation quantities
module radiation_sizes
use precision, only: dp
use cgsphotoconstants, only: ion_freq_HI,& ! HI ionization energy in frequency
sigma_HI_at_ion_freq ! HI cross section at its ionzing frequency
use c2ray_parameters, only: isothermal
implicit none
integer,parameter :: NumFreq = 128 ! Number of integration points in each of the frequency bins
integer,parameter :: NumTau = 2000 ! Number of table points for the optical depth
integer,parameter :: NumBndin1 = 1 ! Number of frequency sub-bins in interval 1
integer,parameter :: NumFreqBnd=NumBndin1 ! Total number of frequency bins
integer,parameter :: NumheatBin=NumBndin1 ! Total number of heating bins
! Optical depths at the entrance of the grid.
! It can be used if radiation enters the simulation volume from the outside.
real(kind=dp) :: boundary_tauHI = 0.0
real(kind=dp), dimension(:), allocatable :: delta_freq ! Frequency width of integration
real(kind=dp), dimension(:), allocatable :: freq_max ! Maximum freqeucny of integration
real(kind=dp), dimension(:), allocatable :: freq_min ! Minimum freqeucny of integration
! Power law fit parameter for frequency range 1:3
real(kind=dp), dimension(:), allocatable :: pl_index_cross_section_HI ! Power law index of cross section of HI
! Cross section of atoms
real(kind=dp), dimension(:), allocatable :: sigma_HI ! Cross section of HI
contains
! set up some scaling factors arrays
subroutine setup_scalingfactors (freq_min_src, freq_max_src)
! Contains the maximum frequency, determined by source input
real(kind=dp),intent(in) :: freq_min_src, freq_max_src
integer :: i_subband
! Allocate size of arrays
allocate(delta_freq(1:NumFreqBnd))
allocate(freq_max(1:NumFreqBnd))
allocate(freq_min(1:NumFreqBnd))
allocate(pl_index_cross_section_HI(1:NumFreqBnd))
allocate(sigma_HI(1:NumFreqBnd))
! Assignment of maximum frequency in the sub-bin partition.
select case (NumBndin1)
case (1)
freq_max(NumBndin1) = freq_max_src
end select
! Assignment of minimum frequency in the sub-bin partition.
freq_min(NumBndin1) = max(ion_freq_HI, freq_min_src)
! calculate the width of frequency sub-bin
do i_subband=1,NumFreqBnd
delta_freq(i_subband) = (freq_max(i_subband)-freq_min(i_subband))/real(NumFreq)
enddo
! Assign value sigma of HI, HeI and HeII at different frequencies
! This value is used to characterise how absorbing the medium is and thus
! is independent of the SED
select case (NumBndin1)
case (1)
sigma_HI(1) = sigma_HI_at_ion_freq
end select
! Assign power-law index of HI, HeI and HeII at different frequencies (about absorption)
select case (NumBndin1)
case (1)
!pl_index_cross_section_HI(1) = 2.761_dp
pl_index_cross_section_HI(1) = 2.8_dp
end select
end subroutine setup_scalingfactors
end module radiation_sizes