cloudy trunk
Loading...
Searching...
No Matches
init_coreload.cpp
Go to the documentation of this file.
1/* This file is part of Cloudy and is copyright (C)1978-2013 by Gary J. Ferland and
2 * others. For conditions of distribution and use see copyright notice in license.txt */
3/*InitCoreload one time initialization of core load, called from cdDrive, this sets
4* minimum set of values needed for the code to start - called after
5* input lines have been read in and checked for VARY or GRID - so
6* known whether single or multiple sims will be run */
7#include "cddefines.h"
8#include "optimize.h"
9#include "parse.h"
10#include "struc.h"
11#include "input.h"
12#include "dense.h"
13#include "hcmap.h"
14#include "h2.h"
15#include "version.h"
16#include "h2.h"
17#include "hextra.h"
18#include "heavy.h"
19#include "grid.h"
20#include "ionbal.h"
21#include "iso.h"
22#include "taulines.h"
23#include "cosmology.h"
24#include "physconst.h"
25#include "broke.h"
26#include "rfield.h"
27
28// //////////////////////////////////////////////////////////////////////////
29//
30//
31// NB DO NOT ADD VARIABLES TO THIS FILE! THE GOAL IS TO REMOVE THIS FILE
32// initialization of variables done one time per coreload should be done in
33// a constructor for the data
34//
35//
36// //////////////////////////////////////////////////////////////////////////
37
38/*InitCoreload one time initialization of core load, called from cdDrive, this sets
39* minimum set of values needed for the code to start - called after
40* input lines have been read in and checked for VARY or GRID - so
41* known whether single or multiple sims will be run */
42void InitCoreload( void )
43{
44 static int nCalled=0;
45 long int nelem;
46
47 DEBUG_ENTRY( "InitCoreload()" );
48
49 /* return if already called */
50 if( nCalled )
51 return;
52
53 ++nCalled;
54
55 rfield.lgMeshSetUp = false;
56
57 hcmap.lgMapOK = true;
58 hcmap.lgMapDone = false;
59
60 // subdirectory delimiter character
61# ifdef _WIN32
62 strcpy( input.chDelimiter, "\\" );
63# else
64 strcpy( input.chDelimiter, "/" );
65# endif
66
67 /* will be reset to positive number when map actually done */
68 hcmap.nMapAlloc = 0;
69 hcmap.nmap = 0;
70 hcmap.lgMapBeingDone = false;
71
72 /* name of output file from optimization run */
73 strncpy( chOptimFileName , "optimal.in" , sizeof( chOptimFileName ) );
74
75 /* number of electrons in valence shell that can Compton recoil ionize */
76 long int nCom[LIMELM] =
77 {
78 1 , 2 , /* K 1s shell */
79 1 , 2 , /* L 2s shell */
80 1 , 2 , 3 , 4 , 5 , 6 , /* L 2p shell */
81 1 , 2 , /* M 3s shell */
82 1 , 2 , 3 , 4 , 5 , 6 , /* M 3p shell */
83 1 , 2 , 3 , 4 , 5 , 6 , 7 , 8 , /* M 3d shell */
84 1 , 2 , /* N 4s shell */
85 1 , 2 /* N 4p shell */
86 };
87
88 for( nelem=ipHYDROGEN; nelem<LIMELM; ++nelem )
89 {
90 ionbal.nCompRecoilElec[nelem] = nCom[nelem];
91 }
92
93 /* list of shells, 1 to 7 */
94 strcpy( Heavy.chShell[0], "1s" );
95 strcpy( Heavy.chShell[1], "2s" );
96 strcpy( Heavy.chShell[2], "2p" );
97 strcpy( Heavy.chShell[3], "3s" );
98 strcpy( Heavy.chShell[4], "3p" );
99 strcpy( Heavy.chShell[5], "3d" );
100 strcpy( Heavy.chShell[6], "4s" );
101
102 /* variables for H-like sequence */
103 /* default number of levels for hydrogen iso sequence */
104 for( nelem=ipHYDROGEN; nelem < LIMELM; ++nelem )
105 {
106 iso_sp[ipH_LIKE][nelem].n_HighestResolved_max = 5;
107 iso_sp[ipH_LIKE][nelem].nCollapsed_max = 2;
108 }
109
110 iso_sp[ipH_LIKE][ipHYDROGEN].n_HighestResolved_max = 10;
111 iso_sp[ipH_LIKE][ipHELIUM].n_HighestResolved_max = 10;
112
113 iso_sp[ipH_LIKE][ipHYDROGEN].nCollapsed_max = 15;
114 iso_sp[ipH_LIKE][ipHELIUM].nCollapsed_max = 15;
115
116 /* variables for He-like sequence */
117 /* "he-like" hydrogen (H-) is treated elsewhere */
118 iso_sp[ipHE_LIKE][ipHYDROGEN].n_HighestResolved_max = -LONG_MAX;
119 iso_sp[ipHE_LIKE][ipHYDROGEN].numLevels_max = -LONG_MAX;
120 iso_sp[ipHE_LIKE][ipHYDROGEN].nCollapsed_max = -LONG_MAX;
121
122 for( nelem=ipHELIUM; nelem < LIMELM; ++nelem )
123 {
124 /* put at least three resolved and 1 collapsed in every element for he-like */
125 iso_sp[ipHE_LIKE][nelem].n_HighestResolved_max = 3;
126 iso_sp[ipHE_LIKE][nelem].nCollapsed_max = 1;
127 }
128
129 iso_sp[ipHE_LIKE][ipHELIUM].n_HighestResolved_max = 6;
130 iso_sp[ipHE_LIKE][ipHELIUM].nCollapsed_max = 20;
131
132 /* And n=5 for these because they are most abundant */
133 iso_sp[ipHE_LIKE][ipCARBON].n_HighestResolved_max = 5;
134 iso_sp[ipHE_LIKE][ipNITROGEN].n_HighestResolved_max = 5;
135 iso_sp[ipHE_LIKE][ipOXYGEN].n_HighestResolved_max = 5;
136 iso_sp[ipHE_LIKE][ipNEON].n_HighestResolved_max = 5;
137 iso_sp[ipHE_LIKE][ipSILICON].n_HighestResolved_max = 5;
138 iso_sp[ipHE_LIKE][ipMAGNESIUM].n_HighestResolved_max = 5;
139 iso_sp[ipHE_LIKE][ipSULPHUR].n_HighestResolved_max = 5;
140 iso_sp[ipHE_LIKE][ipIRON].n_HighestResolved_max = 5;
141 /* also set this, for exercising any possible issues with biggest charge models */
142 iso_sp[ipHE_LIKE][LIMELM-1].n_HighestResolved_max = 5;
143
144 iso_ctrl.chISO[ipH_LIKE] = "H-like ";
145 iso_ctrl.chISO[ipHE_LIKE] = "He-like";
146
147 max_num_levels = 0;
148 for( long ipISO = ipH_LIKE; ipISO < NISO; ipISO++ )
149 {
150 for( nelem=ipISO; nelem < LIMELM; nelem++ )
151 {
152 /* set this to LONG_MAX, reduce to actual number later,
153 * then verify number of levels is not increased after initial coreload */
154 iso_sp[ipISO][nelem].numLevels_malloc = LONG_MAX;
155 iso_update_num_levels( ipISO, nelem );
156 }
157 }
158
159 lgStatesAdded = false;
160 lgLinesAdded = false;
161
162 /* turn element on if this is first call to this routine,
163 * but do not reset otherwise since would clobber how elements were set up */
164 for( nelem=0; nelem < LIMELM; nelem++ )
165 {
166 /* never turn on element if turned off on first iteration */
167 dense.lgElmtOn[nelem] = true;
168
169 /* option to set ionization with element ioniz cmnd,
170 * default (false) is to solve for ionization */
171 dense.lgSetIoniz[nelem] = false;
172
173 // are we using Chianti for this ion?
174 for( int ion=0; ion<=nelem+1; ++ion )
175 {
176 dense.lgIonChiantiOn[nelem][ion] = false;
177 dense.lgIonStoutOn[nelem][ion] = false;
178 dense.maxWN[nelem][ion] = 0.;
179 }
180 }
181
182 /* smallest density to permit in any ion - if smaller then set to zero */
183 dense.density_low_limit = SMALLFLOAT * 1e3;
184 dense.density_low_limit = MAX2( dense.density_low_limit, 1e-50 );
185
186 /* default cosmic ray background */
187 /* >>chng 99 jun 24, slight change to value
188 * quoted by
189 * >>refer cosmic ray ionization rate McKee, C.M., 1999, astro-ph 9901370
190 * this will produce a total
191 * secondary ionization rate of 2.5e-17 s^-1, as tested in
192 * test suite cosmicray.in. If each ionization produces 2.4 eV of heat,
193 * the background heating rate should be 9.6e-29 * n*/
194 /* >>chng 04 jan 26, update cosmic ray ionization rate for H0 to
195 * >>refer cosmic ray ionization Williams, J.P., Bergin, E.A., Caselli, P.,
196 * >>refercon Myers, P.C., & Plume, R. 1998, ApJ, 503, 689,
197 * H0 ionization rate of 2.5e-17 s-1 and a H2 ionization rate twice this
198 * >>chng 04 mar 15, comment said 2.5e-17 which is correct, but code produce 8e-17,
199 * fix back to correct value
200 */
201 /* NB - the rate is derived from the density. these two are related by the secondary
202 * ionization efficiency problem. background_rate is only here to provide the relationship
203 * for predominantly neutral gas. the background_density is the real rate.
204 hextra.background_density = 1.99e-9f;*/
205 /* >>chng 05 apr 16, to get proper ionization rate in ism_set_cr_rate, where
206 * H is forced to be fully atomic, no molecules, density from 1.99 to 2.15 */
207 /* >>chng 02 apr 05, update to
208 * >>refer cosmic ray ionization Indriolo, N., Geballe, T., Oka, T., & McCall, B.J. 2007, ApJ, 671, 1736
209 */
210 hextra.background_density = 2.15e-9f*7.9f;
211 hextra.background_rate = 2.5e-17f*7.9f;
212
213 /* initialization for save files - must call after input commands have
214 * been scanned for grid and vary options. So known if grid or single run
215 * default save is different for these */
216 grid.lgGridDone = false;
217 grid.lgStrictRepeat = false;
218
219 /* these are energy range... if not changed with command, 0. says just use energy limits of mesh */
220 grid.LoEnergy_keV = 0.;
221 grid.HiEnergy_keV = 0.;
222 grid.ipLoEnergy = 0;
223 grid.ipHiEnergy = 0;
224 grid.seqNum = 0;
225
226 for( long i=0; i < LIMPAR; ++i )
227 optimize.lgOptimizeAsLinear[i] = false;
228
229 /* limit on ionization we check for zoning and prtcomment */
230 struc.dr_ionfrac_limit = 1e-3f;
231
233
234 diatoms_init();
235
236 /* initialize cosmological information */
237 cosmology.redshift_current = 0.f;
238 cosmology.redshift_start = 0.f;
239 cosmology.redshift_step = 0.f;
240 cosmology.omega_baryon = 0.04592f;
241 cosmology.omega_rad = 8.23e-5f;
242 cosmology.omega_lambda = 0.7299177f;
243 cosmology.omega_matter = 0.27f;
244 cosmology.omega_k = 0.f;
245 /* the Hubble parameter in 100 km/s/Mpc */
246 cosmology.h = 0.71f;
247 /* the Hubble parameter in km/s/Mpc */
248 cosmology.H_0 = 100.f*cosmology.h;
249 cosmology.lgDo = false;
250
251 // the code is ok at startup; only init here so that code remains broken
252 // in a grid if any single model announces broken
253 broke.lgBroke = false;
254 broke.lgFixit = false;
255 broke.lgCheckit = false;
256
257 return;
258}
t_broke broke
Definition broke.cpp:5
const int ipSILICON
Definition cddefines.h:318
const int ipNITROGEN
Definition cddefines.h:311
const int ipIRON
Definition cddefines.h:330
const int ipCARBON
Definition cddefines.h:310
const int ipOXYGEN
Definition cddefines.h:312
const int LIMELM
Definition cddefines.h:258
const int ipMAGNESIUM
Definition cddefines.h:316
const int NISO
Definition cddefines.h:261
const int ipHELIUM
Definition cddefines.h:306
#define MAX2
Definition cddefines.h:782
const int ipSULPHUR
Definition cddefines.h:320
const int ipHYDROGEN
Definition cddefines.h:305
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
const int ipNEON
Definition cddefines.h:314
t_cosmology cosmology
Definition cosmology.cpp:11
const realnum SMALLFLOAT
Definition cpu.h:191
t_dense dense
Definition dense.cpp:24
t_grid grid
Definition grid.cpp:5
void diatoms_init(void)
Definition h2.cpp:13
t_hcmap hcmap
Definition hcmap.cpp:21
t_Heavy Heavy
Definition heavy.cpp:5
t_hextra hextra
Definition hextra.cpp:5
void InitCoreload(void)
t_input input
Definition input.cpp:12
t_ionbal ionbal
Definition ionbal.cpp:5
t_iso_sp iso_sp[NISO][LIMELM]
Definition iso.cpp:8
long int max_num_levels
Definition iso.cpp:10
t_isoCTRL iso_ctrl
Definition iso.cpp:6
void iso_update_num_levels(long ipISO, long nelem)
const int ipHE_LIKE
Definition iso.h:63
const int ipH_LIKE
Definition iso.h:62
char chOptimFileName[INPUT_LINE_LENGTH]
Definition optimize.cpp:6
t_optimize optimize
Definition optimize.cpp:5
const long LIMPAR
Definition optimize.h:61
void SaveFilesInit(void)
t_rfield rfield
Definition rfield.cpp:8
t_struc struc
Definition struc.cpp:6
bool lgStatesAdded
Definition taulines.cpp:10
bool lgLinesAdded
Definition taulines.cpp:11