1+ # This file is part of AGNI. License is GPL-3.0: https://www.gnu.org/licenses
2+
3+ """
4+ **Run FastChem to solve gas phase chemistry.**
5+
6+ FastChem is a fast and robust chemical equilibrium solver, designed for use in exoplanet
7+ and brown dwarf atmospheres. This module writes the necessary input files for FastChem,
8+ runs the code, and parses the output to update the gas phase composition.
9+ """
110module fastchem
211
312 # Import packages
@@ -6,8 +15,10 @@ module fastchem
615 import DelimitedFiles: readdlm
716
817 # Include local modules
9- using .. atmosphere
10- using .. phys
18+ import .. atmosphere
19+ import .. phys
20+ import .. consts
21+ import .. formulae
1122
1223 # Constants
1324 const SMOOTH_SCALE:: Float64 = 12.0 # smoothing scale for fastchem floor
@@ -93,7 +104,7 @@ module fastchem
93104 @debug " Elements set by user-provided metallicities"
94105
95106 # copy original to calculated; set elem to zero if it was not provided
96- for e in phys . elems_standard
107+ for e in consts . elems_standard
97108 atmos. metal_calc[e] = get (atmos. metal_orig, e, 0.0 )
98109 end
99110 atmos. metal_calc[" H" ] = 1.0
@@ -104,8 +115,8 @@ module fastchem
104115
105116 # Calculate elemental abundances from surface mixing ratios [molecules/m^3]
106117 # assuming ideal gas: N/V = P*x/(Kb*T) , where x is the VMR
107- N_inp_t = zeros (Float64, length (phys . elems_standard)) # total atoms in all gases
108- N_inp_g = zeros (Float64, length (phys . elems_standard)) # atoms in current gas
118+ N_inp_t = zeros (Float64, length (consts . elems_standard)) # total atoms in all gases
119+ N_inp_g = zeros (Float64, length (consts . elems_standard)) # atoms in current gas
109120 # loop over gases
110121 for gas in atmos. gas_names
111122 fill! (N_inp_g, 0.0 )
@@ -116,16 +127,16 @@ module fastchem
116127 end
117128
118129 # count atoms in this gas
119- d = phys . count_atoms (gas)
120- for (i,e) in enumerate (phys . elems_standard)
130+ d = formulae . count_atoms (gas)
131+ for (i,e) in enumerate (consts . elems_standard)
121132 if haskey (d, e)
122133 N_inp_g[i] += d[e] # N_inp_g stores num of atoms in this gas
123134 end
124135 end
125136
126137 # Get gas abundance from original VMR value
127138 # scale number of atoms by the abundance of the gas (p = Ng kB T)
128- N_inp_g *= atmos. gas_vmr[gas][end ] * atmos. p[end ] / (phys . k_B * atmos. tmp[end ])
139+ N_inp_g *= atmos. gas_vmr[gas][end ] * atmos. p[end ] / (consts . k_B * atmos. tmp[end ])
129140
130141 # Add atoms from this gas to total atoms in the mixture
131142 N_inp_t += N_inp_g
@@ -138,7 +149,7 @@ module fastchem
138149 end
139150
140151 # Convert elemental abundances to metallicity number ratios, rel to hydrogen
141- for (i,e) in enumerate (phys . elems_standard)
152+ for (i,e) in enumerate (consts . elems_standard)
142153 atmos. metal_calc[e] = N_inp_t[i]/ N_inp_t[1 ]
143154 end
144155 end
@@ -148,7 +159,7 @@ module fastchem
148159 # for each element `e`, value = log10(N_e/N_H) + 12
149160 open (atmos. fastchem_elem," w" ) do f
150161 write (f," # Elemental abundances file written by AGNI \n " )
151- for e in phys . elems_standard
162+ for e in consts . elems_standard
152163
153164 # skip this element if its abundance is too small
154165 if atmos. metal_calc[e] < 1.0e-30
@@ -160,7 +171,7 @@ module fastchem
160171 write (f, @sprintf (" %s %.3f \n " ,e,log10 (atmos. metal_calc[e]) + 12.0 ))
161172 else
162173 @warn " Got non-finite metallicity for $e - adopting solar value"
163- write (f, @sprintf (" %s %.3f \n " ,e,phys . consts. _solar_metallicity [e]))
174+ write (f, @sprintf (" %s %.3f \n " ,e,consts. solar_metallicity [e]))
164175 end
165176
166177 count_elem_nonzero += 1
@@ -276,10 +287,10 @@ module fastchem
276287
277288 # not stored => search based on atoms
278289 if ! match
279- d_fc = phys . count_atoms (g_fc) # get atoms dict from FC name
290+ d_fc = formulae . count_atoms (g_fc) # get atoms dict from FC name
280291
281292 for g in atmos. gas_names
282- if phys . same_atoms (d_fc, atmos. gas_dat[g]. atoms)
293+ if formulae . same_atoms (d_fc, atmos. gas_dat[g]. atoms)
283294 match = true
284295 g_in = g
285296 break
0 commit comments