Line data Source code
1 : module mod_test_chem
2 :
3 : use chem_lib
4 : use chem_def
5 : use const_def, only: dp
6 : use utils_lib, only: mesa_error
7 :
8 : implicit none
9 :
10 : contains
11 :
12 1 : subroutine do_test
13 : use const_lib, only: const_init
14 : use math_lib, only: math_init
15 : integer :: ierr, i
16 : character(len=32) :: my_mesa_dir
17 1 : my_mesa_dir = '../..'
18 1 : call const_init(my_mesa_dir, ierr)
19 1 : if (ierr /= 0) then
20 0 : write (*, *) 'const_init failed'
21 0 : call mesa_error(__FILE__, __LINE__)
22 : end if
23 1 : call math_init()
24 1 : call chem_init('isotopes.data', ierr)
25 1 : if (ierr /= 0) then
26 0 : write (*, *) 'FATAL ERROR: failed in chem_init'
27 0 : call mesa_error(__FILE__, __LINE__)
28 : end if
29 1 : call do_test_lodders
30 1 : call do_tests
31 51 : do i = 1, 50
32 50 : call do_tests
33 51 : write (*, *) 'done i', i
34 : end do
35 1 : end subroutine do_test
36 :
37 51 : subroutine do_tests
38 :
39 51 : write (*, *) 'call do_test_chem'
40 51 : flush (6)
41 51 : call do_test_chem
42 :
43 51 : write (*, *) 'call do_test_composition_info'
44 51 : flush (6)
45 51 : call do_test_composition_info
46 :
47 51 : write (*, *) 'call do_test_Qtotal'
48 51 : flush (6)
49 51 : call do_test_Qtotal
50 :
51 51 : write (*, *) 'done do_tests'
52 51 : flush (6)
53 :
54 1 : end subroutine do_tests
55 :
56 51 : subroutine do_test_chem
57 51 : real(dp) :: c, n, o, cno
58 : integer :: ic12, in14, io16
59 : include 'formats'
60 :
61 51 : ic12 = get_nuclide_index('c12')
62 51 : in14 = get_nuclide_index('n14')
63 51 : io16 = get_nuclide_index('o16')
64 :
65 51 : write (*, *)
66 51 : write (*, *)
67 51 : write (*, 1) 'chem_W(io16)', chem_isos%W(io16)
68 51 : flush (6)
69 51 : write (*, 1) 'chem_Z(io16)', dble(chem_isos%Z(io16))
70 51 : flush (6)
71 51 : write (*, 1) 'chem_binding_energy(io16)', chem_isos%binding_energy(io16)
72 51 : flush (6)
73 51 : write (*, *)
74 51 : write (*, *) 'chem_name(io16) ', chem_isos%name(io16)
75 51 : flush (6)
76 51 : write (*, *) 'get_nuclide_index("o16") == io16', get_nuclide_index("o16") == io16
77 51 : flush (6)
78 51 : write (*, *)
79 51 : write (*, *) 'chem_element_Name(e_he) ', chem_element_Name(e_he)
80 51 : flush (6)
81 51 : write (*, *) 'chem_get_element_id("he") == e_he', chem_get_element_id("he") == e_he
82 51 : flush (6)
83 51 : write (*, *)
84 51 : write (*, *)
85 51 : write (*, 1) 'Anders & Grevesse 1989 zsol', zsol
86 51 : write (*, 1) 'Anders & Grevesse 1989 yesol', yesol
87 51 : write (*, *)
88 51 : write (*, *)
89 51 : write (*, *) 'cno fraction by mass of Z'
90 51 : write (*, *)
91 51 : flush (6)
92 :
93 51 : c = chem_Xsol('c12')
94 51 : n = chem_Xsol('n14')
95 51 : o = chem_Xsol('o16')
96 : !write(*,1) 'c', c
97 : !write(*,1) 'n', n
98 : !write(*,1) 'o', o
99 51 : cno = c + n + o
100 51 : write (*, 1) 'Anders & Grevesse 1989', cno/zsol
101 51 : write (*, *)
102 51 : cno = GN93_element_zfrac(e_c) + GN93_element_zfrac(e_n) + GN93_element_zfrac(e_o)
103 51 : write (*, 1) 'Grevesse and Noels 1993', cno
104 51 : write (*, *)
105 51 : cno = GS98_element_zfrac(e_c) + GS98_element_zfrac(e_n) + GS98_element_zfrac(e_o)
106 51 : write (*, 1) 'Grevesse and Sauval 1998', cno
107 51 : write (*, *)
108 51 : cno = L03_element_zfrac(e_c) + L03_element_zfrac(e_n) + L03_element_zfrac(e_o)
109 51 : write (*, 1) 'Lodders 2003', cno
110 51 : write (*, *)
111 51 : cno = AGS05_element_zfrac(e_c) + AGS05_element_zfrac(e_n) + AGS05_element_zfrac(e_o)
112 51 : write (*, 1) 'Asplund, Grevesse and Sauval 2005', cno
113 51 : write (*, *)
114 51 : return
115 :
116 : write (*, 1) 'Grevesse and Sauval 1998 C', GS98_element_zfrac(e_c)
117 : write (*, *)
118 : write (*, 1) 'Grevesse and Sauval 1998 N', GS98_element_zfrac(e_n)
119 : write (*, *)
120 : write (*, 1) 'Grevesse and Sauval 1998 O', GS98_element_zfrac(e_o)
121 : write (*, *)
122 : write (*, 1) 'Grevesse and Sauval 1998 Ne', GS98_element_zfrac(e_ne)
123 : write (*, *)
124 : write (*, 1) 'Grevesse and Sauval 1998 Mg', GS98_element_zfrac(e_mg)
125 : write (*, *)
126 : write (*, 1) 'Grevesse and Sauval 1998 Si', GS98_element_zfrac(e_si)
127 : write (*, *)
128 : write (*, 1) 'Grevesse and Sauval 1998 Fe', GS98_element_zfrac(e_fe)
129 : write (*, *)
130 : return
131 :
132 : stop
133 :
134 : end subroutine do_test_chem
135 :
136 51 : subroutine do_test_composition_info
137 : integer, parameter :: num_species = 2
138 : integer, dimension(num_species) :: ids
139 153 : real(dp), dimension(num_species) :: X
140 51 : real(dp) :: xh, xhe, xz, abar, zbar, z2bar, z53bar, ye, mass_correction, sumx
141 357 : real(dp) :: dabar_dx(num_species), dzbar_dx(num_species), dmc_dx(num_species)
142 : character(len=*), parameter :: form1 = '(a,t12,"=",f11.6)'
143 51 : ids(1) = get_nuclide_index('c12')
144 51 : ids(2) = get_nuclide_index('fe56')
145 153 : X = 0.5d0
146 : call composition_info(num_species, ids, X, xh, xhe, xz, &
147 51 : abar, zbar, z2bar, z53bar, ye, mass_correction, sumx, dabar_dx, dzbar_dx, dmc_dx)
148 51 : write (*, *) 'test moments of composition and the mass correction factor'
149 51 : write (*, *) 'for a C12-Fe56 mixture (both 50% by mass)'
150 51 : write (*, *)
151 51 : write (*, form1) 'Abar', abar
152 51 : write (*, form1) 'Zbar', zbar
153 51 : write (*, form1) 'Z2bar', z2bar
154 51 : write (*, form1) 'Ye', ye
155 51 : write (*, form1) 'sum(X*W/A)', mass_correction
156 51 : end subroutine do_test_composition_info
157 :
158 51 : subroutine do_test_Qtotal
159 :
160 : integer, parameter :: num_in = 2, num_out = 3
161 : integer :: reactants(num_in + num_out)
162 51 : real(dp) :: Qtotal
163 :
164 : include 'formats'
165 :
166 51 : write (*, *) 'test Qtotal for 2 he3 => 2 h1 + he4'
167 51 : write (*, *)
168 :
169 51 : reactants(1) = ihe3
170 51 : reactants(2) = ihe3
171 51 : reactants(3) = ihe4
172 51 : reactants(4) = ih1
173 51 : reactants(5) = ih1
174 51 : Qtotal = reaction_Qtotal(num_in, num_out, reactants, chem_isos)
175 :
176 51 : write (*, 1) 'Qtotal', Qtotal ! expect 12.86
177 51 : write (*, *)
178 :
179 51 : end subroutine do_test_Qtotal
180 :
181 1 : subroutine do_test_lodders
182 : integer :: i
183 1 : real(dp) :: percent
184 1 : write (*, *)
185 1 : write (*, '(a,/,72("="))') 'output of solar abundances: compare with Lodders (2003) table'
186 1 : write (*, '(a7,tr3,a11)') 'isotope', '% abundance'
187 287 : do i = 1, size(namsol)
188 286 : percent = lodders03_element_atom_percent(namsol(i))
189 287 : write (*, '(a7,tr3,f11.6)') namsol(i), percent
190 : end do
191 1 : end subroutine do_test_lodders
192 :
193 0 : subroutine write_chem_ids_file
194 : integer, parameter :: iounit = 33
195 : integer :: ierr
196 0 : ierr = 0
197 0 : open (unit=iounit, file=trim('chem_ids.list'), action='write', status='replace', iostat=ierr)
198 0 : if (ierr /= 0) then
199 0 : write (*, *) 'failed to open file for write_chem_ids'
200 0 : call mesa_error(__FILE__, __LINE__)
201 : end if
202 0 : call write_chem_ids(iounit)
203 0 : close (iounit)
204 0 : end subroutine write_chem_ids_file
205 :
206 0 : subroutine write_chem_ids(iounit)
207 : integer, intent(in) :: iounit
208 : integer :: i
209 0 : do i = 1, num_chem_isos
210 0 : write (iounit, '(5x,i5,3x,a5)') i, chem_isos%name(i)
211 : end do
212 0 : write (iounit, *)
213 0 : end subroutine write_chem_ids
214 :
215 0 : subroutine write_element_ids
216 : integer :: i
217 0 : do i = 0, max_el_z
218 0 : write (*, *) el_name(i), el_long_name(i)
219 : end do
220 0 : end subroutine write_element_ids
221 :
222 : end module mod_test_chem
223 :
224 1 : program test_chem
225 1 : use mod_test_chem
226 : implicit none
227 : !call write_element_ids; stop
228 1 : call do_test
229 1 : end program test_chem
|