00001
00002
00003
00004
00005
00006
00007 #include "cddefines.h"
00008 #include "physconst.h"
00009 #include "taulines.h"
00010 #include "thermal.h"
00011 #include "yield.h"
00012 #include "ipoint.h"
00013 #include "ionbal.h"
00014 #include "cddrive.h"
00015 #include "iterations.h"
00016 #include "trace.h"
00017 #include "dense.h"
00018 #include "prt.h"
00019 #include "rt.h"
00020 #include "coolheavy.h"
00021 #include "rfield.h"
00022 #include "phycon.h"
00023 #include "elementnames.h"
00024 #include "iso.h"
00025 #include "hyperfine.h"
00026 #include "hydrogenic.h"
00027 #include "lines_service.h"
00028 #include "atmdat.h"
00029 #include "lines.h"
00030 #include "radius.h"
00031
00032 STATIC void Drive_cdLine( void );
00033
00034 void lines(void)
00035 {
00036 char chLabel[5];
00037 long int i,
00038 ipnt,
00039 nelem;
00040 double BigstExtra,
00041 ExtraCool,
00042 f2, sum;
00043
00044 DEBUG_ENTRY( "lines()" );
00045
00046
00047
00048
00049
00050
00051
00052
00053
00054
00055
00056
00057
00058
00059
00060
00061
00062
00063
00064
00065 if( trace.lgTrace )
00066 {
00067 fprintf( ioQQQ, " lines called\n" );
00068 }
00069
00070
00071 if( trace.lgDrv_cdLine && LineSave.ipass > 0 )
00072 Drive_cdLine();
00073
00074
00075
00076 thermal.power += thermal.htot*radius.dVeffAper;
00077
00078
00079 thermal.FreeFreeTotHeat += CoolHeavy.brems_heat_total*radius.dVeffAper;
00080
00081
00082 rfield.comtot += rfield.cmheat*radius.dVeffAper;
00083 thermal.totcol += thermal.ctot*radius.dVeffAper;
00084
00085
00086 for( nelem=0; nelem<LIMELM; ++nelem )
00087 {
00088 hydro.cintot += iso.RecomInducCool_Rate[ipH_LIKE][nelem]*radius.dVeffAper;
00089 }
00090
00091
00092 LineSave.nsum = 0;
00093 LineSave.nComment = 0;
00094
00095
00096
00097
00098 rt.fracin = 0.5;
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113
00114 PntForLine(0.,"FILL",&i);
00115
00116
00117
00118
00119 t_ADfA::Inst().rec_lines(phycon.te,LineSave.RecCoefCNO);
00120
00121
00122 PutExtra(0.);
00123
00124
00125 linadd(0.f,0,"zero",'i' , "null placeholder");
00126
00127
00128
00129
00130 linadd( 1.e-10 , 1 , "Unit" , 'i' , "unit integration placeholder");
00131 static long int ipOneAng=-1;
00132 if( LineSave.ipass<0 )
00133 ipOneAng = ipoint( RYDLAM );
00134 lindst( 1.e-10 , 1. , "UntD" , ipOneAng , 'i' , false,"unit integration placeholder");
00135
00136
00137 lines_general();
00138
00139
00140 lines_continuum();
00141
00142
00143 lines_grains();
00144
00145
00146 iso_satellite_update();
00147
00148
00149 lines_hydro();
00150
00151
00152 lines_helium();
00153
00154 #if 0
00155
00156
00157
00158 if( iteration > 1 )
00159 {
00160 fprintf( ioQQQ,"ipHi\tipLo\tnu\tlu\tsu\tnl\tll\tsl\tWL\tintens\n" );
00161 for( long ipHi=5; ipHi<= iso.numLevels_local[ipHE_LIKE][ipHELIUM] - iso.nCollapsed_local[ipHE_LIKE][ipHELIUM]; ipHi++ )
00162 {
00163 for( long ipLo=0; ipLo<ipHi; ipLo++ )
00164 {
00165 if( Transitions[ipHE_LIKE][ipHELIUM][ipHi][ipLo].ipCont > 0 )
00166 {
00167 double relint, absint;
00168
00169 if( cdLine("He 1",
00170 Transitions[ipHE_LIKE][ipHELIUM][ipHi][ipLo].WLAng,
00171 &relint, &absint ) )
00172 {
00173
00174
00175 if( Transitions[ipHE_LIKE][ipHELIUM][ipHi][ipLo].WLAng < 1.1E4 &&
00176 Transitions[ipHE_LIKE][ipHELIUM][ipHi][ipLo].WLAng > 3.59E3 &&
00177 ipLo!=3 && ipLo!=4 && relint >= 0.0009 )
00178 {
00179 fprintf( ioQQQ,"%li\t%li\t%li\t%li\t%li\t%li\t%li\t%li\t%e\t%e\n",
00180 ipHi,
00181 ipLo,
00182 StatesElemNEW[ipHELIUM][ipHELIUM-ipHE_LIKE][ipHi].n,
00183 StatesElemNEW[ipHELIUM][ipHELIUM-ipHE_LIKE][ipHi].l,
00184 StatesElemNEW[ipHELIUM][ipHELIUM-ipHE_LIKE][ipHi].S,
00185 StatesElemNEW[ipHELIUM][ipHELIUM-ipHE_LIKE][ipLo].n,
00186 StatesElemNEW[ipHELIUM][ipHELIUM-ipHE_LIKE][ipLo].l,
00187 StatesElemNEW[ipHELIUM][ipHELIUM-ipHE_LIKE][ipLo].S,
00188 Transitions[ipHE_LIKE][ipHELIUM][ipHi][ipLo].WLAng,
00189 relint );
00190 }
00191 }
00192 }
00193 }
00194 }
00195 }
00196 #endif
00197
00198
00199 lines_lv1_li_ne();
00200
00201
00202 lines_lv1_na_ar();
00203
00204
00205 lines_lv1_k_zn();
00206
00207
00208 sum = PrtLineSum();
00209
00210
00211 if( LineSave.ipass > 0 )
00212 {
00213 LineSv[LineSave.nsum].SumLine[0] = 0.;
00214 LineSv[LineSave.nsum].SumLine[1] = 0.;
00215 }
00216
00217 linadd(sum/radius.dVeffAper,0,"Stoy",'i' ,
00218 "Stoy method energy sum ");
00219
00220
00221 i = StuffComment( "recombination" );
00222 linadd( 0., (realnum)i , "####", 'i' ,
00223 "recombination lines");
00224
00225
00226
00227
00228
00229 for( i=0; i < 471; i++ )
00230 {
00231
00232
00233 if( LineSave.ipass <= 0 )
00234 {
00235
00236 strcpy( chLabel, elementnames.chElementSym[(long)(LineSave.RecCoefCNO[0][i])-1] );
00237 strcat( chLabel, elementnames.chIonStage[(long)(LineSave.RecCoefCNO[0][i]-
00238 LineSave.RecCoefCNO[1][i]+1.01)-1] );
00239 }
00240 else
00241 chLabel[0] = ' ';
00242
00243
00244
00245
00246
00247
00248
00249
00250
00251
00252 if( dense.eden < 1e8 )
00253 {
00254 nelem = (long)LineSave.RecCoefCNO[0][i]-1;
00255 long int ion = (long)(LineSave.RecCoefCNO[0][i]-LineSave.RecCoefCNO[1][i]+2)-1;
00256 f2 = LineSave.RecCoefCNO[3][i]*dense.eden*
00257 dense.xIonDense[nelem][ion];
00258
00259
00260 f2 = f2*1.99e-8/LineSave.RecCoefCNO[2][i];
00261 }
00262 else
00263 {
00264 f2 = 0.;
00265 }
00266
00267 PntForLine(LineSave.RecCoefCNO[2][i],chLabel,&ipnt);
00268 lindst(f2,LineSave.RecCoefCNO[2][i],chLabel,ipnt, 'r',true ,
00269 "recombination line");
00270 }
00271
00272
00273 i = StuffComment( "level2 lines" );
00274 linadd( 0., (realnum)i , "####", 'i' ,
00275 "level2 lines");
00276
00277
00278
00279 ExtraCool = 0.;
00280 BigstExtra = 0.;
00281 for( i=0; i < nWindLine; i++ )
00282 {
00283 if( TauLine2[i].Hi->IonStg < TauLine2[i].Hi->nelem+1-NISO )
00284 {
00285 PutLine(&TauLine2[i],"level 2 line");
00286 if( TauLine2[i].Coll.cool > BigstExtra )
00287 {
00288 BigstExtra = TauLine2[i].Coll.cool;
00289 thermal.ipMaxExtra = i+1;
00290 }
00291 ExtraCool += TauLine2[i].Coll.cool;
00292 }
00293 }
00294
00295 thermal.GBarMax = MAX2(thermal.GBarMax,(realnum)(ExtraCool/thermal.ctot));
00296
00297
00298 i = StuffComment( "hyperfine structure" );
00299 linadd( 0., (realnum)i , "####", 'i' ,
00300 "hyperfine structure lines ");
00301
00302
00303 linadd( hyperfine.cooling_total, 0., "hfin", 'i' ,
00304 "total cooling all hyperfine structure lines ");
00305
00306
00307 hyperfine.cooling_max = (realnum)MAX2(hyperfine.cooling_max,hyperfine.cooling_total/thermal.ctot);
00308
00309
00310 for( i=0; i < nHFLines; i++ )
00311 {
00312 PutLine(&HFLines[i],
00313 "hyperfine structure line");
00314 }
00315
00316
00317 i = StuffComment( "database lines" );
00318 linadd( 0., (realnum)i , "####", 'i' ,
00319 "database lines ");
00320
00321
00322 for( i=0; i < linesAdded2; i++ )
00323 {
00324
00325 PutLine(dBaseLines[i].tran, "lines from third-party databases", dBaseLines[i].tran->Hi->chLabel);
00326 }
00327
00328
00329 i = StuffComment( "inner shell" );
00330 linadd( 0., (realnum)i , "####", 'i' ,
00331 "inner shell lines");
00332
00333
00334 for( i=0; i < t_yield::Inst().nlines(); ++i )
00335 {
00336 double xInten =
00337
00338 dense.xIonDense[t_yield::Inst().nelem(i)][t_yield::Inst().ion(i)] *
00339
00340 ionbal.PhotoRate_Shell[t_yield::Inst().nelem(i)][t_yield::Inst().ion(i)][t_yield::Inst().nshell(i)][0]*
00341
00342 t_yield::Inst().yield(i) *
00343
00344 t_yield::Inst().energy(i) * EN1RYD;
00345
00346
00347 if( LineSave.ipass == 0 )
00348 {
00349
00350 strcpy( chLabel , elementnames.chElementSym[t_yield::Inst().nelem(i)] );
00351 strcat( chLabel , elementnames.chIonStage[t_yield::Inst().ion_emit(i)] );
00352 # if 0
00353
00354 if( t_yield::Inst().ion(i) == 0 && t_yield::Inst().nelem(i) == ipIRON )
00355 fprintf(ioQQQ,"DEBUGyeild\t%s\t%.3f\t%.3e\n",
00356
00357 chLabel , t_yield::Inst().energy(i)*EVRYD, t_yield::Inst().yield(i) );
00358 # endif
00359 }
00360
00361
00362 lindst(
00363
00364 xInten,
00365
00366 (realnum)RYDLAM / t_yield::Inst().energy(i),
00367
00368 chLabel ,
00369
00370 t_yield::Inst().ipoint(i),
00371
00372 'r',
00373
00374 true ,
00375 "inner shell line");
00376 }
00377
00378
00379
00380 {
00381 static long nLineSave=-1 , ndLineSave=-1;
00382 if( LineSave.ipass == 0 )
00383 {
00384 nLineSave = LineSave.nsum;
00385 ndLineSave = LineSave.nsum;
00386 }
00387 else if( LineSave.ipass > 0 )
00388 {
00389
00390 if( nLineSave<= 0 || ndLineSave < 0 )
00391 TotalInsanity();
00392
00393
00394
00395
00396 if( nLineSave != LineSave.nsum )
00397 {
00398 fprintf( ioQQQ, "DISASTER number of lines in LineSave.nsum changed between pass 0 and 1 - this is impossible\n" );
00399 fprintf( ioQQQ, "DISASTER LineSave.nsum is %li and nLineSave is %li\n",
00400 LineSave.nsum ,
00401 nLineSave);
00402 ShowMe();
00403 cdEXIT(EXIT_FAILURE);
00404 }
00405 if( ndLineSave != LineSave.nsum )
00406 {
00407 fprintf( ioQQQ, "DISASTER number of lines in LineSave.nsum changed between pass 0 and 1 - this is impossible\n" );
00408 fprintf( ioQQQ, "DISASTER LineSave.nsum is %li and ndLineSave is %li\n",
00409 LineSave.nsum ,
00410 ndLineSave);
00411 ShowMe();
00412 cdEXIT(EXIT_FAILURE);
00413 }
00414 }
00415 }
00416
00417
00418 lines_molecules();
00419
00420 if( trace.lgTrace )
00421 {
00422 fprintf( ioQQQ, " lines returns\n" );
00423 }
00424 return;
00425 }
00426
00427
00428 STATIC void Drive_cdLine( void )
00429 {
00430 long int j;
00431 bool lgMustPrintHeader = true;
00432 double absval , rel;
00433
00434 DEBUG_ENTRY( "Drive_cdLine()" );
00435
00436 for( j=1; j < LineSave.nsum; j++ )
00437 {
00438 if( cdLine( LineSv[j].chALab , LineSv[j].wavelength , &absval , &rel ) <= 0 )
00439 {
00440
00441 if( lgMustPrintHeader )
00442 {
00443 fprintf(ioQQQ,"n\tlab\twl\n");
00444 lgMustPrintHeader = false;
00445 }
00446
00447 fprintf(ioQQQ,"%li\t%s\t%f\n", j, LineSv[j].chALab , LineSv[j].wavelength );
00448 }
00449 }
00450 fprintf( ioQQQ, " Thanks for checking on the cdLine routine!\n" );
00451 cdEXIT(EXIT_FAILURE);
00452 }
00453