cloudy trunk
Loading...
Searching...
No Matches
init_sim_postparse.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/*InitSimPostparse initialization at start of simulation,
4 * called from cloudy after parser,
5 * called one time per sim in grid
6 * not called for every iteration */
7#include "cddefines.h"
8#include "dense.h"
9#include "deuterium.h"
10#include "struc.h"
11#include "rfield.h"
12#include "elementnames.h"
13#include "atoms.h"
14#include "physconst.h"
15#include "iterations.h"
16#include "pressure.h"
17#include "mole.h"
18#include "trace.h"
19#include "radius.h"
20#include "thermal.h"
21#include "heavy.h"
22#include "wind.h"
23#include "init.h"
24#include "iso.h"
25#include "h2.h"
26#include "co.h"
27#include "monitor_results.h"
28#include "taulines.h"
29
30/*InitSimPostparse initialization after parser,
31 * called one time for each simulation in a grid,
32 * only before first iteration
33 * sets initial or zeros values before start of a simulation */
34void InitSimPostparse( void )
35{
36 DEBUG_ENTRY( "InitSimPostparse()" );
37 static double one=1.0;
38
39 struc.nzonePreviousIteration = -1;
40
41 thermal.thist = 0.;
42 thermal.tlowst = 1e20f;
43 thermal.nUnstable = 0;
44 thermal.lgUnstable = false;
45
46 colliders.init();
47
48 mole_global.init();
49 mole_global.zero();
50
51 for( long nelem=ipHYDROGEN; nelem<LIMELM; ++nelem )
52 {
53 for( long ion=0; ion<nelem+2; ++ion )
54 mole.set_location( nelem, ion, &(dense.xIonDense[nelem][ion]) );
55 }
56
57 findspecieslocal("e-")->location = &(dense.eden);
58 findspecieslocal("grn")->location = &(mole.grain_area);
59 findspecieslocal("PHOTON")->location = &one;
60 findspecieslocal("CRPHOT")->location = &one;
61 findspecieslocal("CRP")->location = &one;
62 if( deut.lgElmtOn )
63 {
64 findspecieslocal("D")->location = &(deut.xIonDense[0]);
65 findspecieslocal("D+")->location = &(deut.xIonDense[1]);
66 }
67
68 /* >> chng 06 Nov 24 rjrw: Move reaction definitions here so they can be controlled by parsed commands */
70
71 //mole_cmp_num_in_out_reactions();
72
73 for( diatom_iter diatom = diatoms.begin(); diatom != diatoms.end(); ++diatom )
74 (*diatom)->H2_Reset();
75
76 /* zero out diffuse Lya emission */
77 for( long nelem=ipHYDROGEN; nelem<LIMELM; ++nelem )
78 for( long ion=0; ion<=nelem; ++ion )
79 Heavy.xLyaHeavy[nelem][ion] = 0;
80
81 /* convert STOP RADIUS command to STOP THICKNESS command now that inner radius is known */
82 if( radius.StopRadius[0] > 0. )
83 {
84 for( long j=0; j < iterations.iter_malloc; j++ )
85 {
86 if( radius.StopRadius[j] > radius.Radius )
87 radius.StopThickness[j] = radius.StopRadius[j] - radius.Radius;
88 else
89 {
90 fprintf( ioQQQ, " PROBLEM stopping radius is <= inner radius. Bailing out.\n" );
92 }
93 }
94 }
95
96 /* this term applies centrifugal acceleration if DISK option on wind command */
97 wind.DiskRadius = 0;
98 if( wind.lgDisk )
99 wind.DiskRadius = radius.Radius;
100
101 iterations.lgLastIt = false;
102
103 rfield_opac_zero( 0 , rfield.nupper );
104
105 rfield.lgUSphON = false;
106
107 /* where cloud is optically thick to brems */
108 rfield.ipEnergyBremsThin = 0;
109 rfield.EnergyBremsThin = 0.;
110 rfield.comtot = 0.;
111 rfield.cmcool = 0.;
112 rfield.cinrat = 0.;
113 rfield.extin_mag_B_point = 0.;
114 rfield.extin_mag_V_point = 0.;
115 rfield.extin_mag_B_extended = 0.;
116 rfield.extin_mag_V_extended = 0.;
117 rfield.EnerGammaRay = 7676.;
118
119 for( vector<TransitionList>::iterator it = AllTransitions.begin(); it != AllTransitions.end(); ++it )
120 {
121 for( TransitionList::iterator tr = it->begin(); tr != it->end(); ++tr )
122 {
123 tr->Emis().TauTrack().clear();
124 }
125 }
126
127 /* usually computed in pressure_change, but not called before first zone */
128 wind.AccelGravity = (realnum)(GRAV_CONST*wind.comass*SOLAR_MASS/POW2(radius.Radius)*
129 (1.-wind.DiskRadius/radius.Radius) );
130 if( trace.lgTrace && trace.lgWind )
131 {
132 fprintf(ioQQQ," InitSimPostparse sets AccelGravity %.3e lgDisk?%c\n",
133 wind.AccelGravity ,
134 TorF(wind.lgDisk) );
135 }
136
137 /* pressure related variables */
138 pressure.PresInteg = 0.;
139 pressure.PresIntegElecThin = 0.;
140 pressure.pinzon = 0.;
141
142 /* update iso sequence level information */
143 for( long ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
144 {
145 for( long nelem=ipISO; nelem < LIMELM; ++nelem )
146 {
147 iso_sp[ipISO][nelem].Reset();
148 /* these elements that are turned off */
149 if( nelem>ipHELIUM && !dense.lgElmtOn[nelem] )
150 {
151 iso_sp[ipISO][nelem].numLevels_max = 0;
152 iso_sp[ipISO][nelem].nCollapsed_max = 0;
153 iso_sp[ipISO][nelem].n_HighestResolved_max = 0;
154
155 iso_sp[ipISO][nelem].numLevels_local = 0;
156 iso_sp[ipISO][nelem].nCollapsed_local = 0;
157 iso_sp[ipISO][nelem].n_HighestResolved_local = 0;
158 }
159 else
160 {
161 iso_update_num_levels( ipISO, nelem );
162 /* must always have at least one collapsed level, unless compiling recomb data file. */
163 ASSERT( iso_sp[ipISO][nelem].nCollapsed_max > 0 || iso_ctrl.lgCompileRecomb[ipISO] );
164 }
165 }
166 }
167
168 if( iso_ctrl.lgPrintNumberOfLevels )
169 {
170 fprintf(ioQQQ,"\n\n Number of levels in ions treated by iso sequences.\n");
171 fprintf(ioQQQ," ISO Element hi-n(l-resolved) #(l-resolved) n(collapsed)\n");
172 /* option to print number of levels for each element */
173 for( long ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
174 {
175 for( long nelem=ipISO; nelem<LIMELM; ++nelem )
176 {
177 /* print number of levels */
178 fprintf(ioQQQ," %s %s %4li %4li %4li \n",
179 iso_ctrl.chISO[ipISO] ,
180 elementnames.chElementSym[nelem],
181 iso_sp[ipISO][nelem].n_HighestResolved_max,
182 iso_sp[ipISO][nelem].numLevels_max-iso_sp[ipISO][nelem].nCollapsed_max,
183 iso_sp[ipISO][nelem].nCollapsed_max);
184 }
185 }
186 }
187 atoms.p2nit = 0.;
188 atoms.d5200r = 0.;
189
190 MonitorResults.SumErrorCaseMonitor = 0.;
191 MonitorResults.nSumErrorCaseMonitor = 0;
192
193 return;
194}
t_atoms atoms
Definition atoms.cpp:5
FILE * ioQQQ
Definition cddefines.cpp:7
#define ASSERT(exp)
Definition cddefines.h:578
const int LIMELM
Definition cddefines.h:258
#define POW2
Definition cddefines.h:929
#define EXIT_FAILURE
Definition cddefines.h:140
char TorF(bool l)
Definition cddefines.h:710
#define cdEXIT(FAIL)
Definition cddefines.h:434
const int NISO
Definition cddefines.h:261
const int ipHELIUM
Definition cddefines.h:306
float realnum
Definition cddefines.h:103
const int ipHYDROGEN
Definition cddefines.h:305
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
TransitionProxy::iterator iterator
Definition transition.h:280
double * location
Definition mole.h:349
ColliderList colliders
Definition collision.cpp:57
void rfield_opac_zero(long lo, long ihi)
t_dense dense
Definition dense.cpp:24
t_deuterium deut
Definition deuterium.cpp:8
t_elementnames elementnames
vector< diatomics * > diatoms
Definition h2.cpp:8
vector< diatomics * >::iterator diatom_iter
Definition h2.h:13
t_Heavy Heavy
Definition heavy.cpp:5
void InitSimPostparse(void)
t_iso_sp iso_sp[NISO][LIMELM]
Definition iso.cpp:8
t_isoCTRL iso_ctrl
Definition iso.cpp:6
void iso_update_num_levels(long ipISO, long nelem)
const int ipH_LIKE
Definition iso.h:62
t_iterations iterations
Definition iterations.cpp:5
t_mole_global mole_global
Definition mole.cpp:6
t_mole_local mole
Definition mole.cpp:7
void mole_create_react(void)
molezone * findspecieslocal(const char buf[])
t_monitorresults MonitorResults
UNUSED const double SOLAR_MASS
Definition physconst.h:71
UNUSED const double GRAV_CONST
Definition physconst.h:109
t_pressure pressure
Definition pressure.cpp:5
t_radius radius
Definition radius.cpp:5
t_rfield rfield
Definition rfield.cpp:8
t_struc struc
Definition struc.cpp:6
vector< TransitionList > AllTransitions
Definition taulines.cpp:8
t_thermal thermal
Definition thermal.cpp:5
t_trace trace
Definition trace.cpp:5
Wind wind
Definition wind.cpp:5