00001
00002
00003
00004
00005
00006
00007
00008
00009
00010 #include "cddefines.h"
00011 #include "physconst.h"
00012 #include "dense.h"
00013 #include "coolheavy.h"
00014 #include "taulines.h"
00015 #include "phycon.h"
00016 #include "iso.h"
00017 #include "conv.h"
00018 #include "cooling.h"
00019 #include "trace.h"
00020 #include "hydrogenic.h"
00021 #include "ligbar.h"
00022 #include "cooling.h"
00023 #include "thermal.h"
00024 #include "lines_service.h"
00025 #include "atoms.h"
00026 #include "atomfeii.h"
00027 #include "fe.h"
00028
00029
00030 STATIC void Fe3Lev14(void);
00031
00032
00033 STATIC void Fe4Lev12(void);
00034
00035
00036 STATIC void Fe2_cooling( void )
00037 {
00038 const int NLFE2 = 6;
00039 long int i, j;
00040 int nNegPop;
00041
00042 static double **AulPump,
00043 **CollRate,
00044 **AulEscp,
00045 **col_str ,
00046 **AulDest,
00047 *depart,
00048 *pops,
00049 *destroy,
00050 *create;
00051
00052 static bool lgFirst=true;
00053 bool lgZeroPop;
00054
00055
00056 static double gFe2[NLFE2]={1.,1.,1.,1.,1.,1.};
00057
00058 static double ex[NLFE2]={0.,3.32e4,5.68e4,6.95e4,1.15e5,1.31e5};
00059
00060
00061
00062 static double TUsed = 0.;
00063 static double AbunUsed = 0.;
00064
00065 static long int nZUsed=-1,
00066
00067 nCall=0;
00068
00069 DEBUG_ENTRY( "Fe2_cooling()" );
00070
00071
00072 if( dense.xIonDense[ipIRON][1] == 0. )
00073 {
00074
00075
00076 FeII.Fe2_large_cool = 0.;
00077 FeII.Fe2_large_heat = 0.;
00078 FeII.ddT_Fe2_large_cool = 0.;
00079
00080
00081 FeII.Fe2_UVsimp_cool = 0.;
00082 FeII.ddT_Fe2_UVsimp_cool = 0.;
00083
00084
00085 FeIIIntenZero();
00086 return;
00087 }
00088
00089
00090
00091
00092
00093
00094
00095
00096 if( FeII.lgSlow ||
00097 conv.lgFirstSweepThisZone || conv.lgLastSweepThisZone ||
00098
00099 ( !fp_equal( phycon.te, TUsed ) && fabs(FeII.Fe2_large_cool/thermal.ctot) > 0.002 &&
00100 fabs(dense.xIonDense[ipIRON][1]-AbunUsed)/SDIV(AbunUsed) > 0.002 ) ||
00101 ( !fp_equal( phycon.te, TUsed ) && fabs(FeII.Fe2_large_cool/thermal.ctot) > 0.01) )
00102 {
00103
00104 if( conv.lgFirstSweepThisZone )
00105
00106 nCall = 0;
00107 else
00108
00109
00110 ++nCall;
00111
00112
00113 if( trace.nTrConvg >= 5 )
00114 {
00115 fprintf( ioQQQ, " CoolIron5 calling FeIILevelPops since ");
00116 if( conv.lgFirstSweepThisZone )
00117 {
00118 fprintf( ioQQQ,
00119 "first sweep this zone." );
00120 }
00121 else if( ! fp_equal( phycon.te, TUsed ) )
00122 {
00123 fprintf( ioQQQ,
00124 "temperature changed, old new are %g %g, nCall %li ",
00125 TUsed, phycon.te , nCall);
00126 }
00127 else if( nzone != nZUsed )
00128 {
00129 fprintf( ioQQQ,
00130 "new zone, nCall %li ", nCall );
00131 }
00132 else if( FeII.lgSlow )
00133 {
00134 fprintf( ioQQQ,
00135 "FeII.lgSlow set %li", nCall );
00136 }
00137 else if( conv.lgSearch )
00138 {
00139 fprintf( ioQQQ,
00140 " in search phase %li", nCall );
00141 }
00142 else if( nCall < 2 )
00143 {
00144 fprintf( ioQQQ,
00145 "not second nCall %li " , nCall );
00146 }
00147 else if( ! fp_equal( phycon.te, TUsed ) && FeII.Fe2_large_cool/thermal.ctot > 0.001 )
00148 {
00149 fprintf( ioQQQ,
00150 "temp or cooling changed, new are %g %g nCall %li ",
00151 phycon.te, FeII.Fe2_large_cool, nCall );
00152 }
00153 else
00154 {
00155 fprintf(ioQQQ, "????");
00156 }
00157 fprintf(ioQQQ, "\n");
00158 }
00159
00160
00161 TUsed = phycon.te;
00162 AbunUsed = dense.xIonDense[ipIRON][1];
00163 nZUsed = nzone;
00164
00165
00166 if( FeII.lgPrint )
00167 {
00168 fprintf(ioQQQ,
00169 " FeIILevelPops called zone %4li te %5f abun %10e c(fe/tot):%6f nCall %li\n",
00170 nzone,phycon.te,AbunUsed,FeII.Fe2_large_cool/thermal.ctot,nCall);
00171 }
00172
00173
00174
00175
00176
00177
00178
00179
00180 FeII_RT_Make();
00181 FeIILevelPops();
00182 {
00183 enum{DEBUG_LOC=false};
00184 if( DEBUG_LOC && iteration > 1 && nzone >=4 )
00185 {
00186 fprintf(ioQQQ,"DEBUG1\t%li\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\t%.3e\n",
00187 nzone,
00188 phycon.te,
00189 dense.gas_phase[ipHYDROGEN],
00190 dense.eden,
00191 FeII.Fe2_large_cool ,
00192 FeII.ddT_Fe2_large_cool ,
00193 FeII.Fe2_large_cool/dense.eden/dense.gas_phase[ipHYDROGEN] ,
00194 thermal.ctot );
00195 }
00196 }
00197
00198 if( trace.nTrConvg >= 5 || FeII.lgPrint)
00199 {
00200
00201 fprintf( ioQQQ, " FeIILevelPops5 returned cool=%.2e heat=%.2e derivative=%.2e\n",
00202 FeII.Fe2_large_cool,FeII.Fe2_large_heat ,FeII.ddT_Fe2_large_cool);
00203 }
00204
00205 }
00206 else if( ! fp_equal( dense.xIonDense[ipIRON][1], AbunUsed ) )
00207 {
00208 realnum ratio;
00209
00210
00211
00212 if( trace.nTrConvg >= 5 )
00213 {
00214 fprintf( ioQQQ,
00215 " CoolIron rescaling FeIILevelPops since small change, CFe2=%.2e CTOT=%.2e\n",
00216 FeII.Fe2_large_cool,thermal.ctot);
00217 }
00218 ratio = dense.xIonDense[ipIRON][1]/AbunUsed;
00219 FeII.Fe2_large_cool *= ratio;
00220 FeII.ddT_Fe2_large_cool *= ratio;
00221 FeII.Fe2_large_heat *= ratio;
00222 AbunUsed = dense.xIonDense[ipIRON][1];
00223 }
00224 else
00225 {
00226
00227 if( trace.nTrConvg >= 5 )
00228 {
00229 fprintf( ioQQQ, " CoolIron NOT calling FeIILevelPops\n");
00230 }
00231 }
00232
00233
00234
00235
00236 CoolAdd("Fe 2",0,MAX2(0.,FeII.Fe2_large_cool));
00237
00238
00239 thermal.heating[25][27] = MAX2(0.,FeII.Fe2_large_heat);
00240
00241
00242 if( FeII.Fe2_large_cool > 0. )
00243 {
00244
00245 thermal.dCooldT += 3.*FeII.ddT_Fe2_large_cool;
00246 }
00247
00248 if( trace.lgTrace && trace.lgCoolTr )
00249 {
00250 fprintf( ioQQQ, " Large FeII returns te, cooling, dc=%11.3e%11.3e%11.3e\n",
00251 phycon.te, FeII.Fe2_large_cool, FeII.ddT_Fe2_large_cool );
00252 }
00253
00254
00255 if( !FeII.lgFeIILargeOn )
00256 {
00257
00258
00259
00260
00261
00262
00263 if( lgFirst )
00264 {
00265
00266 lgFirst = false;
00267
00268 pops = (double *)MALLOC( sizeof(double)*(NLFE2) );
00269 create = (double *)MALLOC( sizeof(double)*(NLFE2) );
00270 destroy = (double *)MALLOC( sizeof(double)*(NLFE2) );
00271 depart = (double *)MALLOC( sizeof(double)*(NLFE2) );
00272
00273 AulPump = ((double **)MALLOC((NLFE2)*sizeof(double *)));
00274 CollRate = ((double **)MALLOC((NLFE2)*sizeof(double *)));
00275 AulDest = ((double **)MALLOC((NLFE2)*sizeof(double *)));
00276 AulEscp = ((double **)MALLOC((NLFE2)*sizeof(double *)));
00277 col_str = ((double **)MALLOC((NLFE2)*sizeof(double *)));
00278
00279 for( i=0; i < NLFE2; ++i )
00280 {
00281 AulPump[i] = ((double *)MALLOC((NLFE2)*sizeof(double )));
00282 CollRate[i] = ((double *)MALLOC((NLFE2)*sizeof(double )));
00283 AulDest[i] = ((double *)MALLOC((NLFE2)*sizeof(double )));
00284 AulEscp[i] = ((double *)MALLOC((NLFE2)*sizeof(double )));
00285 col_str[i] = ((double *)MALLOC((NLFE2)*sizeof(double )));
00286 }
00287 }
00288
00289
00290 for( i=0; i < NLFE2; i++ )
00291 {
00292 create[i] = 0.;
00293 destroy[i] = 0.;
00294 for( j=0; j < NLFE2; j++ )
00295 {
00296
00297 col_str[j][i] = 0.;
00298 AulEscp[j][i] = 0.;
00299 AulDest[j][i] = 0.;
00300 AulPump[j][i] = 0.;
00301 }
00302 }
00303
00304
00305 AulEscp[1][0] = 1.;
00306 AulEscp[2][0] = ( TauLines[ipTuv3].Emis().Pesc() + TauLines[ipTuv3].Emis().Pelec_esc())*TauLines[ipTuv3].Emis().Aul();
00307 AulDest[2][0] = TauLines[ipTuv3].Emis().Pdest()*TauLines[ipTuv3].Emis().Aul();
00308 AulPump[0][2] = TauLines[ipTuv3].Emis().pump();
00309
00310 AulEscp[5][0] = (TauLines[ipTFe16].Emis().Pesc() + TauLines[ipTFe16].Emis().Pelec_esc())*TauLines[ipTFe16].Emis().Aul();
00311 AulDest[5][0] = TauLines[ipTFe16].Emis().Pdest()*TauLines[ipTFe16].Emis().Aul();
00312
00313 AulPump[0][5] = TauLines[ipTFe16].Emis().pump();
00314
00315
00316 double PumpLyaFeII = iso_sp[ipH_LIKE][ipHYDROGEN].st[ipH2p].Pop()*
00317 iso_sp[ipH_LIKE][ipHYDROGEN].trans(ipH2p,ipH1s).Emis().Aul()*
00318 hydro.dstfe2lya/SDIV(dense.xIonDense[ipIRON][1]);
00319
00320 AulPump[0][5] += PumpLyaFeII;
00321
00322 AulEscp[2][1] = (TauLines[ipTr48].Emis().Pesc() + TauLines[ipTr48].Emis().Pelec_esc())*TauLines[ipTr48].Emis().Aul();
00323 AulDest[2][1] = TauLines[ipTr48].Emis().Pdest()*TauLines[ipTr48].Emis().Aul();
00324 AulPump[1][2] = TauLines[ipTr48].Emis().pump();
00325
00326 AulEscp[5][1] = (TauLines[ipTFe26].Emis().Pesc() + TauLines[ipTFe26].Emis().Pelec_esc())*TauLines[ipTFe26].Emis().Aul();
00327 AulDest[5][1] = TauLines[ipTFe26].Emis().Pdest()*TauLines[ipTFe26].Emis().Aul();
00328 AulPump[1][5] = TauLines[ipTFe26].Emis().pump();
00329
00330 AulEscp[3][2] = (TauLines[ipTFe34].Emis().Pesc() + TauLines[ipTFe34].Emis().Pelec_esc())*TauLines[ipTFe34].Emis().Aul();
00331 AulDest[3][2] = TauLines[ipTFe34].Emis().Pdest()*TauLines[ipTFe34].Emis().Aul();
00332 AulPump[2][3] = TauLines[ipTFe34].Emis().pump();
00333
00334 AulEscp[4][2] = (TauLines[ipTFe35].Emis().Pesc() + TauLines[ipTFe35].Emis().Pelec_esc())*TauLines[ipTFe35].Emis().Aul();
00335 AulDest[4][2] = TauLines[ipTFe35].Emis().Pdest()*TauLines[ipTFe35].Emis().Aul();
00336 AulPump[2][4] = TauLines[ipTFe35].Emis().pump();
00337
00338 AulEscp[5][3] = (TauLines[ipTFe46].Emis().Pesc() + TauLines[ipTFe46].Emis().Pelec_esc())*TauLines[ipTFe46].Emis().Aul();
00339 AulDest[5][3] = TauLines[ipTFe46].Emis().Pdest()*TauLines[ipTFe46].Emis().Aul();
00340 AulPump[3][5] = TauLines[ipTFe46].Emis().pump();
00341
00342 AulEscp[5][4] = (TauLines[ipTFe56].Emis().Pesc() + TauLines[ipTFe56].Emis().Pelec_esc())*TauLines[ipTFe56].Emis().Aul();
00343 AulDest[5][4] = TauLines[ipTFe56].Emis().Pdest()*TauLines[ipTFe56].Emis().Aul();
00344 AulPump[4][5] = TauLines[ipTFe56].Emis().pump();
00345
00346
00347 col_str[1][0] = 1.;
00348 col_str[2][0] = 12.;
00349 col_str[3][0] = 1.;
00350 col_str[4][0] = 1.;
00351 col_str[5][0] = 12.;
00352 col_str[2][1] = 6.;
00353 col_str[3][1] = 1.;
00354 col_str[4][1] = 1.;
00355 col_str[5][1] = 12.;
00356 col_str[3][2] = 6.;
00357 col_str[4][2] = 12.;
00358 col_str[5][2] = 1.;
00359 col_str[4][3] = 1.;
00360 col_str[5][3] = 12.;
00361 col_str[5][4] = 6.;
00362
00363
00364
00365 atom_levelN(NLFE2,
00366 dense.xIonDense[ipIRON][1],
00367 gFe2,
00368 ex,
00369 'K',
00370 pops,
00371 depart,
00372 &AulEscp ,
00373 &col_str,
00374 &AulDest,
00375 &AulPump,
00376 &CollRate,
00377 create,
00378 destroy,
00379 false,
00380
00381 &FeII.Fe2_UVsimp_cool,
00382 &FeII.ddT_Fe2_UVsimp_cool,
00383 "FeII",
00384 &nNegPop,
00385 &lgZeroPop,
00386 false );
00387
00388
00389 if( nNegPop > 0 )
00390 {
00391 fprintf(ioQQQ," PROBLEM, atom_levelN returned negative population for simple UV FeII.\n");
00392 }
00393
00394 ;
00395 CoolAdd("Fe 2",0,MAX2(0.,FeII.Fe2_UVsimp_cool));
00396 thermal.heating[25][27] = MAX2(0.,-FeII.Fe2_UVsimp_cool);
00397 thermal.dCooldT += FeII.ddT_Fe2_UVsimp_cool;
00398
00399
00400 ASSERT( NLFE2 <= LIMLEVELN );
00401 for( i=0; i < NLFE2; ++i )
00402 {
00403 atoms.PopLevels[i] = pops[i];
00404 atoms.DepLTELevels[i] = depart[i];
00405 }
00406
00407 (*TauLines[ipTuv3].Lo()).Pop() = pops[0];
00408 (*TauLines[ipTuv3].Hi()).Pop() = pops[2];
00409 TauLines[ipTuv3].Emis().PopOpc() = (pops[0] - pops[2]);
00410 TauLines[ipTuv3].Emis().phots() = pops[2]*AulEscp[2][0];
00411 TauLines[ipTuv3].Emis().xIntensity() =
00412 TauLines[ipTuv3].Emis().phots()*TauLines[ipTuv3].EnergyErg();
00413
00414 (*TauLines[ipTr48].Lo()).Pop() = pops[1];
00415 (*TauLines[ipTr48].Hi()).Pop() = pops[2];
00416 TauLines[ipTr48].Emis().PopOpc() = (pops[1] - pops[2]);
00417 TauLines[ipTr48].Emis().phots() = pops[2]*AulEscp[2][1];
00418 TauLines[ipTr48].Emis().xIntensity() =
00419 TauLines[ipTr48].Emis().phots()*TauLines[ipTr48].EnergyErg();
00420
00421 FeII.for7 = pops[1]*AulEscp[1][0]*4.65e-12;
00422
00423 (*TauLines[ipTFe16].Lo()).Pop() = pops[0];
00424 (*TauLines[ipTFe16].Hi()).Pop() = pops[5];
00425
00426
00427
00428
00429 TauLines[ipTFe16].Emis().PopOpc() = (pops[0] - pops[5]);
00430 TauLines[ipTFe16].Emis().phots() = pops[5]*AulEscp[5][0];
00431 TauLines[ipTFe16].Emis().xIntensity() =
00432 TauLines[ipTFe16].Emis().phots()*TauLines[ipTFe16].EnergyErg();
00433
00434 (*TauLines[ipTFe26].Lo()).Pop() = pops[1];
00435 (*TauLines[ipTFe26].Hi()).Pop() = pops[5];
00436 TauLines[ipTFe26].Emis().PopOpc() = (pops[1] - pops[5]);
00437 TauLines[ipTFe26].Emis().phots() = pops[5]*AulEscp[5][1];
00438 TauLines[ipTFe26].Emis().xIntensity() =
00439 TauLines[ipTFe26].Emis().phots()*TauLines[ipTFe26].EnergyErg();
00440
00441 (*TauLines[ipTFe34].Lo()).Pop() = pops[2];
00442 (*TauLines[ipTFe34].Hi()).Pop() = pops[3];
00443 TauLines[ipTFe34].Emis().PopOpc() = (pops[2] - pops[3]);
00444 TauLines[ipTFe34].Emis().phots() = pops[3]*AulEscp[3][2];
00445 TauLines[ipTFe34].Emis().xIntensity() =
00446 TauLines[ipTFe34].Emis().phots()*TauLines[ipTFe34].EnergyErg();
00447
00448 (*TauLines[ipTFe35].Lo()).Pop() = pops[2];
00449 (*TauLines[ipTFe35].Hi()).Pop() = pops[4];
00450 TauLines[ipTFe35].Emis().PopOpc() = (pops[2] - pops[4]);
00451 TauLines[ipTFe35].Emis().phots() = pops[4]*AulEscp[4][2];
00452 TauLines[ipTFe35].Emis().xIntensity() =
00453 TauLines[ipTFe35].Emis().phots()*TauLines[ipTFe35].EnergyErg();
00454
00455 (*TauLines[ipTFe46].Lo()).Pop() = pops[3];
00456 (*TauLines[ipTFe46].Hi()).Pop() = pops[5];
00457 TauLines[ipTFe46].Emis().PopOpc() = (pops[3] - pops[5]);
00458 TauLines[ipTFe46].Emis().phots() = pops[5]*AulEscp[5][3];
00459 TauLines[ipTFe46].Emis().xIntensity() =
00460 TauLines[ipTFe46].Emis().phots()*TauLines[ipTFe46].EnergyErg();
00461
00462 (*TauLines[ipTFe56].Lo()).Pop() = pops[4];
00463 (*TauLines[ipTFe56].Hi()).Pop() = pops[5];
00464 TauLines[ipTFe56].Emis().PopOpc() = (pops[4] - pops[5]);
00465 TauLines[ipTFe56].Emis().phots() = pops[5]*AulEscp[5][4];
00466 TauLines[ipTFe56].Emis().xIntensity() =
00467 TauLines[ipTFe56].Emis().phots()*TauLines[ipTFe56].EnergyErg();
00468
00469
00470
00471
00472 PutCS(10.,TauLines[ipT191]);
00473 atom_level2(TauLines[ipT191]);
00474 }
00475
00476 {
00477
00478 enum{DEBUG_LOC=false};
00479
00480 if( DEBUG_LOC && iteration > 1 && nzone >=4 )
00481 {
00482 fprintf(ioQQQ,"DEBUG2\t%.2e\t%.2e\t%.2e\n",
00483 phycon.te,
00484 FeII.Fe2_large_cool ,
00485 FeII.Fe2_UVsimp_cool );
00486 }
00487 }
00488
00489 return;
00490
00491 }
00492
00493
00494 void CoolIron(void)
00495 {
00496 realnum rate;
00497
00498 DEBUG_ENTRY( "CoolIron()" );
00499
00500
00501
00505
00506
00507 rate = (realnum)(1.2e-7 * dense.eden +
00508
00509
00510 8.0e-10*pow((phycon.te/100.), 0.17 )*dense.xIonDense[ipHYDROGEN][0]);
00511 LineConvRate2CS( TauLines[ipFe1_24m] , rate );
00512
00513 rate = (realnum)(9.3e-8 * dense.eden +
00514
00515
00516 5.3e-10*pow((phycon.te/100.), 0.17 )*dense.xIonDense[ipHYDROGEN][0]);
00517 LineConvRate2CS( TauLines[ipFe1_35m] , rate );
00518
00519 rate = (realnum)(1.2e-7 * dense.eden +
00520
00521
00522 6.9e-10*pow((phycon.te/100.), 0.17 )*dense.xIonDense[ipHYDROGEN][0]);
00523 (*(*TauDummy).Hi()).g() = (*TauLines[ipFe1_35m].Hi()).g();
00524 LineConvRate2CS( *TauDummy , rate );
00525
00526 (*(*TauDummy).Hi()).g() = 0.;
00527
00528 atom_level3(TauLines[ipFe1_24m],TauLines[ipFe1_35m],*TauDummy);
00529
00530
00531
00532
00533 MakeCS(TauLines[ipFeI3884]);
00534 atom_level2(TauLines[ipFeI3884]);
00535
00536
00537 MakeCS(TauLines[ipFeI3729]);
00538 atom_level2(TauLines[ipFeI3729]);
00539
00540
00541 MakeCS(TauLines[ipFeI3457]);
00542 atom_level2(TauLines[ipFeI3457]);
00543
00544
00545 MakeCS(TauLines[ipFeI3021]);
00546 atom_level2(TauLines[ipFeI3021]);
00547
00548
00549 MakeCS(TauLines[ipFeI2966]);
00550 atom_level2(TauLines[ipFeI2966]);
00551
00552
00553 Fe2_cooling();
00554
00555
00556
00557
00558
00559
00560
00561
00562
00563
00564
00565
00566
00567 PutCS(25.,TauLines[ipT1122]);
00568 atom_level2(TauLines[ipT1122]);
00569
00570
00571 Fe3Lev14();
00572
00573
00574 Fe4Lev12();
00575
00576
00577 CoolHeavy.c3892 = atom_pop2(7.4,25.,5.,0.6,3.7e4,dense.xIonDense[ipIRON][4])*
00578 5.11e-12;
00579 CoolAdd("Fe 5",3892,CoolHeavy.c3892);
00580
00581 return;
00582 }
00583
00584
00585 STATIC void Fe4Lev12(void)
00586 {
00587 const int NLFE4 = 12;
00588 bool lgZeroPop;
00589 int nNegPop;
00590 long int i,
00591 j;
00592 static bool lgFirst=true;
00593
00594 double dfe4dt;
00595
00596
00597 static double
00598 **AulEscp,
00599 **col_str,
00600 **AulDest,
00601 depart[NLFE4],
00602 pops[NLFE4],
00603 destroy[NLFE4],
00604 create[NLFE4],
00605 **CollRate,
00606 **AulPump;
00607
00608 static const double Fe4A[NLFE4][NLFE4] = {
00609 {0.,0.,0.,1.e-5,0.,1.368,.89,0.,1.3e-3,1.8e-4,.056,.028},
00610 {0.,0.,2.6e-8,0.,0.,0.,0.,0.,1.7e-7,0.,0.,0.},
00611 {0.,0.,0.,0.,3.5e-7,6.4e-10,0.,0.,6.315e-4,0.,6.7e-7,0.},
00612 {0.,0.,0.,0.,1.1e-6,6.8e-5,8.6e-6,3.4e-10,7.6e-5,1.e-7,5.8e-4,2.8e-4},
00613 {0.,0.,0.,0.,0.,1.5e-5,1.3e-9,0.,7.6e-4,0.,1.1e-6,6.0e-7},
00614 {0.,0.,0.,0.,0.,0.,1.1e-5,1.2e-13,.038,9.9e-7,.022,.018},
00615 {0.,0.,0.,0.,0.,0.,0.,3.7e-5,2.9e-6,.034,3.5e-3,.039},
00616 {0.,0.,0.,0.,0.,0.,0.,0.,0.,.058,3.1e-6,1.4e-3},
00617 {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,1.3e-4,3.1e-14},
00618 {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,1.9e-19,1.0e-5},
00619 {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,1.3e-7},
00620 {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.}
00621 };
00622
00623 static const double gfe4[NLFE4]={6.,12.,10.,6.,8.,6.,4.,2.,8.,2.,6.,4.};
00624
00625
00626
00627
00628
00629
00630 static const double excit_wn[NLFE4]={0.,32245.5,32292.8,32301.2,32305.7,35253.8,
00631 35333.3,35406.6,38779.4,38896.7,38935.1,38938.2};
00632
00633 DEBUG_ENTRY( "Fe4Lev12()" );
00634
00635 if( lgFirst )
00636 {
00637
00638 lgFirst = false;
00639
00640
00641 AulPump = ((double **)MALLOC((NLFE4)*sizeof(double *)));
00642 CollRate = ((double **)MALLOC((NLFE4)*sizeof(double *)));
00643 AulDest = ((double **)MALLOC((NLFE4)*sizeof(double *)));
00644 AulEscp = ((double **)MALLOC((NLFE4)*sizeof(double *)));
00645 col_str = ((double **)MALLOC((NLFE4)*sizeof(double *)));
00646 for( i=0; i < NLFE4; ++i )
00647 {
00648 AulPump[i] = ((double *)MALLOC((NLFE4)*sizeof(double )));
00649 CollRate[i] = ((double *)MALLOC((NLFE4)*sizeof(double )));
00650 AulDest[i] = ((double *)MALLOC((NLFE4)*sizeof(double )));
00651 AulEscp[i] = ((double *)MALLOC((NLFE4)*sizeof(double )));
00652 col_str[i] = ((double *)MALLOC((NLFE4)*sizeof(double )));
00653 }
00654 }
00655
00656
00657 if( dense.xIonDense[ipIRON][3] <= 0. )
00658 {
00659 fe.Fe4CoolTot = 0.;
00660 fe.fe40401 = 0.;
00661 fe.fe42836 = 0.;
00662 fe.fe42829 = 0.;
00663 fe.fe42567 = 0.;
00664 fe.fe41207 = 0.;
00665 fe.fe41206 = 0.;
00666 fe.fe41106 = 0.;
00667 fe.fe41007 = 0.;
00668 fe.fe41008 = 0.;
00669 fe.fe40906 = 0.;
00670 CoolAdd("Fe 4",0,0.);
00671
00672
00673
00674 ASSERT( NLFE4 <= LIMLEVELN);
00675 for( i=0; i < NLFE4; i++ )
00676 {
00677 atoms.PopLevels[i] = 0.;
00678 atoms.DepLTELevels[i] = 1.;
00679 }
00680 return;
00681 }
00682
00683
00684
00685
00686
00687
00688
00689
00690
00691
00692
00693
00694
00695
00696
00697
00698
00699 for( i=0; i < NLFE4; i++ )
00700 {
00701 create[i] = 0.;
00702 destroy[i] = 0.;
00703 for( j=0; j < NLFE4; j++ )
00704 {
00705 col_str[j][i] = 0.;
00706 AulEscp[j][i] = 0.;
00707 AulDest[j][i] = 0.;
00708 AulPump[j][i] = 0.;
00709 }
00710 }
00711
00712
00713 for( long ipHi=1; ipHi < NLFE4; ipHi++ )
00714 {
00715 for( long ipLo=0; ipLo < ipHi; ipLo++ )
00716 {
00717 AulEscp[ipHi][ipLo] = Fe4A[ipLo][ipHi];
00718 col_str[ipHi][ipLo] = Fe4_cs(ipLo, ipHi);
00719 }
00720 }
00721
00722
00723
00724
00725 atom_levelN(NLFE4,
00726 dense.xIonDense[ipIRON][3],
00727 gfe4,
00728 excit_wn,
00729 'w',
00730 pops,
00731 depart,
00732 &AulEscp ,
00733 &col_str ,
00734 &AulDest,
00735 &AulPump,
00736 &CollRate,
00737 create,
00738 destroy,
00739
00740 false,
00741 &fe.Fe4CoolTot,
00742 &dfe4dt,
00743 "FeIV",
00744
00745 &nNegPop,
00746 &lgZeroPop,
00747 false );
00748
00749
00750 ASSERT( NLFE4 <= LIMLEVELN );
00751 for( i=0; i < NLFE4; ++i )
00752 {
00753 atoms.PopLevels[i] = pops[i];
00754 atoms.DepLTELevels[i] = depart[i];
00755 }
00756
00757 if( nNegPop > 0 )
00758 {
00759 fprintf( ioQQQ, " fe4levl2 found negative populations\n" );
00760 ShowMe();
00761 cdEXIT(EXIT_FAILURE);
00762 }
00763
00764 CoolAdd("Fe 4",0,fe.Fe4CoolTot);
00765
00766 thermal.dCooldT += dfe4dt;
00767
00768
00769
00770 fe.fe40401 = (pops[3]*Fe4A[0][3]*(excit_wn[3] - excit_wn[0]) +
00771 pops[4]*Fe4A[0][4]*(excit_wn[4] - excit_wn[0]) )*T1CM*BOLTZMANN;
00772
00773 fe.fe42836 = pops[5]*Fe4A[0][5]*(excit_wn[5] - excit_wn[0])*T1CM*BOLTZMANN;
00774
00775 fe.fe42829 = pops[6]*Fe4A[0][6]*(excit_wn[5] - excit_wn[0])*T1CM*BOLTZMANN;
00776
00777 fe.fe42567 = (pops[10]*Fe4A[0][10]*(excit_wn[10] - excit_wn[0]) +
00778 pops[11]*Fe4A[0][11]*(excit_wn[10] - excit_wn[0]))*T1CM*BOLTZMANN;
00779
00780 fe.fe41207 = pops[11]*Fe4A[6][11]*(excit_wn[11] - excit_wn[6])*T1CM*BOLTZMANN;
00781 fe.fe41206 = pops[11]*Fe4A[5][11]*(excit_wn[11] - excit_wn[5])*T1CM*BOLTZMANN;
00782 fe.fe41106 = pops[10]*Fe4A[5][10]*(excit_wn[10] - excit_wn[5])*T1CM*BOLTZMANN;
00783 fe.fe41007 = pops[9]*Fe4A[6][9]*(excit_wn[9] - excit_wn[6])*T1CM*BOLTZMANN;
00784 fe.fe41008 = pops[9]*Fe4A[7][9]*(excit_wn[9] - excit_wn[7])*T1CM*BOLTZMANN;
00785 fe.fe40906 = pops[8]*Fe4A[5][8]*(excit_wn[8] - excit_wn[5])*T1CM*BOLTZMANN;
00786 return;
00787 }
00788
00789
00790
00791
00792
00793 STATIC void Fe3Lev14(void)
00794 {
00795 bool lgZeroPop;
00796 int nNegPop;
00797 long int i,
00798 j;
00799 static bool lgFirst=true;
00800
00801 double dfe3dt;
00802
00803 long int ihi , ilo;
00804 static double
00805 *depart,
00806 *pops,
00807 *destroy,
00808 *create ,
00809 **AulDest,
00810 **CollRate,
00811 **AulPump,
00812 **AulNet,
00813 **col_str;
00814
00815
00816 static double gfe3[NLFE3]={9.,7.,5.,3.,1.,5.,13.,11.,9.,3.,1.,9.,7.,5.};
00817
00818
00819
00820
00821
00822 static double excit_wn[NLFE3]={
00823 0.0 , 436.2, 738.9, 932.4, 1027.3,
00824 19404.8, 20051.1, 20300.8, 20481.9, 20688.4,
00825 21208.5, 21462.2, 21699.9, 21857.2 };
00826
00827 DEBUG_ENTRY( "Fe3Lev14()" );
00828
00829 if( lgFirst )
00830 {
00831
00832 lgFirst = false;
00833
00834
00835 depart = ((double *)MALLOC((NLFE3)*sizeof(double)));
00836 pops = ((double *)MALLOC((NLFE3)*sizeof(double)));
00837 destroy = ((double *)MALLOC((NLFE3)*sizeof(double)));
00838 create = ((double *)MALLOC((NLFE3)*sizeof(double)));
00839
00840 fe.Fe3_wl = ((double **)MALLOC((NLFE3)*sizeof(double *)));
00841 fe.Fe3_emiss = ((double **)MALLOC((NLFE3)*sizeof(double *)));
00842 AulNet = ((double **)MALLOC((NLFE3)*sizeof(double *)));
00843 col_str = ((double **)MALLOC((NLFE3)*sizeof(double *)));
00844 AulPump = ((double **)MALLOC((NLFE3)*sizeof(double *)));
00845 CollRate = ((double **)MALLOC((NLFE3)*sizeof(double *)));
00846 AulDest = ((double **)MALLOC((NLFE3)*sizeof(double *)));
00847 for( i=0; i < NLFE3; ++i )
00848 {
00849 fe.Fe3_wl[i] = ((double *)MALLOC((NLFE3)*sizeof(double )));
00850 fe.Fe3_emiss[i] = ((double *)MALLOC((NLFE3)*sizeof(double )));
00851 AulNet[i] = ((double *)MALLOC((NLFE3)*sizeof(double )));
00852 col_str[i] = ((double *)MALLOC((NLFE3)*sizeof(double )));
00853 AulPump[i] = ((double *)MALLOC((NLFE3)*sizeof(double )));
00854 CollRate[i] = ((double *)MALLOC((NLFE3)*sizeof(double )));
00855 AulDest[i] = ((double *)MALLOC((NLFE3)*sizeof(double )));
00856 }
00857
00858
00859 for( i=0; i < NLFE3; ++i )
00860 {
00861 create[i] = 0.;
00862 destroy[i] = 0.;
00863 for( j=0; j < NLFE3; ++j )
00864 {
00865 AulNet[i][j] = 0.;
00866 col_str[i][j] = 0.;
00867 CollRate[i][j] = 0.;
00868 AulDest[i][j] = 0.;
00869 AulPump[i][j] = 0.;
00870 fe.Fe3_wl[i][j] = 0.;
00871 fe.Fe3_emiss[i][j] = 0.;
00872 }
00873 }
00874
00875
00876 for( ihi=1; ihi < NLFE3; ++ihi )
00877 {
00878 for( ilo=0; ilo < ihi; ++ilo )
00879 {
00880 fe.Fe3_wl[ihi][ilo] = 1e8/(excit_wn[ihi]-excit_wn[ilo]) /
00881 RefIndex( (excit_wn[ihi]-excit_wn[ilo]) );
00882 }
00883 }
00884
00885
00886
00887 AulNet[1][0] = 2.8e-3;
00888 AulNet[7][0] = 4.9e-6;
00889 AulNet[8][0] = 5.7e-3;
00890 AulNet[11][0] = 4.5e-1;
00891 AulNet[12][0] = 4.2e-2;
00892
00893 AulNet[2][1] = 1.8e-3;
00894 AulNet[5][1] = 4.2e-1;
00895 AulNet[8][1] = 1.0e-3;
00896 AulNet[11][1] = 8.4e-2;
00897 AulNet[12][1] = 2.5e-1;
00898 AulNet[13][1] = 2.7e-2;
00899
00900 AulNet[3][2] = 7.0e-4;
00901 AulNet[5][2] = 5.1e-5;
00902 AulNet[9][2] = 5.4e-1;
00903 AulNet[12][2] = 8.5e-2;
00904 AulNet[13][2] = 9.8e-2;
00905
00906 AulNet[4][3] = 1.4e-4;
00907 AulNet[5][3] = 3.9e-2;
00908 AulNet[9][3] = 4.1e-5;
00909 AulNet[10][3] = 7.0e-1;
00910 AulNet[13][3] = 4.7e-2;
00911
00912 AulNet[9][4] = 9.3e-2;
00913
00914 AulNet[9][5] = 4.7e-2;
00915 AulNet[12][5] = 2.5e-6;
00916 AulNet[13][5] = 1.7e-5;
00917
00918 AulNet[7][6] = 2.7e-4;
00919
00920 AulNet[8][7] = 1.2e-4;
00921 AulNet[11][7] = 6.6e-4;
00922
00923 AulNet[11][8] = 1.6e-3;
00924 AulNet[12][8] = 7.8e-4;
00925
00926 AulNet[10][9] = 8.4e-3;
00927 AulNet[13][9] = 2.8e-7;
00928
00929 AulNet[12][11] = 3.0e-4;
00930
00931 AulNet[13][12] = 1.4e-4;
00932
00933 for( int ipHi = 1; ipHi < NLFE3; ipHi++)
00934 {
00935 for( int ipLo = 0; ipLo < ipHi; ipLo++)
00936 {
00937 col_str[ipHi][ipLo] = Fe3_cs(ipLo,ipHi);
00938 }
00939 }
00940 }
00941
00942
00943 if( dense.xIonDense[ipIRON][2] <= 0. )
00944 {
00945 CoolAdd("Fe 3",0,0.);
00946
00947 fe.Fe3CoolTot = 0.;
00948 for( ihi=1; ihi < NLFE3; ++ihi )
00949 {
00950 for( ilo=0; ilo < ihi; ++ilo )
00951 {
00952 fe.Fe3_emiss[ihi][ilo] = 0.;
00953 }
00954 }
00955
00956
00957 ASSERT( NLFE3 <= LIMLEVELN);
00958 for( i=0; i < NLFE3; i++ )
00959 {
00960 atoms.PopLevels[i] = 0.;
00961 atoms.DepLTELevels[i] = 1.;
00962 }
00963 return;
00964 }
00965
00966
00967 atom_levelN(
00968
00969 NLFE3,
00970
00971 dense.xIonDense[ipIRON][2],
00972
00973 gfe3,
00974
00975 excit_wn,
00976 'w',
00977
00978 pops,
00979
00980 depart,
00981
00982 &AulNet ,
00983
00984 &col_str ,
00985
00986 &AulDest,
00987
00988 &AulPump,
00989
00990 &CollRate,
00991
00992 create,
00993
00994 destroy,
00995
00996 false,
00997 &fe.Fe3CoolTot,
00998 &dfe3dt,
00999 "Fe 3",
01000 &nNegPop,
01001 &lgZeroPop,
01002 false );
01003
01004
01005 ASSERT( NLFE3 <= LIMLEVELN );
01006 for( i=0; i < NLFE3; ++i )
01007 {
01008 atoms.PopLevels[i] = pops[i];
01009 atoms.DepLTELevels[i] = depart[i];
01010 }
01011
01012 if( nNegPop > 0 )
01013 {
01014 fprintf( ioQQQ, " Fe3Lev14 found negative populations\n" );
01015 ShowMe();
01016 cdEXIT(EXIT_FAILURE);
01017 }
01018
01019
01020 CoolAdd("Fe 3",0,fe.Fe3CoolTot);
01021
01022 thermal.dCooldT += dfe3dt;
01023
01024
01025 for( ihi=1; ihi < NLFE3; ++ihi )
01026 {
01027 for( ilo=0; ilo < ihi; ++ilo )
01028 {
01029
01030 fe.Fe3_emiss[ihi][ilo] = pops[ihi]*AulNet[ihi][ilo]*(excit_wn[ihi] - excit_wn[ilo])*T1CM*BOLTZMANN;
01031 }
01032 }
01033 return;
01034 }
01035
01036 double Fe3_cs(long ipLo, long ipHi)
01037 {
01038 static double col_str[14][14];
01039
01040 static double lgOneTimeMustInit=true;
01041 if( lgOneTimeMustInit )
01042 {
01043 lgOneTimeMustInit = false;
01044
01046
01047 col_str[1][0] = 2.92;
01048 col_str[2][0] = 1.24;
01049 col_str[3][0] = 0.595;
01050 col_str[4][0] = 0.180;
01051 col_str[5][0] = 0.580;
01052 col_str[6][0] = 1.34;
01053 col_str[7][0] = 0.489;
01054 col_str[8][0] = 0.0926;
01055 col_str[9][0] = 0.165;
01056 col_str[10][0] = 0.0213;
01057 col_str[11][0] = 1.07;
01058 col_str[12][0] = 0.435;
01059 col_str[13][0] = 0.157;
01060
01061 col_str[2][1] = 2.06;
01062 col_str[3][1] = 0.799;
01063 col_str[4][1] = 0.225;
01064 col_str[5][1] = 0.335;
01065 col_str[6][1] = 0.555;
01066 col_str[7][1] = 0.609;
01067 col_str[8][1] = 0.367;
01068 col_str[9][1] = 0.195;
01069 col_str[10][1] = 0.0698;
01070 col_str[11][1] = 0.538;
01071 col_str[12][1] = 0.484;
01072 col_str[13][1] = 0.285;
01073
01074 col_str[3][2] = 1.29;
01075 col_str[4][2] = 0.312;
01076 col_str[5][2] = 0.173;
01077 col_str[6][2] = 0.178;
01078 col_str[7][2] = 0.430;
01079 col_str[8][2] = 0.486;
01080 col_str[9][2] = 0.179;
01081 col_str[10][2] = 0.0741;
01082 col_str[11][2] = 0.249;
01083 col_str[12][2] = 0.362;
01084 col_str[13][2] = 0.324;
01085
01086 col_str[4][3] = 0.493;
01087 col_str[5][3] = 0.0767;
01088 col_str[6][3] = 0.0348;
01089 col_str[7][3] = 0.223;
01090 col_str[8][3] = 0.401;
01091 col_str[9][3] = 0.126;
01092 col_str[10][3] = 0.0528;
01093 col_str[11][3] = 0.101;
01094 col_str[12][3] = 0.207;
01095 col_str[13][3] = 0.253;
01096
01097 col_str[5][4] = 0.0211;
01098 col_str[6][4] = 0.00122;
01099 col_str[7][4] = 0.0653;
01100 col_str[8][4] = 0.154;
01101 col_str[9][4] = 0.0453;
01102 col_str[10][4] = 0.0189;
01103 col_str[11][4] = 0.0265;
01104 col_str[12][4] = 0.0654;
01105 col_str[13][4] = 0.0950;
01106
01107 col_str[6][5] = 0.403;
01108 col_str[7][5] = 0.213;
01109 col_str[8][5] = 0.0939;
01110 col_str[9][5] = 1.10;
01111 col_str[10][5] = 0.282;
01112 col_str[11][5] = 0.942;
01113 col_str[12][5] = 0.768;
01114 col_str[13][5] = 0.579;
01115
01116 col_str[7][6] = 2.84;
01117 col_str[8][6] = 0.379;
01118 col_str[9][6] = 0.0876;
01119 col_str[10][6] = 0.00807;
01120 col_str[11][6] = 1.85;
01121 col_str[12][6] = 0.667;
01122 col_str[13][6] = 0.0905;
01123
01124 col_str[8][7] = 3.07;
01125 col_str[9][7] = 0.167;
01126 col_str[10][7] = 0.0526;
01127 col_str[11][7] = 0.814;
01128 col_str[12][7] = 0.837;
01129 col_str[13][7] = 0.626;
01130
01131 col_str[9][8] = 0.181;
01132 col_str[10][8] = 0.0854;
01133 col_str[11][8] = 0.180;
01134 col_str[12][8] = 0.778;
01135 col_str[13][8] = 0.941;
01136
01137 col_str[10][9] = 0.377;
01138 col_str[11][9] = 0.603;
01139 col_str[12][9] = 0.472;
01140 col_str[13][9] = 0.302;
01141
01142 col_str[11][10] = 0.216;
01143 col_str[12][10] = 0.137;
01144 col_str[13][10] = 0.106;
01145
01146 col_str[12][11] = 1.25;
01147 col_str[13][11] = 0.292;
01148
01149 col_str[13][12] = 1.10;
01150 }
01151
01152 ASSERT( ipHi > ipLo );
01153 double CollisionStrength = col_str[ipHi][ipLo];
01154 ASSERT( CollisionStrength >0. );
01155
01156 return( CollisionStrength );
01157 }
01158
01159 double Fe4_cs(long ipLo, long ipHi)
01160 {
01161 const int NLFE4=12;
01162
01163
01164
01165 static const double Fe4CS[NLFE4][NLFE4] = {
01166 {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.},
01167 {0.98,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.},
01168 {0.8167,3.72,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.},
01169 {0.49,0.0475,0.330,0.,0.,0.,0.,0.,0.,0.,0.,0.},
01170 {0.6533,0.473,2.26,1.64,0.,0.,0.,0.,0.,0.,0.,0.},
01171 {0.45,0.686,0.446,0.106,0.254,0.,0.,0.,0.,0.,0.,0.},
01172 {0.30,0.392,0.152,0.269,0.199,0.605,0.,0.,0.,0.,0.,0.},
01173 {0.15,0.0207,0.190,0.0857,0.166,0.195,0.327,0.,0.,0.,0.,0.},
01174 {0.512,1.23,0.733,0.174,0.398,0.623,0.335,0.102,0.,0.,0.,0.},
01175 {0.128,0.0583,0.185,0.200,0.188,0.0835,0.127,0.0498,0.0787,0.,0.,0.},
01176 {0.384,0.578,0.534,0.363,0.417,0.396,0.210,0.171,0.810,0.101,0.,0.},
01177 {0.256,0.234,0.306,0.318,0.403,0.209,0.195,0.112,0.195,0.458,0.727,0.}
01178 };
01179
01180 ASSERT( ipHi > ipLo );
01181 double CollisionStrength = Fe4CS[ipHi][ipLo];
01182 ASSERT( CollisionStrength >0. );
01183
01184 return( CollisionStrength );
01185 }
01186
01187 double Fe5_cs(long ipLo, long ipHi)
01188 {
01189 const int NLFE5 = 14;
01190 static double col_str[NLFE5][NLFE5];
01191
01192
01193
01194 static double lgOneTimeMustInit=true;
01195 if( lgOneTimeMustInit )
01196 {
01197 lgOneTimeMustInit = false;
01198 for( int i = 0;i < NLFE5;i++)
01199 {
01200 for( int j = 0;j < NLFE5;j++)
01201 {
01202 col_str[i][j] = 1.;
01203 }
01204 }
01205
01206 col_str[10][3] = 1.4;
01207 col_str[7][2] = 1.1;
01208 col_str[13][4] = 3.7;
01209 col_str[12][3] = 3.7;
01210 col_str[11][2] = 2.0;
01211 }
01212
01213 ASSERT( ipHi > ipLo );
01214 double CollisionStrength = col_str[ipHi][ipLo];
01215 ASSERT( CollisionStrength >0. );
01216
01217 return( CollisionStrength );
01218 }