Line data Source code
1 : ! ***********************************************************************
2 : !
3 : ! Copyright (C) 2011-2020 Bill Paxton & The MESA Team
4 : !
5 : ! This program is free software: you can redistribute it and/or modify
6 : ! it under the terms of the GNU Lesser General Public License
7 : ! as published by the Free Software Foundation,
8 : ! either version 3 of the License, or (at your option) any later version.
9 : !
10 : ! This program is distributed in the hope that it will be useful,
11 : ! but WITHOUT ANY WARRANTY; without even the implied warranty of
12 : ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
13 : ! See the GNU Lesser General Public License for more details.
14 : !
15 : ! You should have received a copy of the GNU Lesser General Public License
16 : ! along with this program. If not, see <https://www.gnu.org/licenses/>.
17 : !
18 : ! ***********************************************************************
19 :
20 : module test_screen
21 :
22 : use const_def, only: dp
23 : use rates_lib
24 : use rates_def
25 : use utils_lib, only: mesa_error
26 :
27 : implicit none
28 :
29 : contains
30 :
31 2 : subroutine do_test_screen
32 : use chem_def
33 : use chem_lib
34 : use const_lib, only: const_init
35 : use math_lib
36 :
37 : integer, parameter :: num_isos = 8, max_z_to_cache = 12
38 : integer :: chem_id(num_isos), i1, i2, ierr
39 2 : integer, pointer :: net_iso(:)
40 52 : real(dp) :: xin(num_isos), y(num_isos), iso_z(num_isos), xz, abar, zbar, z2bar, z53bar, ye, sumx, &
41 36 : dabar_dx(num_isos), dzbar_dx(num_isos), temp, den, logT, logRho, &
42 36 : sc1a, sc1adt, sc1add, xh, xhe, dmc_dx(num_isos), iso_z158(num_isos)
43 : type(Screen_Info) :: sc
44 2 : real(dp) :: zs13, zhat, zhat2, lzav, aznut, zs13inv, mass_correction !approx_abar, approx_zbar
45 : integer :: i
46 : integer :: h1, he3, he4, c12, n14, o16, ne20, mg24
47 : character(len=32) :: my_mesa_dir
48 :
49 : include 'formats'
50 :
51 2 : ierr = 0
52 2 : my_mesa_dir = '../..'
53 2 : call const_init(my_mesa_dir, ierr)
54 2 : if (ierr /= 0) then
55 0 : write (*, *) 'const_init failed'
56 0 : call mesa_error(__FILE__, __LINE__)
57 : end if
58 2 : call math_init()
59 :
60 2 : call chem_init('isotopes.data', ierr)
61 2 : if (ierr /= 0) then
62 0 : write (*, *) 'chem_init failed'
63 0 : call mesa_error(__FILE__, __LINE__)
64 : end if
65 :
66 2 : h1 = 1
67 2 : he3 = 2
68 2 : he4 = 3
69 2 : c12 = 4
70 2 : n14 = 5
71 2 : o16 = 6
72 2 : ne20 = 7
73 2 : mg24 = 8
74 :
75 2 : allocate (net_iso(num_chem_isos))
76 :
77 15714 : net_iso = 0
78 :
79 2 : net_iso(ih1) = h1; chem_id(h1) = ih1
80 2 : net_iso(ihe3) = he3; chem_id(he3) = ihe3
81 2 : net_iso(ihe4) = he4; chem_id(he4) = ihe4
82 2 : net_iso(ic12) = c12; chem_id(c12) = ic12
83 2 : net_iso(in14) = n14; chem_id(n14) = in14
84 2 : net_iso(io16) = o16; chem_id(o16) = io16
85 2 : net_iso(ine20) = ne20; chem_id(ne20) = ine20
86 2 : net_iso(img24) = mg24; chem_id(mg24) = img24
87 :
88 2 : logT = 7.7110722845770692D+00
89 2 : logRho = 4.5306372623742392D+00
90 :
91 2 : xin(net_iso(ih1)) = 9.1649493293186402D-22
92 2 : xin(net_iso(ihe3)) = 1.8583723770506327D-24
93 2 : xin(net_iso(ihe4)) = 9.8958688392029037D-01
94 2 : xin(net_iso(ic12)) = 7.3034226840994307D-04
95 2 : xin(net_iso(in14)) = 6.2231918940828723D-03
96 2 : xin(net_iso(io16)) = 3.6335428720509150D-04
97 2 : xin(net_iso(ine20)) = 1.0527812894109640D-03
98 2 : xin(net_iso(img24)) = 2.0434463406007555D-03
99 :
100 2 : i1 = ihe4
101 2 : i2 = ic12
102 :
103 : if (.false.) then ! TESTING
104 :
105 : xin = 0
106 : xin(net_iso(ih1)) = 0.72d0
107 : xin(net_iso(ihe4)) = 0.26d0
108 : xin(net_iso(in14)) = 0.02d0
109 :
110 : i1 = ih1
111 : i2 = in14
112 :
113 : write (*, 1) 'sum(xin)', sum(xin(:))
114 :
115 : logT = 7d0
116 : logRho = 1d0
117 :
118 : end if
119 :
120 : call composition_info( &
121 : num_isos, chem_id, xin, xh, xhe, xz, abar, zbar, z2bar, z53bar, &
122 2 : ye, mass_correction, sumx, dabar_dx, dzbar_dx, dmc_dx)
123 :
124 18 : iso_z(:) = chem_isos%Z(chem_id(:))
125 :
126 18 : do i = 1, num_isos
127 18 : iso_z158(i) = pow(real(chem_isos%Z(chem_id(i)), kind=dp), 1.58d0)
128 : end do
129 18 : y(:) = xin(:)/chem_isos%Z_plus_N(chem_id(:))
130 :
131 2 : call do1(salpeter_screening, ierr)
132 2 : if (ierr /= 0) call mesa_error(__FILE__, __LINE__)
133 :
134 2 : call do1(extended_screening, ierr)
135 2 : if (ierr /= 0) call mesa_error(__FILE__, __LINE__)
136 :
137 2 : call do1(chugunov_screening, ierr)
138 2 : if (ierr /= 0) call mesa_error(__FILE__, __LINE__)
139 :
140 2 : call do1(no_screening, ierr)
141 2 : if (ierr /= 0) call mesa_error(__FILE__, __LINE__)
142 :
143 2 : deallocate (net_iso)
144 :
145 8 : write (*, *) 'done'
146 :
147 : contains
148 :
149 8 : subroutine do1(sc_mode, ierr)
150 2 : use math_lib
151 : integer, intent(in) :: sc_mode
152 : integer, intent(out) :: ierr
153 : character(len=64) :: sc_str
154 : include 'formats'
155 8 : call screening_option_str(sc_mode, sc_str, ierr)
156 8 : if (ierr /= 0) return
157 8 : write (*, *) trim(sc_str)
158 :
159 8 : temp = exp10(logT)
160 8 : den = exp10(logRho)
161 :
162 : call screen_init_AZ_info( &
163 0 : chem_isos%W(i1), dble(chem_isos%Z(i1)), &
164 0 : chem_isos%W(i2), dble(chem_isos%Z(i2)), &
165 : zs13, &
166 : zhat, zhat2, lzav, aznut, zs13inv, &
167 8 : ierr)
168 8 : if (ierr /= 0) return
169 :
170 : ! set the context for the reactions
171 : call screen_set_context( &
172 : sc, temp, den, logT, logRho, zbar, abar, z2bar, &
173 8 : sc_mode, num_isos, y, iso_z158)
174 :
175 : call screen_pair( &
176 : sc, &
177 0 : chem_isos%W(i1), dble(chem_isos%Z(i1)), &
178 0 : chem_isos%W(i2), dble(chem_isos%Z(i2)), &
179 : sc_mode, &
180 : zs13, zhat, zhat2, lzav, aznut, zs13inv, rattab_tlo, &
181 8 : sc1a, sc1adt, sc1add, ierr)
182 8 : if (ierr /= 0) return
183 :
184 8 : write (*, 1) 'logT = ', logT
185 8 : write (*, 1) 'logRho = ', logRho
186 8 : write (*, 1) 'zbar = ', zbar
187 8 : write (*, 1) 'abar = ', abar
188 8 : write (*, 1) 'z2bar = ', z2bar
189 8 : write (*, 1) 'sc1a = ', sc1a
190 8 : write (*, 1) 'sc1adt = ', sc1adt
191 8 : write (*, 1) 'sc1add = ', sc1add
192 8 : write (*, '(A)')
193 :
194 8 : end subroutine do1
195 :
196 : end subroutine do_test_screen
197 :
198 : end module test_screen
199 :
|