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