cloudy
trunk
|
00001 /* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and 00002 * others. For conditions of distribution and use see copyright notice in license.txt */ 00003 /*AtomSeqBoron compute cooling from 5-level boron sequence model atom */ 00004 #include "cddefines.h" 00005 #include "cooling.h" 00006 #include "thermal.h" 00007 #include "dense.h" 00008 #include "atoms.h" 00009 00010 void AtomSeqBoron( 00011 /* indices for all lines are on the C scale since they will be stuffed into 00012 * C arrays. so, t10 refers to the 2-1 transition */ 00013 transition * t10, 00014 transition * t20, 00015 transition * t30, 00016 transition * t21, 00017 transition * t31, 00018 transition * t41, 00019 double cs40, 00020 double cs32, 00021 double cs42, 00022 double cs43, 00023 /* pump rate s-1 due to UV permitted lines */ 00024 double pump_rate , 00025 /* string used to identify calling program in case of error */ 00026 const char *chLabel 00027 ) 00028 { 00029 00030 /* this routine has three possible returns: 00031 * abundance is zero 00032 * too cool for full 5-level atom, but still do ground term 00033 * full solution */ 00034 00035 /* boron sequence is now a five level atom */ 00036 # define N_SEQ_BORON 5 00037 static double 00038 **AulEscp , 00039 **col_str , 00040 **AulDest, 00041 /* AulPump[low][high] is rate (s^-1) from lower to upper level */ 00042 **AulPump, 00043 **CollRate, 00044 *pops, 00045 *create, 00046 *destroy, 00047 *depart, 00048 /* statistical weight */ 00049 *stat , 00050 /* excitation energies in kelvin */ 00051 *excit; 00052 00053 double b_cooling, 00054 dCoolDT; 00055 double EnrLU, EnrUL; 00056 realnum abundan; 00057 00058 static bool lgFirst=true, 00059 lgZeroPop; 00060 int i , j; 00061 int 00062 /* flag to signal negative level populations */ 00063 lgNegPop; 00064 /* flag to turn on debug print in atom_levelN */ 00065 bool lgDeBug; 00066 00067 DEBUG_ENTRY( "AtomSeqBoron()" ); 00068 00069 if( lgFirst ) 00070 { 00071 /* will never do this again */ 00072 lgFirst = false; 00073 /* allocate the 1D arrays*/ 00074 excit = (double *)MALLOC( sizeof(double)*(N_SEQ_BORON) ); 00075 stat = (double *)MALLOC( sizeof(double)*(N_SEQ_BORON) ); 00076 pops = (double *)MALLOC( sizeof(double)*(N_SEQ_BORON) ); 00077 create = (double *)MALLOC( sizeof(double)*(N_SEQ_BORON) ); 00078 destroy = (double *)MALLOC( sizeof(double)*(N_SEQ_BORON) ); 00079 depart = (double *)MALLOC( sizeof(double)*(N_SEQ_BORON) ); 00080 /* create space for the 2D arrays */ 00081 AulPump = ((double **)MALLOC((N_SEQ_BORON)*sizeof(double *))); 00082 CollRate = ((double **)MALLOC((N_SEQ_BORON)*sizeof(double *))); 00083 AulDest = ((double **)MALLOC((N_SEQ_BORON)*sizeof(double *))); 00084 AulEscp = ((double **)MALLOC((N_SEQ_BORON)*sizeof(double *))); 00085 col_str = ((double **)MALLOC((N_SEQ_BORON)*sizeof(double *))); 00086 for( i=0; i<(N_SEQ_BORON); ++i ) 00087 { 00088 AulPump[i] = ((double *)MALLOC((N_SEQ_BORON)*sizeof(double ))); 00089 CollRate[i] = ((double *)MALLOC((N_SEQ_BORON)*sizeof(double ))); 00090 AulDest[i] = ((double *)MALLOC((N_SEQ_BORON)*sizeof(double ))); 00091 AulEscp[i] = ((double *)MALLOC((N_SEQ_BORON)*sizeof(double ))); 00092 col_str[i] = ((double *)MALLOC((N_SEQ_BORON)*sizeof(double ))); 00093 } 00094 } 00095 00096 /* total abundance of this species */ 00097 abundan = dense.xIonDense[ t10->Hi->nelem -1][t10->Hi->IonStg-1]; 00098 00100 if( abundan <= 0. ) 00101 { 00102 /* this branch, no abundance of ion */ 00103 t10->Lo->Pop = 0.; 00104 t20->Lo->Pop = 0.; 00105 t30->Lo->Pop = 0.; 00106 t21->Lo->Pop = 0.; 00107 t31->Lo->Pop = 0.; 00108 t41->Lo->Pop = 0.; 00109 00110 t10->Emis->PopOpc = 0.; 00111 t20->Emis->PopOpc = 0.; 00112 t30->Emis->PopOpc = 0.; 00113 t21->Emis->PopOpc = 0.; 00114 t31->Emis->PopOpc = 0.; 00115 t41->Emis->PopOpc = 0.; 00116 00117 t10->Hi->Pop = 0.; 00118 t20->Hi->Pop = 0.; 00119 t30->Hi->Pop = 0.; 00120 t21->Hi->Pop = 0.; 00121 t31->Hi->Pop = 0.; 00122 t41->Hi->Pop = 0.; 00123 00124 t10->Emis->xIntensity = 0.; 00125 t20->Emis->xIntensity = 0.; 00126 t30->Emis->xIntensity = 0.; 00127 t21->Emis->xIntensity = 0.; 00128 t31->Emis->xIntensity = 0.; 00129 t41->Emis->xIntensity = 0.; 00130 00131 t10->Coll.cool = 0.; 00132 t20->Coll.cool = 0.; 00133 t30->Coll.cool = 0.; 00134 t21->Coll.cool = 0.; 00135 t31->Coll.cool = 0.; 00136 t41->Coll.cool = 0.; 00137 00138 t10->Emis->phots = 0.; 00139 t20->Emis->phots = 0.; 00140 t30->Emis->phots = 0.; 00141 t21->Emis->phots = 0.; 00142 t31->Emis->phots = 0.; 00143 t41->Emis->phots = 0.; 00144 00145 t10->Emis->ColOvTot = 0.; 00146 t20->Emis->ColOvTot = 0.; 00147 t30->Emis->ColOvTot = 0.; 00148 t21->Emis->ColOvTot = 0.; 00149 t31->Emis->ColOvTot = 0.; 00150 t41->Emis->ColOvTot = 0.; 00151 00152 t10->Coll.heat = 0.; 00153 t20->Coll.heat = 0.; 00154 t30->Coll.heat = 0.; 00155 t21->Coll.heat = 0.; 00156 t31->Coll.heat = 0.; 00157 t41->Coll.heat = 0.; 00158 00159 CoolAdd( chLabel, t10->WLAng , 0.); 00160 CoolAdd( chLabel, t20->WLAng , 0.); 00161 CoolAdd( chLabel, t30->WLAng , 0.); 00162 CoolAdd( chLabel, t21->WLAng , 0.); 00163 CoolAdd( chLabel, t31->WLAng , 0.); 00164 CoolAdd( chLabel, t41->WLAng , 0.); 00165 00166 /* level populations */ 00167 /* LIMLEVELN is the dimension of the atoms vectors */ 00168 ASSERT( N_SEQ_BORON <= LIMLEVELN); 00169 for( i=0; i < N_SEQ_BORON; i++ ) 00170 { 00171 atoms.PopLevels[i] = 0.; 00172 atoms.DepLTELevels[i] = 1.; 00173 } 00174 return; 00175 } 00176 00177 ASSERT( t10->Coll.cs > 0.); 00178 ASSERT( t20->Coll.cs > 0.); 00179 ASSERT( t30->Coll.cs > 0.); 00180 ASSERT( t21->Coll.cs > 0.); 00181 ASSERT( t31->Coll.cs > 0.); 00182 ASSERT( t41->Coll.cs > 0.); 00183 ASSERT( cs40>0.); 00184 ASSERT( cs32>0.); 00185 ASSERT( cs42>0.); 00186 ASSERT( cs43>0.); 00187 00188 /* all elements are used, and must be set to zero if zero */ 00189 for( i=0; i < N_SEQ_BORON; i++ ) 00190 { 00191 create[i] = 0.; 00192 destroy[i] = 0.; 00193 for( j=0; j < N_SEQ_BORON; j++ ) 00194 { 00195 /*data[j][i] = -1e33;*/ 00196 AulEscp[j][i] = 0.; 00197 AulDest[j][i] = 0.; 00198 AulPump[j][i] = 0.; 00199 col_str[j][i] = 0.; 00200 } 00201 } 00202 00203 /* statistical weights */ 00204 stat[0] = t10->Lo->g; 00205 stat[1] = t10->Hi->g; 00206 stat[2] = t20->Hi->g; 00207 stat[3] = t30->Hi->g; 00208 stat[4] = t41->Hi->g; 00209 ASSERT( stat[0]>0. && stat[1]>0. &&stat[2]>0. &&stat[3]>0. &&stat[4]>0.); 00210 ASSERT( fabs(t10->Lo->g/2.-1.) < FLT_EPSILON); 00211 ASSERT( fabs(t10->Hi->g/4.-1.) < FLT_EPSILON); 00212 ASSERT( fabs(t20->Lo->g/2.-1.) < FLT_EPSILON); 00213 ASSERT( fabs(t20->Hi->g/2.-1.) < FLT_EPSILON); 00214 ASSERT( fabs(t30->Lo->g/2.-1.) < FLT_EPSILON); 00215 ASSERT( fabs(t30->Hi->g/4.-1.) < FLT_EPSILON); 00216 ASSERT( fabs(t21->Lo->g/4.-1.) < FLT_EPSILON); 00217 ASSERT( fabs(t21->Hi->g/2.-1.) < FLT_EPSILON); 00218 ASSERT( fabs(t31->Lo->g/4.-1.) < FLT_EPSILON); 00219 ASSERT( fabs(t31->Hi->g/4.-1.) < FLT_EPSILON); 00220 ASSERT( fabs(t41->Lo->g/4.-1.) < FLT_EPSILON); 00221 ASSERT( fabs(t41->Hi->g/6.-1.) < FLT_EPSILON); 00222 00223 /* excitation energy of each level relative to ground, in Kelvin */ 00224 excit[0] = 0.; 00225 excit[1] = t10->EnergyK; 00226 excit[2] = t20->EnergyK; 00227 excit[3] = t30->EnergyK; 00228 excit[4] = t41->EnergyK + t10->EnergyK; 00229 ASSERT( excit[1]>0. &&excit[2]>0. &&excit[3]>0. &&excit[4]>0.); 00230 00231 /* fill in Einstein As, collision strengths, pumping rates */ 00232 AulEscp[1][0] = t10->Emis->Aul*(t10->Emis->Pesc + t10->Emis->Pelec_esc); 00233 AulDest[1][0] = t10->Emis->Aul*t10->Emis->Pdest; 00234 col_str[1][0] = t10->Coll.cs; 00235 AulPump[0][1] = t10->Emis->pump; 00236 00237 /* add FUV pump transitions to this pump rate */ 00238 AulPump[0][1] += pump_rate; 00239 00240 AulEscp[2][0] = t20->Emis->Aul*(t20->Emis->Pesc + t20->Emis->Pelec_esc); 00241 AulDest[2][0] = t20->Emis->Aul*t20->Emis->Pdest; 00242 col_str[2][0] = t20->Coll.cs; 00243 AulPump[0][2] = t20->Emis->pump; 00244 00245 AulEscp[3][0] = t30->Emis->Aul*(t30->Emis->Pesc + t30->Emis->Pelec_esc); 00246 AulDest[3][0] = t30->Emis->Aul*t30->Emis->Pdest; 00247 col_str[3][0] = t30->Coll.cs; 00248 AulPump[0][3] = t30->Emis->pump; 00249 00250 AulEscp[4][0] = 1e-8;/* made up trans prob */ 00251 AulDest[4][0] = 0.; 00252 col_str[4][0] = cs40; 00253 AulPump[0][4] = 0.; 00254 00255 AulEscp[2][1] = t21->Emis->Aul*(t21->Emis->Pesc + t21->Emis->Pelec_esc); 00256 AulDest[2][1] = t21->Emis->Aul*t21->Emis->Pdest; 00257 col_str[2][1] = t21->Coll.cs; 00258 AulPump[1][2] = t21->Emis->pump; 00259 00260 AulEscp[3][1] = t31->Emis->Aul*(t31->Emis->Pesc + t31->Emis->Pelec_esc); 00261 AulDest[3][1] = t31->Emis->Aul*t31->Emis->Pdest; 00262 col_str[3][1] = t31->Coll.cs; 00263 AulPump[1][3] = t31->Emis->pump; 00264 00265 AulEscp[4][1] = t41->Emis->Aul*(t41->Emis->Pesc + t41->Emis->Pelec_esc); 00266 AulDest[4][1] = t41->Emis->Aul*t41->Emis->Pdest; 00267 col_str[4][1] = t41->Coll.cs; 00268 AulPump[1][4] = t41->Emis->pump; 00269 00270 AulEscp[3][2] = 1e-8;/* made up trans prob */ 00271 AulDest[3][2] = 0.; 00272 col_str[3][2] = cs32; 00273 AulPump[2][3] = 0.; 00274 00275 AulEscp[4][2] = 1e-8;/* made up trans prob */ 00276 AulDest[4][2] = 0.; 00277 col_str[4][2] = cs42; 00278 AulPump[2][4] = 0.; 00279 00280 AulEscp[4][3] = 1e-8;/* made up trans prob */ 00281 AulDest[4][3] = 0.; 00282 col_str[4][3] = cs43; 00283 AulPump[3][4] = 0.; 00284 00285 lgDeBug = false; 00286 00287 /* lgNegPop positive if negative pops occurred, negative if too cold */ 00288 atom_levelN(N_SEQ_BORON, 00289 abundan, 00290 stat, 00291 excit, 00292 'K', 00293 pops, 00294 depart, 00295 &AulEscp, 00296 &col_str, 00297 &AulDest, 00298 &AulPump, 00299 &CollRate, 00300 create, 00301 destroy, 00302 false,/* say atom_levelN should evaluate coll rates from cs */ 00303 &b_cooling, 00304 &dCoolDT, 00305 chLabel, 00306 &lgNegPop, 00307 &lgZeroPop, 00308 lgDeBug );/* option to print stuff - set to true for debug printout */ 00309 00310 /* this returned lgNegPop < 0 in case were gas too cool to excite highest level. 00311 * in this case we still want to evaluate the cooling due to the split ground term 00312 * in this sequence */ 00313 if( lgNegPop < 0 ) 00314 { 00315 /* it was too cool for atom_levelN to do whole atom, but let's still do the 00316 * split ground term - this routine adds cooling and its temperature 00317 * derivative to total cooling and dD/dT */ 00318 atom_level2(t10); 00319 00320 /* now put in pop for absorption into higher levels, for line optical depths, 00321 * PopLevels was evaluated by atom_level2 */ 00322 t20->Lo->Pop = atoms.PopLevels[0]; 00323 t30->Lo->Pop = atoms.PopLevels[0]; 00324 t21->Lo->Pop = atoms.PopLevels[1]; 00325 t31->Lo->Pop = atoms.PopLevels[1]; 00326 t41->Lo->Pop = atoms.PopLevels[1]; 00327 00328 t20->Emis->PopOpc = atoms.PopLevels[0]; 00329 t30->Emis->PopOpc = atoms.PopLevels[0]; 00330 t21->Emis->PopOpc = atoms.PopLevels[1]; 00331 t31->Emis->PopOpc = atoms.PopLevels[1]; 00332 t41->Emis->PopOpc = atoms.PopLevels[1]; 00333 00334 t20->Hi->Pop = 0.; 00335 t30->Hi->Pop = 0.; 00336 t21->Hi->Pop = 0.; 00337 t31->Hi->Pop = 0.; 00338 t41->Hi->Pop = 0.; 00339 00340 t20->Emis->xIntensity = 0.; 00341 t30->Emis->xIntensity = 0.; 00342 t21->Emis->xIntensity = 0.; 00343 t31->Emis->xIntensity = 0.; 00344 t41->Emis->xIntensity = 0.; 00345 00346 t20->Coll.cool = 0.; 00347 t30->Coll.cool = 0.; 00348 t21->Coll.cool = 0.; 00349 t31->Coll.cool = 0.; 00350 t41->Coll.cool = 0.; 00351 00352 t20->Emis->phots = 0.; 00353 t30->Emis->phots = 0.; 00354 t21->Emis->phots = 0.; 00355 t31->Emis->phots = 0.; 00356 t41->Emis->phots = 0.; 00357 00358 t20->Emis->ColOvTot = 0.; 00359 t30->Emis->ColOvTot = 0.; 00360 t21->Emis->ColOvTot = 0.; 00361 t31->Emis->ColOvTot = 0.; 00362 t41->Emis->ColOvTot = 0.; 00363 00364 t20->Coll.heat = 0.; 00365 t30->Coll.heat = 0.; 00366 t21->Coll.heat = 0.; 00367 t31->Coll.heat = 0.; 00368 t41->Coll.heat = 0.; 00369 00370 CoolAdd( chLabel, t20->WLAng , 0.); 00371 CoolAdd( chLabel, t30->WLAng , 0.); 00372 CoolAdd( chLabel, t21->WLAng , 0.); 00373 CoolAdd( chLabel, t31->WLAng , 0.); 00374 CoolAdd( chLabel, t41->WLAng , 0.); 00375 00376 /* level populations */ 00377 /* LIMLEVELN is the dimension of the atoms vectors */ 00378 ASSERT( N_SEQ_BORON <= LIMLEVELN); 00379 for( i=2; i < N_SEQ_BORON; i++ ) 00380 { 00381 atoms.PopLevels[i] = 0.; 00382 atoms.DepLTELevels[i] = 1.; 00383 } 00384 } 00385 00386 else 00387 { 00388 /* atom_levelN did not evaluate PopLevels, so save pops here */ 00389 /* LIMLEVELN is the dimension of the atoms vectors */ 00390 ASSERT( N_SEQ_BORON <= LIMLEVELN); 00391 for( i=0; i< N_SEQ_BORON; ++i ) 00392 { 00393 atoms.PopLevels[i] = pops[i]; 00394 atoms.DepLTELevels[i] = depart[i]; 00395 } 00396 /* this branch, we have a full valid solution */ 00397 t10->Lo->Pop = pops[0]; 00398 t20->Lo->Pop = pops[0]; 00399 t30->Lo->Pop = pops[0]; 00400 t21->Lo->Pop = pops[1]; 00401 t31->Lo->Pop = pops[1]; 00402 t41->Lo->Pop = pops[1]; 00403 00404 t10->Emis->PopOpc = (pops[0] - pops[1]*t10->Lo->g/t10->Hi->g); 00405 t20->Emis->PopOpc = (pops[0] - pops[2]*t20->Lo->g/t20->Hi->g); 00406 t30->Emis->PopOpc = (pops[0] - pops[3]*t30->Lo->g/t30->Hi->g); 00407 t21->Emis->PopOpc = (pops[1] - pops[2]*t21->Lo->g/t21->Hi->g); 00408 t31->Emis->PopOpc = (pops[1] - pops[3]*t31->Lo->g/t31->Hi->g); 00409 t41->Emis->PopOpc = (pops[1] - pops[4]*t41->Lo->g/t41->Hi->g); 00410 00411 t10->Hi->Pop = pops[1]; 00412 t20->Hi->Pop = pops[2]; 00413 t30->Hi->Pop = pops[3]; 00414 t21->Hi->Pop = pops[2]; 00415 t31->Hi->Pop = pops[3]; 00416 t41->Hi->Pop = pops[4]; 00417 00418 t10->Emis->phots = t10->Emis->Aul*(t10->Emis->Pesc + t10->Emis->Pelec_esc)*pops[1]; 00419 t20->Emis->phots = t20->Emis->Aul*(t20->Emis->Pesc + t20->Emis->Pelec_esc)*pops[2]; 00420 t30->Emis->phots = t30->Emis->Aul*(t30->Emis->Pesc + t30->Emis->Pelec_esc)*pops[3]; 00421 t21->Emis->phots = t21->Emis->Aul*(t21->Emis->Pesc + t21->Emis->Pelec_esc)*pops[2]; 00422 t31->Emis->phots = t31->Emis->Aul*(t31->Emis->Pesc + t31->Emis->Pelec_esc)*pops[3]; 00423 t41->Emis->phots = t41->Emis->Aul*(t41->Emis->Pesc + t41->Emis->Pelec_esc)*pops[4]; 00424 00425 t10->Emis->xIntensity = t10->Emis->phots*t10->EnergyErg; 00426 t20->Emis->xIntensity = t20->Emis->phots*t20->EnergyErg; 00427 t30->Emis->xIntensity = t30->Emis->phots*t30->EnergyErg; 00428 t21->Emis->xIntensity = t21->Emis->phots*t21->EnergyErg; 00429 t31->Emis->xIntensity = t31->Emis->phots*t31->EnergyErg; 00430 t41->Emis->xIntensity = t41->Emis->phots*t41->EnergyErg; 00431 00432 /* ratio of collisional to total excitation */ 00433 t10->Emis->ColOvTot = CollRate[0][1]/(CollRate[0][1]+t10->Emis->pump); 00434 t20->Emis->ColOvTot = CollRate[0][2]/(CollRate[0][2]+t20->Emis->pump); 00435 t30->Emis->ColOvTot = CollRate[0][3]/(CollRate[0][3]+t30->Emis->pump); 00436 t21->Emis->ColOvTot = CollRate[1][2]/(CollRate[1][2]+t21->Emis->pump); 00437 t31->Emis->ColOvTot = CollRate[1][3]/(CollRate[1][3]+t31->Emis->pump); 00438 t41->Emis->ColOvTot = CollRate[1][4]/(CollRate[1][4]+t41->Emis->pump); 00439 00440 /* derivative of cooling function */ 00441 thermal.dCooldT += dCoolDT; 00442 00443 /* two cases - collisionally excited (usual case) or 00444 * radiatively excited - in which case line can be a heat source 00445 * following are correct heat exchange, they will mix to get correct derivative 00446 * the sum of heat-cool will add up to EnrUL - EnrLU - this is a trick to 00447 * keep stable solution by effectively dividing up heating and cooling, 00448 * so that negative cooling does not occur */ 00449 00450 EnrLU = t10->Lo->Pop*CollRate[0][1]*t10->EnergyErg; 00451 EnrUL = t10->Hi->Pop*CollRate[1][0]*t10->EnergyErg; 00452 /* energy exchange due to this level 00453 * net cooling due to excitation minus part of de-excit */ 00454 t10->Coll.cool = EnrLU - EnrUL*t10->Emis->ColOvTot; 00455 /* net heating is remainder */ 00456 t10->Coll.heat = EnrUL*(1. - t10->Emis->ColOvTot); 00457 /* add to cooling */ 00458 CoolAdd( chLabel, t10->WLAng , t10->Coll.cool); 00459 /* derivative of cooling function */ 00460 thermal.dCooldT += t10->Coll.cool * (t10->EnergyK * thermal.tsq1 - thermal.halfte ); 00461 00462 EnrLU = t20->Lo->Pop*CollRate[0][2]*t20->EnergyErg; 00463 EnrUL = t20->Hi->Pop*CollRate[2][0]*t20->EnergyErg; 00464 t20->Coll.cool = EnrLU - EnrUL*t20->Emis->ColOvTot; 00465 t20->Coll.heat = EnrUL*(1. - t20->Emis->ColOvTot); 00466 /* add to cooling */ 00467 CoolAdd( chLabel, t20->WLAng , t20->Coll.cool); 00468 thermal.dCooldT += t20->Coll.cool * (t20->EnergyK * thermal.tsq1 - thermal.halfte ); 00469 00470 EnrLU = t30->Lo->Pop*CollRate[0][3]*t30->EnergyErg; 00471 EnrUL = t30->Hi->Pop*CollRate[3][0]*t30->EnergyErg; 00472 t30->Coll.cool = EnrLU - EnrUL*t30->Emis->ColOvTot; 00473 t30->Coll.heat = EnrUL*(1. - t30->Emis->ColOvTot); 00474 /* add to cooling */ 00475 CoolAdd( chLabel, t30->WLAng , t30->Coll.cool); 00476 thermal.dCooldT += t30->Coll.cool * (t30->EnergyK * thermal.tsq1 - thermal.halfte ); 00477 00478 EnrLU = t21->Lo->Pop*CollRate[1][2]*t21->EnergyErg; 00479 EnrUL = t21->Hi->Pop*CollRate[2][1]*t21->EnergyErg; 00480 t21->Coll.cool = EnrLU - EnrUL*t21->Emis->ColOvTot; 00481 t21->Coll.heat = EnrUL*(1. - t21->Emis->ColOvTot); 00482 /* add to cooling */ 00483 CoolAdd( chLabel, t21->WLAng , t21->Coll.cool); 00484 /* use of 20 in intentional in following - that is Boltzmann factor */ 00485 thermal.dCooldT += t21->Coll.cool * (t20->EnergyK * thermal.tsq1 - thermal.halfte ); 00486 00487 EnrLU = t31->Lo->Pop*CollRate[1][3]*t31->EnergyErg; 00488 EnrUL = t31->Hi->Pop*CollRate[3][1]*t31->EnergyErg; 00489 t31->Coll.cool = EnrLU - EnrUL*t31->Emis->ColOvTot; 00490 t31->Coll.heat = EnrUL*(1. - t31->Emis->ColOvTot); 00491 /* add to cooling */ 00492 CoolAdd( chLabel, t31->WLAng , t31->Coll.cool); 00493 /* use of 30 in intentional in following - that is Boltzmann factor */ 00494 thermal.dCooldT += t31->Coll.cool * (t30->EnergyK * thermal.tsq1 - thermal.halfte ); 00495 00496 EnrLU = t41->Lo->Pop*CollRate[1][4]*t41->EnergyErg; 00497 EnrUL = t41->Hi->Pop*CollRate[4][1]*t41->EnergyErg; 00498 t41->Coll.cool = EnrLU - EnrUL*t41->Emis->ColOvTot; 00499 t41->Coll.heat = EnrUL*(1. - t41->Emis->ColOvTot); 00500 /* add to cooling */ 00501 CoolAdd( chLabel, t41->WLAng , t41->Coll.cool); 00502 /* use of 31 in intentional in following - that is Boltzmann factor (no 40 here) */ 00503 thermal.dCooldT += t41->Coll.cool * (t41->EnergyK * thermal.tsq1 - thermal.halfte ); 00504 } 00505 00506 return; 00507 }