00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #include "cddefines.h"
00012 #include "taulines.h"
00013 #include "opacity.h"
00014 #include "dense.h"
00015 #include "iso.h"
00016 #include "hmi.h"
00017 #include "h2.h"
00018 #include "rfield.h"
00019 #include "conv.h"
00020 #include "rt.h"
00021 #include "atomfeii.h"
00022 #include "heavy.h"
00023 #include "he.h"
00024 #include "trace.h"
00025
00026
00027 static int nOTS_Line_type = -1;
00028 static int nOTS1=-1 , nOTS2=-1;
00029
00030 STATIC void RT_OTS_AddCont(
00031
00032 realnum ots,
00033
00034 long int ip);
00035
00036
00037
00038 void RT_OTS(void)
00039 {
00040 long int
00041 ipla,
00042 ipISO ,
00043 nelem,
00044 n;
00045 realnum
00046 difflya,
00047 esc,
00048 ots;
00049
00050
00051
00052 const realnum BOWEN = 0.5f;
00053 long int ipHi,
00054 ipLo;
00055
00056 double bwnfac;
00057 double ots660;
00058 realnum cont_phot_destroyed;
00059 double save_lya_dest,
00060 save_he2lya_dest;
00061
00062 double save_he2rec_dest;
00063
00064
00065
00066
00067
00068 DEBUG_ENTRY( "RT_OTS()" );
00069
00070 for( long i=0; i < rfield.nflux; i++ )
00071 {
00072 rfield.otslin[i] = 0.;
00073 rfield.otscon[i] = 0.;
00074 }
00075
00076
00077
00078
00079
00080
00081 nOTS_Line_type = 0;
00082 nelem = ipHELIUM;
00083 if( dense.lgElmtOn[nelem] )
00084 {
00085
00086
00087 bwnfac = BOWEN * MAX2(0.f,1.f- Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->Pesc -
00088 Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->Pelec_esc );
00089
00090
00091
00092 ots660 = Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->Aul*
00093 StatesElemNEW[nelem][nelem-ipH_LIKE][ipH2p].Pop*
00094
00095 Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->Pdest *BOWEN*2.0;
00096
00097
00098 if( ots660 > SMALLFLOAT )
00099 RT_OTS_AddLine(ots660 , he.ip660 );
00100
00101
00102
00103 Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->Pdest *= (realnum)bwnfac;
00104 ASSERT( Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->Pdest >= 0. );
00105 {
00106
00107
00108 enum {DEBUG_LOC=false};
00109 if( DEBUG_LOC )
00110 {
00111 fprintf(ioQQQ,"DEBUG HeII Bowen nzone %li bwnfac:%.2e bwnfac esc:%.2e ots660 %.2e\n",
00112 nzone,
00113 bwnfac ,
00114 bwnfac/BOWEN ,
00115 ots660 );
00116 }
00117 }
00118 }
00119
00120 else
00121 {
00122 bwnfac = 1.;
00123 }
00124
00125
00126 save_lya_dest = Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->Pdest;
00127
00128
00129
00130
00131 Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->Pdest *= rfield.lgLyaOTS;
00132
00133
00134 save_he2lya_dest = Transitions[ipH_LIKE][ipHELIUM][ipH2p][ipH1s].Emis->Pdest;
00135 Transitions[ipH_LIKE][ipHELIUM][ipH2p][ipH1s].Emis->Pdest *= rfield.lgHeIIOTS;
00136
00137
00138 save_he2rec_dest = iso.RadRecomb[ipH_LIKE][ipHELIUM][ipH1s][ipRecRad];
00139 iso.RadRecomb[ipH_LIKE][ipHELIUM][ipH1s][ipRecRad] *= rfield.lgHeIIOTS;
00140
00141 nOTS_Line_type = 1;
00142
00143
00144
00145
00146 for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
00147 {
00148 for( nelem=ipISO; nelem < LIMELM; nelem++ )
00149 {
00150 nOTS1 = ipISO;
00151 nOTS2 = nelem;
00152
00156 if( (dense.IonHigh[nelem] >= nelem+1-ipISO ) )
00157 {
00158
00159
00160
00161
00162 for( ipHi=1; ipHi < iso.numLevels_local[ipISO][nelem]; ipHi++ )
00163 {
00164 for( ipLo=0; ipLo < ipHi; ipLo++ )
00165 {
00166
00167
00168
00169 if( Transitions[ipISO][nelem][ipHi][ipLo].ipCont < 1 ||
00170 Transitions[ipISO][nelem][ipHi][ipLo].Emis->Pdest<= DEST0 )
00171 continue;
00172
00173
00174 Transitions[ipISO][nelem][ipHi][ipLo].Emis->ots =
00175 Transitions[ipISO][nelem][ipHi][ipLo].Emis->Aul*
00176 StatesElemNEW[nelem][nelem-ipISO][ipHi].Pop*
00177 Transitions[ipISO][nelem][ipHi][ipLo].Emis->Pdest;
00178
00179 ASSERT( Transitions[ipISO][nelem][ipHi][ipLo].Emis->ots >= 0. );
00180
00181
00182
00183
00184
00185
00186
00187 if( Transitions[ipISO][nelem][ipHi][ipLo].Emis->ots > SMALLFLOAT )
00188 RT_OTS_AddLine(Transitions[ipISO][nelem][ipHi][ipLo].Emis->ots,
00189 Transitions[ipISO][nelem][ipHi][ipLo].ipCont );
00190 }
00191 }
00192 {
00193
00194
00195
00196 enum {DEBUG_LOC=false};
00197
00198 if( DEBUG_LOC )
00199 {
00200 long ip;
00201 if( ipISO==0 && nelem==0 && nzone>500 )
00202 {
00203 ipHi = 2;
00204 ipLo = 0;
00205 ip = Transitions[ipISO][nelem][ipHi][ipLo].ipCont;
00206 fprintf(ioQQQ,"DEBUG hlyaots\t%.2f\tEdenTrue\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\n",
00207 fnzone,
00208 dense.EdenTrue ,
00209 Transitions[ipISO][nelem][ipHi][ipLo].Emis->ots,
00210 opac.opacity_abs[ip-1],
00211 StatesElemNEW[nelem][nelem-ipISO][ipHi].Pop,
00212 dense.xIonDense[nelem][nelem+1-ipISO],
00213 Transitions[ipISO][nelem][ipHi][ipLo].Emis->Pdest,
00214 rfield.otslin[ip-1]);
00215 }
00216 }
00217 }
00218
00219
00220
00221
00222
00223
00224
00225
00226
00227 for( n=0; n < iso.numLevels_local[ipISO][nelem]; n++ )
00228 {
00229 cont_phot_destroyed = (realnum)(iso.RadRecomb[ipISO][nelem][n][ipRecRad]*
00230 (1. - iso.RadRecomb[ipISO][nelem][n][ipRecEsc])*dense.eden*
00231 dense.xIonDense[nelem][nelem+1-ipISO]);
00232 ASSERT( cont_phot_destroyed >= 0. );
00233
00234
00235 RT_OTS_AddCont(cont_phot_destroyed,iso.ipIsoLevNIonCon[ipISO][nelem][n]);
00236
00237 {
00238
00239 enum {DEBUG_LOC=false};
00240
00241 if( DEBUG_LOC && nzone > 400 && nelem==0 && n==2 )
00242 {
00243 long ip = iso.ipIsoLevNIonCon[ipISO][nelem][n]-1;
00244 fprintf(ioQQQ,"hotsdebugg\t%.3f\t%li\th con ots\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\n",
00245 fnzone,
00246 n ,
00247 StatesElemNEW[ipHYDROGEN][ipHYDROGEN-ipISO][2].Pop,
00248 cont_phot_destroyed,
00249 cont_phot_destroyed/opac.opacity_abs[ip],
00250 rfield.otscon[ip] ,
00251 opac.opacity_abs[ip] ,
00252 hmi.Hmolec[ipMHm] ,
00253 hmi.HMinus_photo_rate);
00254 }
00255 }
00256 }
00257 }
00258 }
00259 }
00260
00261 {
00262 enum {DEBUG_LOC=false};
00263 if( DEBUG_LOC )
00264 {
00265 nelem = 0;
00266 fprintf(ioQQQ,"hotsdebugg %li \t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\n",
00267 nzone,
00268 rfield.otscon[iso.ipIsoLevNIonCon[ipH_LIKE][nelem][0]-1],
00269 rfield.otscon[iso.ipIsoLevNIonCon[ipH_LIKE][nelem][1]-1],
00270 rfield.otscon[iso.ipIsoLevNIonCon[ipH_LIKE][nelem][3]-1],
00271 rfield.otscon[iso.ipIsoLevNIonCon[ipH_LIKE][nelem][4]-1],
00272 rfield.otscon[iso.ipIsoLevNIonCon[ipH_LIKE][nelem][5]-1],
00273 rfield.otscon[iso.ipIsoLevNIonCon[ipH_LIKE][nelem][6]-1],
00274 opac.opacity_abs[iso.ipIsoLevNIonCon[ipH_LIKE][nelem][6]-1]);
00275 }
00276 }
00277
00278
00279 Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->Pdest = (realnum)save_lya_dest;
00280 Transitions[ipH_LIKE][ipHELIUM][ipH2p][ipH1s].Emis->Pdest = (realnum)save_he2lya_dest;
00281 iso.RadRecomb[ipH_LIKE][ipHELIUM][ipH1s][ipRecRad] = save_he2rec_dest;
00282
00283 nelem = ipHELIUM;
00284 if( dense.lgElmtOn[nelem] && bwnfac > 0. )
00285 {
00286
00287 Transitions[ipH_LIKE][ipHELIUM][ipH2p][ipH1s].Emis->Pdest /= (realnum)bwnfac;
00288 }
00289
00290 if( trace.lgTrace )
00291 {
00292 fprintf(ioQQQ," RT_OTS Pdest %.2e ots rate %.2e in otslin %.2e con opac %.2e\n",
00293 Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->Pdest , Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->ots ,
00294 rfield.otslin[Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].ipCont-1] ,
00295 opac.opacity_abs[Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].ipCont-1]
00296 );
00297 }
00298
00299 nOTS_Line_type = 2;
00300
00301 for( nelem=NISO; nelem < LIMELM; nelem++ )
00302 {
00303 long int ion;
00304
00305
00306 for( ion=0; ion < nelem+1-NISO; ion++ )
00307 {
00308 if( dense.xIonDense[nelem][ion+1] > 0. )
00309 {
00310 nOTS1 = nelem;
00311 nOTS2 = ion;
00312
00313 ipla = Heavy.ipLyHeavy[nelem][ion];
00314 ASSERT( ipla>0 );
00315 esc = opac.ExpmTau[ipla-1];
00316
00317 difflya = Heavy.xLyaHeavy[nelem][ion]*dense.xIonDense[nelem][ion+1];
00318
00319
00320 ots = difflya*MAX2(0.f,1.f-esc);
00321
00322
00323 ASSERT( ots >= 0.);
00324
00325
00326
00327
00328 if( ots > SMALLFLOAT )
00329 RT_OTS_AddLine(ots,ipla);
00330
00331
00332 ipla = Heavy.ipBalHeavy[nelem][ion];
00333 esc = opac.ExpmTau[ipla-1];
00334
00335 difflya = Heavy.xLyaHeavy[nelem][ion]*dense.xIonDense[nelem][ion+1];
00336
00337
00338 ots = difflya*MAX2(0.f,1.f-esc);
00339 ASSERT( ots >= 0.);
00340
00341
00342 if( ots > SMALLFLOAT )
00343 RT_OTS_AddLine(ots,ipla);
00344 }
00345 }
00346 }
00347
00348 nOTS_Line_type = 3;
00349
00350 FeII_OTS();
00351
00352 nOTS_Line_type = 4;
00353
00354 for( nOTS1=1; nOTS1 < nLevel1; nOTS1++ )
00355 {
00356 TauLines[nOTS1].Emis->ots = TauLines[nOTS1].Hi->Pop * TauLines[nOTS1].Emis->Aul * TauLines[nOTS1].Emis->Pdest;
00357 if( TauLines[nOTS1].Emis->ots > SMALLFLOAT )
00358 RT_OTS_AddLine( TauLines[nOTS1].Emis->ots , TauLines[nOTS1].ipCont);
00359 }
00360
00361 nOTS_Line_type = 5;
00362
00363 for( nOTS1=0; nOTS1 < nWindLine; nOTS1++ )
00364 {
00365 if( TauLine2[nOTS1].Hi->IonStg < TauLine2[nOTS1].Hi->nelem+1-NISO )
00366 {
00367 TauLine2[nOTS1].Emis->ots = TauLine2[nOTS1].Hi->Pop * TauLine2[nOTS1].Emis->Aul * TauLine2[nOTS1].Emis->Pdest;
00368 if( TauLine2[nOTS1].Emis->ots > SMALLFLOAT )
00369 RT_OTS_AddLine( TauLine2[nOTS1].Emis->ots , TauLine2[nOTS1].ipCont);
00370 }
00371 }
00372
00373 nOTS_Line_type = 6;
00374 for( long ipSpecies=0; ipSpecies<nSpecies; ipSpecies++ )
00375 {
00376 if( Species[ipSpecies].lgActive )
00377 {
00378 for( long ipHi=1; ipHi<Species[ipSpecies].numLevels_local; ipHi++ )
00379 {
00380 transition *t = &dBaseTrans[ipSpecies][ipHi][0];
00381 for(long ipLo = 0; ipLo < ipHi; ipLo++ )
00382 {
00383 if( t->ipCont > 0 )
00384 {
00385 t->Emis->ots = t->Hi->Pop * t->Emis->Aul * t->Emis->Pdest;
00386 RT_OTS_AddLine( t->Emis->ots , t->ipCont);
00387 }
00388 ++t;
00389 }
00390 }
00391 }
00392 }
00393
00394 nOTS_Line_type = 7;
00395
00396 H2_RT_OTS();
00397
00398 return;
00399 }
00400
00401
00402
00403 void RT_OTS_AddLine(double ots,
00404
00405 long int ip )
00406 {
00407
00408 DEBUG_ENTRY( "RT_OTS_AddLine()" );
00409
00410
00411
00412
00413
00414 if( ip==0 || ip > rfield.nflux )
00415 {
00416 return;
00417 }
00418
00419
00420 ASSERT( ots >= 0. );
00421
00422 ASSERT( ip > 0 );
00423
00424
00425
00426
00427 if( opac.opacity_abs[ip-1] > 0. )
00428 {
00429 rfield.otslin[ip-1] += (realnum)(ots/opac.opacity_abs[ip-1]);
00430 }
00431
00432 {
00433 enum {DEBUG_LOC=false};
00434 if( DEBUG_LOC && ip== 2363 )
00435 {
00436 fprintf(ioQQQ,"DEBUG ots, opc, otsr %.3e\t%.3e\t%.3e\t",
00437 ots ,
00438 opac.opacity_abs[ip-1],
00439 ots/opac.opacity_abs[ip-1] );
00440 fprintf(ioQQQ,"iteration %li type %i %i %i \n",
00441 iteration,
00442 nOTS_Line_type,
00443 nOTS1,nOTS2 );
00444 }
00445 }
00446 return;
00447 }
00448
00449
00450
00451
00452 STATIC void RT_OTS_AddCont(
00453
00454 realnum ots,
00455
00456 long int ip)
00457 {
00458
00459 DEBUG_ENTRY( "RT_OTS_AddCont()" );
00460
00461
00462
00463
00464
00465
00466
00467 if( ip > rfield.nflux )
00468 {
00469 return;
00470 }
00471
00472 ASSERT( ip > 0 );
00473 ASSERT( ots >= 0. );
00474 ASSERT( ip <= rfield.nupper );
00475
00476
00477
00478
00479 if( opac.opacity_abs[ip-1] > 0. )
00480 {
00481 rfield.otscon[ip-1] += (realnum)(ots/opac.opacity_abs[ip-1]);
00482 }
00483 return;
00484 }
00485
00486
00487
00488
00489 void RT_OTS_Update(double *SumOTS)
00490 {
00491 long int i;
00492
00493 DEBUG_ENTRY( "RT_OTS_Update()" );
00494
00495 *SumOTS = 0.;
00496
00497
00498 if( rfield.lgKillOTSLine )
00499 {
00500 for( i=0; i < rfield.nflux; i++ )
00501 {
00502 rfield.otslin[i] = 0.;
00503 }
00504 }
00505
00506
00507 *SumOTS = 0.;
00508
00509 for( i=0; i < rfield.nflux; ++i )
00510 {
00511 double CurrentInverseOpacity = 1./MAX2( SMALLDOUBLE , opac.opacity_abs[i] );
00512
00513
00514
00515 rfield.ConOTS_local_OTS_rate[i] = (realnum)((double)rfield.ConOTS_local_photons[i]*CurrentInverseOpacity);
00516
00517
00518 *SumOTS += (rfield.otscon[i] + rfield.otslin[i])*opac.opacity_abs[i];
00519
00520 rfield.SummedDif[i] = rfield.otscon[i] + rfield.otslin[i] + rfield.outlin_noplot[i]+
00521 rfield.ConInterOut[i]*rfield.lgOutOnly + rfield.outlin[0][i] +
00522 rfield.ConOTS_local_OTS_rate[i];
00523
00524 rfield.SummedCon[i] = rfield.flux[0][i] + rfield.SummedDif[i];
00525 rfield.SummedOcc[i] = rfield.SummedCon[i]*rfield.convoc[i];
00526 }
00527
00528
00529
00530 rfield.flux_accum[rfield.nflux-1] = 0.;
00531 for( i=1; i < rfield.nflux; i++ )
00532 {
00533 rfield.flux_accum[rfield.nflux-i-1] = rfield.flux_accum[rfield.nflux-i] +
00534 rfield.SummedCon[rfield.nflux-i-1];
00535 }
00536
00537
00538
00539 ASSERT( rfield.ipPlasma > 0 );
00540
00541
00542 for( i=0; i < rfield.ipPlasma-1; i++ )
00543 {
00544 rfield.otscon[i] = 0.;
00545 rfield.ConOTS_local_OTS_rate[i] = 0.;
00546 rfield.ConOTS_local_photons[i] = 0.;
00547 rfield.otslin[i] = 0.;
00548 rfield.SummedDif[i] = 0.;
00549 rfield.OccNumbBremsCont[i] = 0.;
00550 rfield.SummedCon[i] = 0.;
00551 rfield.SummedOcc[i] = 0.;
00552 rfield.ConInterOut[i] = 0.;
00553 }
00554
00555
00556 if( rfield.ipEnergyBremsThin > 0 )
00557 {
00558 for( i=rfield.ipPlasma-1; i < rfield.nflux; i++ )
00559 {
00560
00561
00562
00563
00564 realnum factor = MIN2(1.f,rfield.anu2[MAX2(0,rfield.ipEnergyBremsThin-1)] / rfield.anu2[i]);
00565
00566
00567 rfield.OccNumbBremsCont[i] = (realnum)(1./(1./SDIV(rfield.ContBoltz[i]) - 1.)) * factor;
00568 }
00569 }
00570 return;
00571 }
00572
00573
00574
00575
00576
00577 void RT_OTS_Zero( void )
00578 {
00579 long int i;
00580
00581 DEBUG_ENTRY( "RT_OTS_Zero()" );
00582
00583
00584
00585 for( i=0; i <= rfield.nflux; i++ )
00586 {
00587 rfield.SummedDif[i] = 0.;
00588
00589 rfield.otscon[i] = 0.;
00590 rfield.otslin[i] = 0.;
00591 rfield.trans_coef_zone[i] = 1.f;
00592 rfield.trans_coef_total[i] = 1.f;
00593
00594 rfield.ConInterOut[i] = 0.;
00595 rfield.outlin[0][i] = 0.;
00596 rfield.outlin_noplot[i] = 0.;
00597 rfield.SummedDif[i] = 0.;
00598
00599 rfield.SummedCon[i] = rfield.flux[0][i];
00600 rfield.SummedOcc[i] = rfield.SummedCon[i]*rfield.convoc[i];
00601 rfield.ConOTS_local_photons[i] = 0.;
00602 rfield.ConOTS_local_OTS_rate[i] = 0.;
00603 }
00604 return;
00605 }
00606
00607
00608
00609
00610 void RT_OTS_ChkSum(long int ipPnt)
00611 {
00612 static long int nInsane=0;
00613 long int i;
00614 double phisig;
00615 const int LIM_INSAME_PRT = 30;
00616
00617 DEBUG_ENTRY( "RT_OTS_ChkSum()" );
00618
00619
00620
00621
00622 for( i=rfield.ipEnergyBremsThin; i < rfield.nflux; i++ )
00623 {
00624 phisig = rfield.otscon[i] + rfield.otslin[i] + rfield.ConInterOut[i]*rfield.lgOutOnly +
00625 rfield.outlin[0][i]+
00626 rfield.outlin_noplot[i]+
00627 rfield.ConOTS_local_OTS_rate[i];
00628
00629
00630 if( phisig > 0. && rfield.SummedDif[i]> 0.)
00631 {
00632 if( fabs(rfield.SummedDif[i]/phisig-1.) > 1e-3 )
00633 {
00634 ++nInsane;
00635
00636
00637 if( nInsane < LIM_INSAME_PRT )
00638 {
00639 fprintf( ioQQQ, " PROBLEM RT_OTS_ChkSum insane SummedDif at energy %.5e error= %.2e i=%4ld\n",
00640 rfield.anu[i], rfield.SummedDif[i]/phisig - 1., i );
00641 fprintf( ioQQQ, " SummedDif, sum are%11.4e%11.4e\n",
00642 rfield.SummedDif[i], phisig );
00643 fprintf( ioQQQ, " otscon otslin ConInterOut outlin are%11.4e%11.4e%11.4e%11.4e\n",
00644 rfield.otscon[i], rfield.otslin[i]+rfield.outlin_noplot[i], rfield.ConInterOut[i],
00645 rfield.outlin[0][i]+rfield.outlin_noplot[i] );
00646 fprintf( ioQQQ, " line continuum here are %4.4s %4.4s\n",
00647 rfield.chLineLabel[i], rfield.chContLabel[i] );
00648 }
00649 }
00650 }
00651
00652 phisig += rfield.flux[0][i];
00653
00654
00655 if( phisig > 0. && rfield.SummedDif[i]> 0. )
00656 {
00657 if( fabs(rfield.SummedCon[i]/phisig-1.) > 1e-3 )
00658 {
00659 ++nInsane;
00660
00661
00662 if( nInsane < LIM_INSAME_PRT )
00663 {
00664 fprintf( ioQQQ, " PROBLEM RT_OTS_ChkSum %3ld, insane SummedCon at energy %.5e error=%.2e i=%ld\n",
00665 ipPnt, rfield.anu[i], rfield.SummedCon[i]/phisig - 1., i );
00666 fprintf( ioQQQ, " SummedCon, sum are %.4e %.4e\n",
00667 rfield.SummedCon[i], phisig );
00668 fprintf( ioQQQ, " otscon otslin ConInterOut outlin flux are%.4e %.4e %.4e %.4e %.4e\n",
00669 rfield.otscon[i], rfield.otslin[i]+rfield.outlin_noplot[i], rfield.ConInterOut[i],
00670 rfield.outlin[0][i]+rfield.outlin_noplot[i], rfield.flux[0][i] );
00671 fprintf( ioQQQ, " line continuum here are %s %s\n",
00672 rfield.chLineLabel[i], rfield.chContLabel[i]
00673 );
00674 }
00675 }
00676 }
00677 }
00678
00679 if( nInsane > 0 )
00680 {
00681 fprintf( ioQQQ, " PROBLEM RT_OTS_ChkSum too much insanity to continue.\n");
00682
00683 TotalInsanity();
00684 }
00685 return;
00686 }
00687
00688
00689
00690
00691 void RT_OTS_PrtRate(
00692
00693 double weak ,
00694
00695 int chFlag )
00696 {
00697 long int i;
00698
00699 DEBUG_ENTRY( "RT_OTS_PrtRate()" );
00700
00701
00702 ASSERT( chFlag=='l' || chFlag=='c' || chFlag=='b' );
00703
00704
00705
00706
00707
00708
00709
00710
00711
00712
00713
00714
00715
00716 if( chFlag == 'c' || chFlag == 'b' )
00717 {
00718 fprintf( ioQQQ, " DEBUG OTSCON array, anu, otscon, opac, OTS*opac limit:%.2e zone:%.2f IonConv?%c\n",
00719 weak,fnzone ,TorF(conv.lgConvIoniz) );
00720
00721 for( i=0; i < rfield.nupper; i++ )
00722 {
00723 if( rfield.otscon[i]*opac.opacity_abs[i] > weak )
00724 {
00725 fprintf( ioQQQ, " %4ld%12.4e%12.4e%12.4e%12.4e %s \n",
00726 i,
00727 rfield.anu[i],
00728 rfield.otscon[i],
00729 opac.opacity_abs[i],
00730 rfield.otscon[i]*opac.opacity_abs[i],
00731 rfield.chContLabel[i]);
00732
00733 }
00734 }
00735 }
00736
00737
00738
00739
00740 if( chFlag == 'l' || chFlag == 'b' )
00741 {
00742 fprintf( ioQQQ, "DEBUG density He %.2e He+2 %.2e O+2 %.2e\n",
00743 dense.gas_phase[ipHELIUM] , dense.xIonDense[ipHELIUM][2],
00744 dense.xIonDense[ipOXYGEN][2] );
00745 fprintf( ioQQQ, " DEBUG OTSLIN array, anu, otslin, opac, OTS*opac Lab nLine limit:%.2e zone:%.2f IonConv?%c\n",
00746 weak,fnzone,TorF(conv.lgConvIoniz) );
00747
00748 for( i=0; i < rfield.nupper; i++ )
00749 {
00750 if( rfield.otslin[i]*opac.opacity_abs[i] > weak )
00751 {
00752 fprintf( ioQQQ, " %4ld%12.4e%12.4e%12.4e%12.4e %s %3li\n",
00753 i,
00754 rfield.anu[i],
00755 rfield.otslin[i],
00756 opac.opacity_abs[i],
00757 rfield.otslin[i]*opac.opacity_abs[i],
00758 rfield.chLineLabel[i] ,
00759 rfield.line_count[i] );
00760 }
00761 }
00762 }
00763 return;
00764 }