00001
00002
00003
00004
00005
00006 #include "cddefines.h"
00007 #include "physconst.h"
00008 #include "cddrive.h"
00009 #include "lines_service.h"
00010 #include "iso.h"
00011 #include "continuum.h"
00012 #include "stopcalc.h"
00013 #include "hyperfine.h"
00014 #include "dense.h"
00015 #include "grainvar.h"
00016 #include "version.h"
00017 #include "rt.h"
00018 #include "he.h"
00019 #include "ionbal.h"
00020 #include "taulines.h"
00021 #include "hydrogenic.h"
00022 #include "lines.h"
00023 #include "trace.h"
00024 #include "hcmap.h"
00025 #include "hmi.h"
00026 #include "save.h"
00027 #include "h2.h"
00028 #include "conv.h"
00029 #include "dynamics.h"
00030 #include "opacity.h"
00031 #include "geometry.h"
00032 #include "elementnames.h"
00033 #include "ca.h"
00034 #include "broke.h"
00035 #include "pressure.h"
00036 #include "mole.h"
00037 #include "atoms.h"
00038 #include "abund.h"
00039 #include "rfield.h"
00040 #include "colden.h"
00041 #include "phycon.h"
00042 #include "timesc.h"
00043 #include "hextra.h"
00044 #include "radius.h"
00045 #include "iterations.h"
00046 #include "fudgec.h"
00047 #include "called.h"
00048 #include "magnetic.h"
00049 #include "wind.h"
00050 #include "secondaries.h"
00051 #include "struc.h"
00052 #include "oxy.h"
00053 #include "input.h"
00054 #include "thermal.h"
00055 #include "atmdat.h"
00056 #include "warnings.h"
00057
00058
00059 STATIC void chkCaHeps(double *totwid);
00060
00061
00062 STATIC void prt_smooth_predictions(void);
00063
00064 void PrtComment(void)
00065 {
00066 char chLbl[11],
00067 chLine[INPUT_LINE_LENGTH];
00068
00069 bool lgAbort_flag,
00070 lgThick,
00071 lgLots_of_moles;
00072
00073 long int dum1,
00074 dum2,
00075 i,
00076 imas,
00077 ipLo ,
00078 ipHi ,
00079 ipISO,
00080 nelem,
00081 isav,
00082 j,
00083 nc,
00084 nline,
00085 nn,
00086 nneg,
00087 ns,
00088 nw;
00089
00090 double big_ion_jump,
00091 absint ,
00092 aj,
00093 alpha,
00094 big,
00095 bigm,
00096 comfrc,
00097 differ,
00098 error,
00099 flur,
00100 freqn,
00101 rate,
00102 ratio,
00103 rel,
00104 small,
00105 tauneg,
00106 ts ,
00107 HBeta,
00108 relfl ,
00109 relhm,
00110 fedest,
00111 GetHeat,
00112 SumNeg,
00113 c4363,
00114 t4363,
00115 totwid;
00116
00117 double VolComputed , VolExpected , ConComputed , ConExpected;
00118
00119 bool lgLotsSolids;
00120
00121 DEBUG_ENTRY( "PrtComment()" );
00122
00123 if( 0 && lgAbort )
00124 {
00125 return;
00126 }
00127
00128 if( trace.lgTrace )
00129 {
00130 fprintf( ioQQQ, " PrtComment called.\n" );
00131 }
00132
00133
00134
00135
00136
00137
00138 iterations.lgIterAgain = false;
00139
00140
00141 wcnint();
00142
00143 if( t_version::Inst().nBetaVer > 0 )
00144 {
00145 sprintf( chLine,
00146 " !This is beta test version %ld and is intended for testing only.",
00147 t_version::Inst().nBetaVer );
00148 bangin(chLine);
00149 }
00150
00151
00152
00153
00154 if( broke.lgFixit && !t_version::Inst().lgRelease )
00155 {
00156 sprintf( chLine, " !The code needs to be fixed - search for fixit()." );
00157 bangin(chLine);
00158 }
00159
00160
00161
00162
00163 if( broke.lgCheckit && !t_version::Inst().lgRelease )
00164 {
00165 sprintf( chLine, " !New code needs to be reviewed - search for CodeReview()." );
00166 bangin(chLine);
00167 }
00168
00169
00170 if( conv.lgBadStop )
00171 {
00172
00173 sprintf( warnings.chRgcln[0], " W-Calculation stopped because %s Iteration%3ld of %ld",
00174 StopCalc.chReasonStop, iteration, iterations.itermx + 1 );
00175 sprintf( chLine, " W-Calculation stopped because %s",
00176 StopCalc.chReasonStop );
00177 warnin(chLine);
00178 sprintf( chLine, " W-This was not intended." );
00179 warnin(chLine);
00180 }
00181 else
00182 {
00183
00184 if( (conv.lgAutoIt && iteration != iterations.itermx + 1) &&
00185 iteration > 2 )
00186 {
00187 sprintf( warnings.chRgcln[0],
00188 " Calculation stopped because %s Iteration %ld of %ld, not converged due to %s",
00189 StopCalc.chReasonStop,
00190 iteration,
00191 iterations.itermx + 1,
00192 conv.chNotConverged );
00193 }
00194 else
00195 {
00196 sprintf( warnings.chRgcln[0],
00197 " Calculation stopped because %s Iteration %ld of %ld",
00198 StopCalc.chReasonStop, iteration, iterations.itermx + 1 );
00199 }
00200 }
00201
00202
00203
00204 if( (!geometry.lgZoneSet) && geometry.lgZoneTrp )
00205 {
00206 conv.lgBadStop = true;
00207 sprintf( chLine,
00208 " W-Calculation stopped because default number of zones reached. Was this intended???" );
00209 warnin(chLine);
00210 sprintf( chLine,
00211 " W-Default limit can be increased while retaining this check with the SET NEND command." );
00212 warnin(chLine);
00213 }
00214
00215
00216
00217
00218 if( radius.lgDrMinUsed || radius.lgdR2Small )
00219 {
00220 conv.lgBadStop = true;
00221 sprintf( chLine,
00222 " W-Calculation stopped zone thickness became too thin. This was not intended." );
00223 warnin(chLine);
00224 sprintf( chLine,
00225 " W-The most likely reason was an uncontrolled oscillation." );
00226 warnin(chLine);
00227 ShowMe();
00228 }
00229
00230 if( radius.lgdR2Small )
00231 {
00232 sprintf( chLine,
00233 " W-This happened because the globule scale became very small relative to the depth." );
00234 warnin(chLine);
00235 sprintf( chLine,
00236 " W-This problem is described in Hazy." );
00237 warnin(chLine);
00238 }
00239
00240
00241
00242 ASSERT( nzone < struc.nzlim );
00243
00244 if( struc.testr[nzone-1] == 0. )
00245 struc.testr[nzone-1] = struc.testr[nzone-2];
00246
00247 if( struc.ednstr[nzone-1] == 0. )
00248 struc.ednstr[nzone-1] = struc.ednstr[nzone-2];
00249
00250
00251 rel = radius.depth/radius.rinner;
00252 if( rel < 0.1 )
00253 {
00254 sprintf( warnings.chRgcln[1], " The geometry is plane-parallel." );
00255 }
00256 else if( rel >= 0.1 && rel < 3. )
00257 {
00258 sprintf( warnings.chRgcln[1], " The geometry is a thick shell." );
00259 }
00260 else
00261 {
00262 sprintf( warnings.chRgcln[1], " The geometry is spherical." );
00263 }
00264
00265
00266
00267
00268
00269
00270
00271
00272 if( broke.lgBroke )
00273 {
00274 sprintf( chLine, " W-The code is broken - search for broken()." );
00275 warnin(chLine);
00276 }
00277
00278
00279 if( dense.lgEdenBad )
00280 {
00281 if( dense.nzEdenBad == nzone )
00282 {
00283 sprintf( chLine, " C-The assumed electron density was incorrect for the last zone." );
00284 caunin(chLine);
00285 sprintf( chLine, " C-Did a temperature discontinuity occur??" );
00286 caunin(chLine);
00287 }
00288 else
00289 {
00290 sprintf( chLine, " W-The assumed electron density was incorrect during the calculation. This is bad." );
00291 warnin(chLine);
00292 ShowMe();
00293 }
00294 }
00295
00296 if( lgAbort )
00297 {
00298 sprintf( chLine, " W-The calculation aborted. Something REALLY went wrong!" );
00299 warnin(chLine);
00300 }
00301
00302
00303 if( hcmap.lgMapDone && !hcmap.lgMapOK )
00304 {
00305 sprintf( chLine, " !The thermal map had changes in slope - check map output." );
00306 bangin(chLine);
00307 }
00308
00309
00310
00311 if( fudgec.nfudge > 0 || fudgec.lgFudgeUsed )
00312 {
00313 sprintf( chLine, " !Fudge factors were used or were checked. Why?" );
00314 bangin(chLine);
00315 }
00316
00317 if( dense.gas_phase[ipHYDROGEN] > 1.1e13 )
00318 {
00319 if( dense.gas_phase[ipHYDROGEN] > 1e15 )
00320 {
00321 sprintf( chLine, " C-Density greater than 10**15, heavy elements are very uncertain." );
00322 caunin(chLine);
00323 }
00324 else
00325 {
00326 sprintf( chLine, " C-Density greater than 10**13" );
00327 caunin(chLine);
00328 }
00329 }
00330
00331
00332 if( cdLine("Pump",4861.36f,&relfl,&absint)<=0 )
00333 {
00334 fprintf( ioQQQ, " PROBLEM Did not find Pump H-beta, set to unity\n" );
00335 relfl = 1.;
00336 absint = 1.;
00337 }
00338
00339
00340
00341
00342 if( cdLine( "H 1",4861.36f,&HBeta,&absint)<=0 )
00343 {
00344 fprintf( ioQQQ, " NOTE Did not find H 1 H-beta - set intensity to unity, "
00345 "will not check on importance of H 1 pumping.\n" );
00346 HBeta = 1.;
00347 absint = 1.;
00348 }
00349 else
00350 {
00351
00352 if( HBeta>SMALLFLOAT )
00353 {
00354 flur = relfl/HBeta;
00355 if( flur > 0.1 )
00356 {
00357 sprintf( chLine, " !Continuum fluorescent production of H-beta was very important." );
00358 bangin(chLine);
00359 }
00360 else if(flur > 0.01 )
00361 {
00362 sprintf( chLine, " Continuum fluorescent production of H-beta was significant." );
00363 notein(chLine);
00364 }
00365 }
00366 }
00367
00368
00369 if( iteration > 1 && conv.lgAutoIt && iteration<iterations.itermx)
00370 {
00371 sprintf( chLine , " Iteration not converged because %s.",
00372 conv.chNotConverged );
00373 notein(chLine);
00374 }
00375
00376
00377 if( wind.lgBallistic() && ((!wind.lgWindOK) || wind.windv < 1e6) )
00378 {
00379 sprintf( chLine, " C-Wind velocity below sonic point; solution is not valid." );
00380 caunin(chLine);
00381 }
00382
00383
00384 if( !wind.lgStatic() )
00385 {
00386 rel = wind.emdot/(wind.windv*dense.gas_phase[ipHYDROGEN])/radius.r1r0sq;
00387 if( fabs(1.-rel)> 0.02 )
00388 {
00389 sprintf( chLine, " C-Wind mass flux error is %g%%",fabs(1.-rel)*100. );
00390 caunin(chLine);
00391 fprintf(ioQQQ,"DEBUG emdot\t%.3e\t%.3e\t%.3e\t%.3e\n",
00392 wind.emdot , wind.windv*dense.gas_phase[ipHYDROGEN],wind.windv,dense.gas_phase[ipHYDROGEN]);
00393 }
00394 }
00395
00396
00397 if( nzone >= struc.nzlim )
00398 {
00399 TotalInsanity();
00400 }
00401
00402
00403 if( !lgConserveEnergy() )
00404 {
00405 ShowMe();
00406 lgAbort = true;
00407 fprintf(ioQQQ,"\n\n Enable per zone energy conservation check by setting "
00408 "CHECK_ENERGY_EVERY_ZONE=true in cloudy.cpp, recompile, then rerun.\n");
00409 }
00410
00411
00412
00413 if( hextra.cryden*magnetic.lgB > 0. )
00414 {
00415 sprintf( chLine,
00416 " !Magnetic field & cosmic rays both present. Their interactions are not treated." );
00417 bangin(chLine);
00418 }
00419
00420
00421 if( hextra.cryden== 0. && StopCalc.TempLoStopZone < phycon.TEMP_STOP_DEFAULT)
00422 {
00423 sprintf( chLine,
00424 " !Background cosmic rays are not included - is this physical? It affects the chemistry." );
00425 bangin(chLine);
00426 }
00427
00428
00429 if( hextra.cryden > 0. && (colden.colden[ipCOL_H0]/10. + colden.colden[ipCOL_Hp]) > 1e23 )
00430 {
00431 sprintf( chLine,
00432 " C-Model is thick to cosmic rays, which are on." );
00433 caunin(chLine);
00434 }
00435
00436
00437 if( hextra.cryden == 0. && iso.gamnc[ipH_LIKE][ipHYDROGEN][ipH1s] < 1e-17 )
00438 {
00439 sprintf( chLine,
00440 " !Ionization rate fell below background cosmic ray ionization rate. Should this be added too?" );
00441 bangin(chLine);
00442 sprintf( chLine,
00443 " ! Use the COSMIC RAY BACKGROUND command." );
00444 bangin(chLine);
00445 }
00446
00447
00448 if( lgTestCodeCalled )
00449 {
00450 sprintf( chLine, " !Test code is in place." );
00451 bangin(chLine);
00452 }
00453
00454
00455 if( rfield.lgComUndr )
00456 {
00457 sprintf( chLine,
00458 " !Compton cooling rate underflows to zero. Is this important?" );
00459 bangin(chLine);
00460 }
00461
00462
00463 if( input.lgUnderscoreFound )
00464 {
00465 sprintf( chLine,
00466 " !Some input lines contained underscores, these were changed to spaces." );
00467 bangin(chLine);
00468 }
00469
00470
00471 if( input.lgBracketFound )
00472 {
00473 sprintf( chLine,
00474 " !Some input lines contained [ or ], these were changed to spaces." );
00475 bangin(chLine);
00476 }
00477
00478
00479 if( rfield.lgHionRad )
00480 {
00481 sprintf( chLine,
00482 " !There is no hydrogen-ionizing radiation. Was this intended?" );
00483 bangin(chLine);
00484 }
00485
00486
00487 if( thermal.nUnstable > 0 )
00488 {
00489 sprintf( chLine,
00490 " Derivative of net cooling negative and so possibly thermally unstable in%4ld zones.",
00491 thermal.nUnstable );
00492 notein(chLine);
00493 }
00494
00495
00496 if( nzone > 1 &&
00497 (realnum)(thermal.nUnstable)/(realnum)(nzone) > 0.25 )
00498 {
00499 sprintf( chLine,
00500 " !A large fraction of the zones were possibly thermally unstable,%4ld out of%4ld",
00501 thermal.nUnstable, nzone );
00502 bangin(chLine);
00503 }
00504
00505
00506 if( thermal.CoolHeatMax > 0.2 )
00507 {
00508 sprintf( chLine,
00509 " !Negative cooling reached %6.1f%% of the local heating, due to %4.4s %.1f",
00510 thermal.CoolHeatMax*100., thermal.chCoolHeatMax, thermal.wlCoolHeatMax );
00511 bangin(chLine);
00512 }
00513 else if( thermal.CoolHeatMax > 0.05 )
00514 {
00515 sprintf( chLine,
00516 " Negative cooling reached %6.1f%% of the local heating, due to %4.4s %.2f",
00517 thermal.CoolHeatMax*100., thermal.chCoolHeatMax, thermal.wlCoolHeatMax );
00518 notein(chLine);
00519 }
00520
00521
00522 if( dynamics.HeatMax > 0.05 )
00523 {
00524 sprintf( chLine,
00525 " !Advection heating reached %.2f%% of the local heating.",
00526 dynamics.HeatMax*100. );
00527 bangin(chLine);
00528 }
00529 else if( dynamics.HeatMax > 0.005 )
00530 {
00531 sprintf( chLine,
00532 " Advection heating reached %.2f%% of the local heating.",
00533 dynamics.HeatMax*100. );
00534 notein(chLine);
00535 }
00536
00537
00538 if( dynamics.CoolMax > 0.05 )
00539 {
00540 sprintf( chLine,
00541 " !Advection cooling reached %.2f%% of the local cooling.",
00542 dynamics.CoolMax*100. );
00543 bangin(chLine);
00544 }
00545 else if( dynamics.CoolMax > 0.005 )
00546 {
00547 sprintf( chLine,
00548 " Advection cooling reached %.2f%% of the local heating.",
00549 dynamics.CoolMax*100. );
00550 notein(chLine);
00551 }
00552
00553
00554
00555 if( dynamics.lgTimeDependentStatic && dynamics.lgRecom )
00556 {
00557 if( rfield.uh > 1. )
00558 {
00559 sprintf( chLine,
00560 " W-Time dependent ionization front cannot now handle strong-R cases - the ionization parameter is too large." );
00561 warnin(chLine);
00562 }
00563 else if( rfield.uh > 0.1 )
00564 {
00565 sprintf( chLine,
00566 " C-Time dependent ionization front cannot now handle strong-R cases - the ionization parameter is too large." );
00567 caunin(chLine);
00568 }
00569 }
00570
00571
00572 if( hydro.HCollIonMax > 0.10 )
00573 {
00574 sprintf( chLine,
00575 " !Thermal collisional ionization of H reached %.2f%% of the local ionization rate.",
00576 hydro.HCollIonMax*100. );
00577 bangin(chLine);
00578 }
00579 else if( hydro.HCollIonMax > 0.005 )
00580 {
00581 sprintf( chLine,
00582 " Thermal collisional ionization of H reached %.2f%% of the local ionization rate.",
00583 hydro.HCollIonMax*100. );
00584 notein(chLine);
00585 }
00586
00587
00588 if( !atmdat.lgHCaseBOK[1][ipHYDROGEN] )
00589 {
00590 sprintf( chLine,
00591 " Te-ne bounds of Case B lookup table exceeded, H I Case B line intensities set to zero." );
00592 notein(chLine);
00593 }
00594 if( !atmdat.lgHCaseBOK[1][ipHELIUM] )
00595 {
00596 sprintf( chLine,
00597 " Te-ne bounds of Case B lookup table exceeded, He II Case B line intensities set to zero." );
00598 notein(chLine);
00599 }
00600
00601 if( dense.EdenMax>1e8 )
00602 {
00603 sprintf( chLine,
00604 " !The high electron density makes the Nussbaumer/Storey CNO recombination predictions unreliable." );
00605 bangin(chLine);
00606 }
00607
00608
00609 if( secondaries.SecHIonMax > 0.10 )
00610 {
00611 sprintf( chLine,
00612 " !Suprathermal collisional ionization of H reached %.2f%% of the local H ionization rate.",
00613 secondaries.SecHIonMax*100. );
00614 bangin(chLine);
00615 }
00616 else if( secondaries.SecHIonMax > 0.005 )
00617 {
00618 sprintf( chLine,
00619 " Suprathermal collisional ionization of H reached %.2f%% of the local H ionization rate.",
00620 secondaries.SecHIonMax*100. );
00621 notein(chLine);
00622 }
00623
00624
00625 if( hmi.HeatH2DexcMax > 0.05 )
00626 {
00627 sprintf( chLine,
00628 " !H2 vib deexec heating reached %.2f%% of the local heating.",
00629 hmi.HeatH2DexcMax*100. );
00630 bangin(chLine);
00631 }
00632 else if( hmi.HeatH2DexcMax > 0.005 )
00633 {
00634 sprintf( chLine,
00635 " H2 vib deexec heating reached %.2f%% of the local heating.",
00636 hmi.HeatH2DexcMax*100. );
00637 notein(chLine);
00638 }
00639
00640
00641 if( hmi.CoolH2DexcMax > 0.05 )
00642 {
00643 sprintf( chLine,
00644 " !H2 deexec cooling reached %.2f%% of the local heating.",
00645 hmi.CoolH2DexcMax*100. );
00646 bangin(chLine);
00647 }
00648 else if( hmi.CoolH2DexcMax > 0.005 )
00649 {
00650 sprintf( chLine,
00651 " H2 deexec cooling reached %.2f%% of the local heating.",
00652 hmi.CoolH2DexcMax*100. );
00653 notein(chLine);
00654 }
00655
00656
00657 if( atmdat.HIonFracMax > 0.10 )
00658 {
00659 sprintf( chLine,
00660 " !Charge transfer ionization of H reached %.1f%% of the local H ionization rate.",
00661 atmdat.HIonFracMax*100. );
00662 bangin(chLine);
00663 }
00664 else if( atmdat.HIonFracMax > 0.005 )
00665 {
00666 sprintf( chLine,
00667 " Charge transfer ionization of H reached %.2f%% of the local H ionization rate.",
00668 atmdat.HIonFracMax*100. );
00669 notein(chLine);
00670 }
00671
00672
00673 if( atmdat.HCharHeatMax > 0.05 )
00674 {
00675 sprintf( chLine,
00676 " !Charge transfer heating reached %.2f%% of the local heating.",
00677 atmdat.HCharHeatMax*100. );
00678 bangin(chLine);
00679 }
00680 else if( atmdat.HCharHeatMax > 0.005 )
00681 {
00682 sprintf( chLine,
00683 " Charge transfer heating reached %.2f%% of the local heating.",
00684 atmdat.HCharHeatMax*100. );
00685 notein(chLine);
00686 }
00687
00688 if( atmdat.HCharCoolMax > 0.05 )
00689 {
00690 sprintf( chLine,
00691 " !Charge transfer cooling reached %.2f%% of the local heating.",
00692 atmdat.HCharCoolMax*100. );
00693 bangin(chLine);
00694 }
00695 else if( atmdat.HCharCoolMax > 0.005 )
00696 {
00697 sprintf( chLine,
00698 " Charge transfer cooling reached %.2f%% of the local heating.",
00699 atmdat.HCharCoolMax*100. );
00700 notein(chLine);
00701 }
00702
00703
00704 if( atoms.xMg2Max > 0.1 )
00705 {
00706 sprintf( chLine,
00707 " !Photoionization of upper level of Mg II 2798 reached %.1f%% of the total Mg+ photo rate.",
00708 atoms.xMg2Max*100. );
00709 bangin(chLine);
00710 }
00711 else if( atoms.xMg2Max > 0.01 )
00712 {
00713 sprintf( chLine,
00714 " Photoionization of upper level of Mg II 2798 reached %.1f%% of the total Mg+ photo rate.",
00715 atoms.xMg2Max*100. );
00716 notein(chLine);
00717 }
00718
00719
00720 if( oxy.poimax > 0.1 )
00721 {
00722 sprintf( chLine,
00723 " !Photoionization of upper levels of [O I] reached %.1f%% of the total O destruction rate.",
00724 oxy.poimax*100. );
00725 bangin(chLine);
00726 }
00727 else if( oxy.poimax > 0.01 )
00728 {
00729 sprintf( chLine,
00730 " Photoionization of upper levels of [O I] reached %.1f%% of the total O destruction rate.",
00731 oxy.poimax*100. );
00732 notein(chLine);
00733 }
00734
00735
00736 if( (oxy.poiii2Max + oxy.poiii3Max) > 0.1 )
00737 {
00738 sprintf( chLine,
00739 " !Photoionization of upper levels of [O III] reached %.1f%% of the total O++ photo rate.",
00740 (oxy.poiii2Max + oxy.poiii3Max)*100. );
00741 bangin(chLine);
00742 }
00743 else if( (oxy.poiii2Max + oxy.poiii3Max) > 0.01 )
00744 {
00745 sprintf( chLine,
00746 " Photoionization of upper levels of [O III] reached %.1f%% of the total O++ photo rate.",
00747 (oxy.poiii2Max + oxy.poiii3Max)*100. );
00748 notein(chLine);
00749 }
00750
00751
00752 if( he.frac_he0dest_23S > 0.1 )
00753 {
00754 sprintf( chLine,
00755 " !Destruction of He 2TriS reached %.1f%% of the total He0 dest rate"
00756 " at zone %li, %.1f%% of that was photoionization.",
00757 he.frac_he0dest_23S*100,
00758 he.nzone,
00759 he.frac_he0dest_23S_photo*100. );
00760 bangin(chLine);
00761 }
00762 else if( he.frac_he0dest_23S > 0.01 )
00763 {
00764 sprintf( chLine,
00765 " Destruction of He 2TriS reached %.1f%% of the total He0 dest rate"
00766 " at zone %li, %.1f%% of that was photoionization.",
00767 he.frac_he0dest_23S*100,
00768 he.nzone,
00769 he.frac_he0dest_23S_photo*100. );
00770 notein(chLine);
00771 }
00772
00773
00774 for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
00775 {
00776 if( !iso.lgCritDensLMix[ipISO] && dense.lgElmtOn[ipISO] )
00777 {
00778 sprintf( chLine,
00779 " The density is too low to l-mix the lowest %s I collapsed level. "
00780 " More resolved levels are needed for accurate line ratios.",
00781 elementnames.chElementSym[ipISO]);
00782 notein(chLine);
00783 }
00784 }
00785
00786 for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
00787 {
00788
00789 for( nelem=ipISO; nelem<LIMELM; ++nelem )
00790 {
00791 if( iso.lgLevelsLowered[ipISO][nelem] && dense.lgElmtOn[nelem] )
00792 {
00793 sprintf( chLine, " !Continuum was lowered into model %s-like %s due to high density. Highest n is %li",
00794 elementnames.chElementSym[ipISO],
00795 elementnames.chElementSym[nelem],
00796 iso.n_HighestResolved_local[ipISO][nelem]+iso.nCollapsed_local[ipISO][nelem]);
00797 bangin(chLine);
00798 }
00799 else if( iso.lgLevelsEverLowered[ipISO][nelem] && dense.lgElmtOn[nelem] )
00800 {
00801 sprintf( chLine, " !Continuum was lowered into model %s-like %s due to high density at SOME point but NOT at the last zone.",
00802 elementnames.chElementSym[ipISO],
00803 elementnames.chElementNameShort[nelem]);
00804 bangin(chLine);
00805 }
00806 }
00807 }
00808
00809
00810 if( !rfield.lgMMok )
00811 {
00812 sprintf( chLine,
00813 " C-Continuum not defined in extreme infrared - Compton scat, grain heating, not treated properly?" );
00814 caunin(chLine);
00815 }
00816
00817 if( !rfield.lgHPhtOK )
00818 {
00819 sprintf( chLine,
00820 " C-Continuum not defined at photon energies which ionize excited states of H, important for H- and ff heating." );
00821 caunin(chLine);
00822 }
00823
00824 if( !rfield.lgXRayOK )
00825 {
00826 sprintf( chLine,
00827 " C-Continuum not defined at X-Ray energies - Compton scattering and Auger ionization wrong?" );
00828 caunin(chLine);
00829 }
00830
00831 if( !rfield.lgGamrOK )
00832 {
00833 sprintf( chLine,
00834 " C-Continuum not defined at gamma-ray energies - pair production and Compton scattering OK?" );
00835 caunin(chLine);
00836 }
00837
00838 if( continuum.lgCon0 )
00839 {
00840 sprintf( chLine, " C-Continuum zero at some energies." );
00841 caunin(chLine);
00842 }
00843
00844 if( continuum.lgCoStarInterpolationCaution )
00845 {
00846 sprintf( chLine , " C-CoStarInterpolate interpolated between non-adjoining tracks, this may not be valid." );
00847 caunin(chLine);
00848 }
00849
00850 if( rfield.lgOcc1Hi )
00851 {
00852 sprintf( chLine,
00853 " !The continuum occupation number at 1 Ryd is greater than unity." );
00854 bangin(chLine);
00855 }
00856
00857
00858 if( radius.lgDR2Big )
00859 {
00860 sprintf( chLine,
00861 " C-The thickness of the first zone was set larger than optimal by a SET DR command." );
00862 caunin(chLine);
00863
00864
00865 if( nzone<=1 )
00866 sprintf( chLine,
00867 " C-Consider using the STOP THICKNESS command instead." );
00868 caunin(chLine);
00869 }
00870
00871
00872 if( cdLine("TOTL",4363,&t4363,&absint)<=0 )
00873 {
00874 fprintf( ioQQQ, " PrtComment could not find total O III 4363 with cdLine.\n" );
00875 ShowMe();
00876 cdEXIT(EXIT_FAILURE);
00877 }
00878
00879 if( cdLine("Coll",4363,&c4363,&absint)<=0 )
00880 {
00881 fprintf( ioQQQ, " PrtComment could not find collisional O III 4363 with cdLine.\n" );
00882 ShowMe();
00883 cdEXIT(EXIT_FAILURE);
00884 }
00885
00886
00887 if( HBeta > 1e-35 )
00888 {
00889
00890 if( t4363/HBeta > 1e-5 && dense.gas_phase[ipHYDROGEN] < 1e8 )
00891 {
00892 ratio = (t4363 - c4363)/t4363;
00893 if( ratio > 0.01 )
00894 {
00895 sprintf( chLine,
00896 " !Non-collisional excitation of [O III] 4363 reached %.2f%% of the total.",
00897 ratio*100. );
00898 bangin(chLine);
00899 }
00900 else if( ratio > 0.001 )
00901 {
00902 sprintf( chLine,
00903 " Non-collisional excitation of [O III] 4363 reached %.2f%% of the total.",
00904 ratio*100. );
00905 notein(chLine);
00906 }
00907 }
00908 }
00909
00910
00911 if( rfield.lgPlasNu )
00912 {
00913 sprintf( chLine,
00914 " !The largest plasma frequency was %.2e Ryd = %.2e micron The continuum is set to 0 below this.",
00915 rfield.plsfrqmax,
00916
00917 RYDLAM/rfield.plsfrqmax/1e4);
00918 bangin(chLine);
00919 }
00920
00921 if( rfield.occmax > 0.1 )
00922 {
00923 if( rfield.occmnu > 1e-4 )
00924 {
00925 sprintf( chLine,
00926 " !The largest continuum occupation number was %.3e at %.3e Ryd.",
00927 rfield.occmax, rfield.occmnu );
00928 bangin(chLine);
00929 }
00930 else
00931 {
00932
00933
00934 sprintf( chLine,
00935 " The largest continuum occupation number was %.3e at %.3e Ryd.",
00936 rfield.occmax, rfield.occmnu );
00937 notein(chLine);
00938 }
00939 }
00940
00941 if( rfield.occmax > 1e4 && rfield.occ1nu > 0. )
00942 {
00943
00944 if( rfield.occ1nu < 0.0912 )
00945 {
00946 sprintf( chLine,
00947 " The continuum occupation number fell below 1 at %.3e microns.",
00948 0.0912/rfield.occ1nu );
00949 notein(chLine);
00950 }
00951 else if( rfield.occ1nu < 1. )
00952 {
00953 sprintf( chLine,
00954 " The continuum occupation number fell below 1 at %.3e Angstroms.",
00955 912./rfield.occ1nu );
00956 notein(chLine);
00957 }
00958 else
00959 {
00960 sprintf( chLine,
00961 " The continuum occupation number fell below 1 at %.3e Ryd.",
00962 rfield.occ1nu );
00963 notein(chLine);
00964 }
00965 }
00966
00967 if( rfield.tbrmax > 1e3 )
00968 {
00969 sprintf( chLine,
00970 " !The largest continuum brightness temperature was %.3eK at %.3e Ryd.",
00971 rfield.tbrmax, rfield.tbrmnu );
00972 bangin(chLine);
00973 }
00974
00975 if( rfield.tbrmax > 1e4 )
00976 {
00977
00978 if( rfield.tbr4nu < 0.0912 )
00979 {
00980 sprintf( chLine,
00981 " The continuum brightness temperature fell below 10,000K at %.3e microns.",
00982 0.0912/rfield.tbr4nu );
00983 notein(chLine);
00984 }
00985 else if( rfield.tbr4nu < 1. )
00986 {
00987 sprintf( chLine,
00988 " The continuum brightness temperature fell below 10,000K at %.3e Angstroms.",
00989 912./rfield.tbr4nu );
00990 notein(chLine);
00991 }
00992 else
00993 {
00994 sprintf( chLine,
00995 " The continuum brightness temperature fell below 10,000K at %.3e Ryd.",
00996 rfield.tbr4nu );
00997 notein(chLine);
00998 }
00999 }
01000
01001
01002 if( DoppVel.TurbVel > 0. && strcmp(dense.chDenseLaw,"CPRE") == 0 )
01003 {
01004 sprintf( chLine,
01005 " !Both constant pressure and turbulence makes no physical sense?" );
01006 bangin(chLine);
01007 }
01008
01009
01010 if( geometry.FillFac < 1. && strcmp(dense.chDenseLaw,"CPRE") == 0 )
01011 {
01012 sprintf( chLine,
01013 " !Both constant pressure and a filling factor makes no physical sense?" );
01014 bangin(chLine);
01015 }
01016
01017
01018 if( gv.lgDustOn() && abund.lgAbnSolar )
01019 {
01020 sprintf( chLine,
01021 " !Grains are present, but the gas phase abundances were left at the solar default. This is not physical." );
01022 bangin(chLine);
01023 }
01024
01025
01026 if( abund.lgDepln && !gv.lgDustOn() )
01027 {
01028 sprintf( chLine,
01029 " !Grains are not present, but the gas phase abundances were depleted. This is not physical." );
01030 bangin(chLine);
01031 }
01032
01033 if( gv.lgDustOn() )
01034 {
01035 long nBin=0L, nFail=0L;
01036 for( size_t nd=0; nd < gv.bin.size(); nd++ )
01037 {
01038 if( gv.bin[nd]->QHeatFailures > 0L )
01039 {
01040 ++nBin;
01041 nFail += gv.bin[nd]->QHeatFailures;
01042 }
01043 }
01044 if( nFail > 0 )
01045 {
01046 sprintf( chLine,
01047 " !The grain quantum heating treatment failed to converge %ld time(s) in %ld bin(s).", nFail, nBin );
01048 bangin(chLine);
01049 }
01050 }
01051
01052 #if 0
01053
01054
01056 if( gv.lgDustOn() )
01057 {
01058 bool lgPAHsPresent_and_constant = false;
01059 for( size_t nd=0; nd < gv.bin.size(); nd++ )
01060 {
01061 lgPAHsPresent_and_constant = lgPAHsPresent_and_constant ||
01062
01063 (gv.bin[nd]->lgPAHsInIonizedRegion );
01064 }
01065 if( lgPAHsPresent_and_constant )
01066 {
01067 sprintf( chLine,
01068 " C-PAH's were present in the ionized region, this has never been observed in H II regions." );
01069 caunin(chLine);
01070 }
01071 }
01072 #endif
01073
01074
01075 if( thermal.lgTemperatureConstant && thermal.ConstTemp*1.0001 < phycon.TEnerDen )
01076 {
01077 sprintf( chLine,
01078 " C-The continuum energy density temperature (%g K)"
01079 " is greater than the gas kinetic temperature (%g K).",
01080 phycon.TEnerDen , thermal.ConstTemp);
01081 caunin(chLine);
01082 sprintf( chLine, " C-This is unphysical." );
01083 caunin(chLine);
01084 }
01085
01086
01087 if( !gv.lgDustOn() && phycon.TEnerDen < 800. )
01088 {
01089 sprintf( chLine,
01090 " Grains were not present but might survive in this environment (energy density temperature was %.2eK)",
01091 phycon.TEnerDen );
01092 notein(chLine);
01093 }
01094
01095
01096 AgeCheck();
01097
01098
01099
01100 chkCaHeps(&totwid);
01101 if( totwid > 121. )
01102 {
01103 sprintf( chLine, " H-eps and Ca H overlap." );
01104 notein(chLine);
01105 }
01106
01107
01108 if( !phycon.lgPhysOK )
01109 {
01110 sprintf( chLine, " !A physical process has been disabled." );
01111 bangin(chLine);
01112 }
01113
01114
01115 if( dense.gas_phase[ipHYDROGEN] < 1e8 )
01116 {
01117 if( oxy.r5007Max > 0.0263f )
01118 {
01119 sprintf( chLine,
01120 " !Photoionization of upper level of [O III] 5007 reached %.2e%% of the radiative lifetime.",
01121 oxy.r5007Max*100. );
01122 bangin(chLine);
01123 }
01124 else if( oxy.r5007Max > 0.0263f/10.f )
01125 {
01126 sprintf( chLine,
01127 " Photoionization of upper level of [O III] 5007 reached %.2e%% of the radiative lifetime.",
01128 oxy.r5007Max*100. );
01129 notein(chLine);
01130 }
01131 if( oxy.r4363Max > 1.78f )
01132 {
01133 sprintf( chLine,
01134 " !Photoionization of upper level of [O III] 4363 reached %.2e%% of the radiative lifetime.",
01135 oxy.r4363Max*100. );
01136 bangin(chLine);
01137 }
01138 else if( oxy.r4363Max > 1.78f/10.f )
01139 {
01140 sprintf( chLine,
01141 " Photoionization of upper level of [O III] 4363 reached %.2e%% of the radiative lifetime.",
01142 oxy.r4363Max*100. );
01143 notein(chLine);
01144 }
01145 }
01146
01147
01148
01149 error = fabs(thermal.power-thermal.totcol)/SDIV((thermal.power + thermal.totcol)/2.);
01150 if( thermal.lgTemperatureConstant )
01151 {
01152 if( error > 0.05 )
01153 {
01154 sprintf( chLine,
01155 " !Heating - cooling mismatch =%5.1f%%. Caused by constant temperature assumption. ",
01156 error*100. );
01157 bangin(chLine);
01158 }
01159 }
01160
01161 else
01162 {
01163 if( error > 0.05 && error < 0.2 )
01164 {
01165 sprintf( chLine, " C-Heating - cooling mismatch =%.1f%%. What\'s wrong?",
01166 error*100. );
01167 caunin(chLine);
01168 }
01169 else if( error >= 0.2 )
01170 {
01171 sprintf( chLine, " W-Heating - cooling mismatch =%.2e%%. What\'s wrong????",
01172 error*100. );
01173 warnin(chLine);
01174 }
01175 }
01176
01177
01178 if( ca.Ca2RmLya > 0.01 )
01179 {
01180 sprintf( chLine,
01181 " Photoionization of Ca+ 2D level by Ly-alpha reached %6.1f%% of the total rate out.",
01182 ca.Ca2RmLya*100. );
01183 notein(chLine);
01184 }
01185
01186
01187 if( hydro.nLyaHot > 0 )
01188 {
01189 if( hydro.TLyaMax/hydro.TeLyaMax > 1.05 )
01190 {
01191 sprintf( chLine,
01192 " !The excitation temp of Lya exceeded the electron temp, largest value was %.2eK (gas temp there was %.2eK, zone%4ld)",
01193 hydro.TLyaMax, hydro.TeLyaMax, hydro.nZTLaMax );
01194 bangin(chLine);
01195 }
01196 }
01197
01198
01199
01200
01201
01202 if( cdLine("Line",0,&SumNeg,&absint)<=0 )
01203 {
01204 fprintf( ioQQQ, " did not get sumneg cdLine\n" );
01205 ShowMe();
01206 cdEXIT(EXIT_FAILURE);
01207 }
01208
01209
01210 if( cdLine("TotH",0,&GetHeat,&absint)<=0 )
01211 {
01212 fprintf( ioQQQ, " did not get GetHeat cdLine\n" );
01213 ShowMe();
01214 cdEXIT(EXIT_FAILURE);
01215 }
01216
01217 if( GetHeat > 0. )
01218 {
01219 SumNeg /= GetHeat;
01220 if( SumNeg > 0.1 )
01221 {
01222 sprintf( chLine,
01223 " !Line absorption heating reached %.2f%% of the global heating.",
01224 SumNeg*100. );
01225 bangin(chLine);
01226 }
01227 else if( SumNeg > 0.01 )
01228 {
01229 sprintf( chLine,
01230 " Line absorption heating reached %.2f%% of the global heating.",
01231 SumNeg*100. );
01232 notein(chLine);
01233 }
01234 }
01235
01236
01237 if( input.lgSetNoBuffering )
01238 {
01239 sprintf( chLine,
01240 " !NO BUFFERING command was entered - this increases exec time by LARGE amounts.");
01241 bangin(chLine);
01242 }
01243
01244
01245 if( thermal.GBarMax > 0.1 )
01246 {
01247 ASSERT( thermal.ipMaxExtra > 0 );
01248 strcpy( chLbl, chLineLbl(&TauLine2[thermal.ipMaxExtra-1]) );
01249
01250 sprintf( chLine,
01251 " !G-bar cooling lines reached %.2f%% of the local cooling. Line=%.10s",
01252 thermal.GBarMax*100., chLbl );
01253 bangin(chLine);
01254 }
01255
01256 else if( thermal.GBarMax > 0.01 )
01257 {
01258 strcpy( chLbl, chLineLbl(&TauLine2[thermal.ipMaxExtra-1]) );
01259
01260 sprintf( chLine,
01261 " G-bar cooling lines reached %.2f%% of the local cooling. Line=%.10s",
01262 thermal.GBarMax*100., chLbl );
01263 notein(chLine);
01264 }
01265
01266
01267 if( hyperfine.cooling_max > 0.1 )
01268 {
01269 sprintf( chLine,
01270 " !Hyperfine structure line cooling reached %.2f%% of the local cooling.",
01271 hyperfine.cooling_max*100.);
01272 bangin(chLine);
01273 }
01274
01275 else if( hyperfine.cooling_max > 0.01 )
01276 {
01277 sprintf( chLine,
01278 " Hyperfine structure line cooling reached %.2f%% of the local cooling.",
01279 hyperfine.cooling_max*100. );
01280 notein(chLine);
01281 }
01282
01283
01284
01285 if( thermal.HeatLineMax > 0.1 )
01286 {
01287 if( thermal.levlmax == 1 )
01288 {
01289
01290
01291
01292
01293 strcpy( chLbl, chLineLbl(&TauLines[thermal.ipHeatlmax] ) );
01294 }
01295 else if( thermal.levlmax == 2 )
01296 {
01297
01298 strcpy( chLbl, chLineLbl(&TauLine2[thermal.ipHeatlmax]) );
01299 }
01300 else if( thermal.levlmax == 3 )
01301 {
01302
01303 strcpy( chLbl, chLineLbl(&HFLines[thermal.ipHeatlmax]) );
01304 }
01305 else if( thermal.levlmax == 4 )
01306 {
01307
01308 strcpy( chLbl, chLineLbl(dBaseLines[thermal.ipHeatlmax].tran) );
01309 }
01310 else
01311 {
01312 fprintf( ioQQQ, " PROBLEM DISASTER PrtComment has insane levlmax,=%5ld\n",
01313 thermal.levlmax );
01314 }
01315 sprintf( chLine,
01316 " !Line absorption heating reached %.2f%% of the local heating - largest by level%2ld line %.10s",
01317 thermal.HeatLineMax*100., thermal.levlmax, chLbl );
01318 bangin(chLine);
01319 }
01320
01321 else if( thermal.HeatLineMax > 0.01 )
01322 {
01323 sprintf( chLine,
01324 " Line absorption heating reached %.2f%% of the local heating.",
01325 thermal.HeatLineMax*100. );
01326 notein(chLine);
01327 }
01328
01329 if( ionbal.CompHeating_Max > 0.05 )
01330 {
01331 sprintf( chLine,
01332 " !Bound Compton heating reached %.2f%% of the local heating.",
01333 ionbal.CompHeating_Max*100. );
01334 bangin(chLine);
01335 }
01336 else if( ionbal.CompHeating_Max > 0.01 )
01337 {
01338 sprintf( chLine,
01339 " Bound Compton heating reached %.2f%% of the local heating.",
01340 ionbal.CompHeating_Max*100. );
01341 notein(chLine);
01342 }
01343
01344
01345 for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
01346 {
01347 for( nelem=ipISO; nelem<LIMELM; ++nelem )
01348 {
01349 if( dense.lgElmtOn[nelem] )
01350 {
01351
01352 long int nmax = iso.numLevels_local[ipISO][nelem];
01353
01354
01355 for( ipHi=1; ipHi < nmax - 1; ++ipHi )
01356 {
01357 for( ipLo=0; ipLo < ipHi; ++ipLo )
01358 {
01359 if( Transitions[ipISO][nelem][ipHi][ipLo].Emis->Aul <= iso.SmallA )
01360 continue;
01361
01362
01363 if( Transitions[ipISO][nelem][ipHi][ipLo].Emis->TauIn < -0.1 )
01364 {
01365 sprintf( chLine,
01366 " !Some iso-structure lines mased: %s-like %s, line %li-%li had optical depth %.2e",
01367 elementnames.chElementSym[ipISO],
01368 elementnames.chElementNameShort[nelem],
01369 ipHi , ipLo ,
01370 Transitions[ipISO][nelem][ipHi][ipLo].Emis->TauIn );
01371 bangin(chLine);
01372 }
01373 }
01374 }
01375 }
01376 }
01377 }
01378
01379 if( dense.gas_phase[ipHYDROGEN] < 1e7 )
01380 {
01381
01382 lgThick = false;
01383 tauneg = 0.;
01384 alpha = 0.;
01385
01386 for( i=3; i <= nLevel1; i++ )
01387 {
01388
01389 if( TauLines[i].EnergyWN < 10000. )
01390 {
01391 if( TauLines[i].Emis->TauIn > 1. )
01392 {
01393 lgThick = true;
01394 alpha = MAX2(alpha,(double)TauLines[i].Emis->TauIn);
01395 }
01396 else if( TauLines[i].Emis->TauIn < (realnum)tauneg )
01397 {
01398 tauneg = TauLines[i].Emis->TauIn;
01399 strcpy( chLbl, chLineLbl(&TauLines[i]) );
01400 }
01401 }
01402 }
01403
01404 if( lgThick )
01405 {
01406 sprintf( chLine,
01407 " !Some infrared fine structure lines are optically thick: largest tau was %.2e",
01408 alpha );
01409 bangin(chLine);
01410 }
01411
01412 if( tauneg < -0.01 )
01413 {
01414 sprintf( chLine,
01415 " !Some fine structure lines mased: line %s had optical depth %.2e",
01416 chLbl, tauneg );
01417 bangin(chLine);
01418 }
01419 }
01420
01421
01422 tauneg = 0.;
01423 alpha = 0.;
01424 for( i=1; i <= nLevel1; i++ )
01425 {
01426
01427 if( TauLines[i].EnergyWN >= 10000. )
01428 {
01429 if( TauLines[i].Emis->TauIn < (realnum)tauneg )
01430 {
01431 tauneg = TauLines[i].Emis->TauIn;
01432 strcpy( chLbl, chLineLbl(&TauLines[i]) );
01433 }
01434 }
01435 }
01436
01437
01438 if( tauneg < -0.01 )
01439 {
01440 sprintf( chLine,
01441 " !Some level1 lines mased: most negative ion and tau were: %s %.2e",
01442 chLbl, tauneg );
01443 bangin(chLine);
01444 }
01445
01446
01447
01448
01449 if( geometry.lgStatic && iterations.lgLastIt && (iteration == 1) &&
01450 !geometry.lgStaticNoIt)
01451 {
01452 sprintf( chLine, " C-I must iterate when SPHERE STATIC is set." );
01453 caunin(chLine);
01454 iterations.lgIterAgain = true;
01455 }
01456
01457
01458 if( save.lgPunContinuum && iteration == 1 && iterations.lgLastIt)
01459 {
01460 sprintf( chLine, " C-I must iterate when save continuum output is done." );
01461 caunin(chLine);
01462 iterations.lgIterAgain = true;
01463 }
01464
01466
01467 if( iso.TwoNu_induc_dn_max[ipH_LIKE][ipHYDROGEN] > 1. )
01468 {
01469 sprintf( chLine, " !Rate of induced H 2-photon emission reached %.2e s^-1",
01470 iso.TwoNu_induc_dn_max[ipH_LIKE][ipHYDROGEN] );
01471 bangin(chLine);
01472 }
01473
01474 else if( iso.TwoNu_induc_dn_max[ipH_LIKE][ipHYDROGEN] > 0.01 )
01475 {
01476 sprintf( chLine, " Rate of induced H 2-photon emission reached %.2e s^-1",
01477 iso.TwoNu_induc_dn_max[ipH_LIKE][ipHYDROGEN] );
01478 notein(chLine);
01479 }
01480
01481
01482 if( hydro.FracInd > 0.01 )
01483 {
01484 sprintf( chLine,
01485 " Induced recombination was %5.1f%% of the total for H level%3ld",
01486 hydro.FracInd*100., hydro.ndclev );
01487 notein(chLine);
01488 }
01489
01490 if( hydro.fbul > 0.01 )
01491 {
01492 sprintf( chLine,
01493 " Stimulated emission was%6.1f%% of the total for H transition%3ld -%3ld",
01494 hydro.fbul*100., hydro.nbul + 1, hydro.nbul );
01495 notein(chLine);
01496 }
01497
01498
01499
01500 if( cdLine("Fe 2",1216,&fedest,&absint)<=0 )
01501 {
01502 fprintf( ioQQQ, " Did not find Fe II Lya\n" );
01503 ShowMe();
01504 cdEXIT(EXIT_FAILURE);
01505 }
01506
01507
01508 if( cdLine("TOTL",1215.68f,&relhm,&absint)<=0 )
01509 {
01510 fprintf( ioQQQ, " Did not find Lya\n" );
01511 ShowMe();
01512 cdEXIT(EXIT_FAILURE);
01513 }
01514
01515 if( relhm > 0. )
01516 {
01517 ratio = fedest/(fedest + relhm);
01518 if( ratio > 0.1 )
01519 {
01520 sprintf( chLine, " !Fe II destruction of Ly-a removed %.1f%% of the line.",
01521 ratio *100.);
01522 bangin(chLine);
01523 }
01524 else if( ratio > 0.01 )
01525 {
01526 sprintf( chLine, " Fe II destruction of Ly-a removed %.1f%% of the line.",
01527 ratio );
01528 notein(chLine);
01529 }
01530 }
01531
01532 if( cdLine("H-CT",6563,&relhm,&absint)<=0 )
01533 {
01534 fprintf( ioQQQ, " Comment did not find H-CT H-alpha\n" );
01535 ShowMe();
01536 cdEXIT(EXIT_FAILURE);
01537 }
01538
01539 if( HBeta > 0. )
01540 {
01541 if( relhm/HBeta > 0.01 )
01542 {
01543 sprintf( chLine,
01544 " !Mutual neutralization production of H-alpha was significant." );
01545 bangin(chLine);
01546 }
01547 }
01548
01549
01550 if( hydro.lgHiPop2 )
01551 {
01552 sprintf( chLine,
01553 " The population of H n=2 reached %.2e relative to the ground state.",
01554 hydro.pop2mx );
01555 notein(chLine);
01556 }
01557
01558
01559 for( ipISO=ipH_LIKE; ipISO<=ipHE_LIKE; ++ipISO )
01560 {
01561 for( nelem=0; nelem < LIMELM; nelem++ )
01562 {
01563 if( iso.CaseBCheck[ipISO][nelem] > 1.5 )
01564 {
01565 sprintf( chLine,
01566 " Ratio of computed diffuse emission to case B reached %g for iso %li element %li",
01567 iso.CaseBCheck[ipISO][nelem] , ipISO , nelem+1 );
01568 notein(chLine);
01569 }
01570 }
01571 }
01572
01573
01574 if( thermal.thist > 1e9 )
01575 {
01576
01577
01578
01579 if( thermal.thist > 1.0001e10 )
01580 {
01581 sprintf( chLine, " W-Electrons were relativistic; High TE=%.2e",
01582 thermal.thist );
01583 warnin(chLine);
01584 }
01585 else
01586 {
01587 sprintf( chLine, " C-Electrons were mildly relativistic; High TE=%.2e",
01588 thermal.thist );
01589 caunin(chLine);
01590 }
01591 }
01592
01593
01594 rate = timesc.TimeErode*2e-26;
01595 if( rate > 1e-35 )
01596 {
01597
01598
01599
01600 ts = (1./rate)/3e7;
01601 if( ts < 1e3 )
01602 {
01603 sprintf( chLine, " !Timescale-photoerosion of Fe=%.2e yr",
01604 ts );
01605 bangin(chLine);
01606 }
01607 else if( ts < 1e9 )
01608 {
01609 sprintf( chLine, " Timescale-photoerosion of Fe=%.2e yr",
01610 ts );
01611 notein(chLine);
01612 }
01613 }
01614
01615
01616 comfrc = rfield.comtot/SDIV(thermal.power);
01617 if( comfrc > 0.01 )
01618 {
01619 sprintf( chLine, " Compton heating was %5.1f%% of the total.",
01620 comfrc*100. );
01621 notein(chLine);
01622 }
01623
01624
01625 if( comfrc > 0.01 && rfield.cinrat > 0.05 )
01626 {
01627 sprintf( chLine,
01628 " !Induced Compton heating was %.2e of the total Compton heating.",
01629 rfield.cinrat );
01630 bangin(chLine);
01631 }
01632
01633
01634 if( timesc.tcmptn > 5e17 )
01635 {
01636 if( comfrc > 0.05 )
01637 {
01638 sprintf( chLine,
01639 " C-Compton cooling is significant and the equilibrium timescale (%.2e s) is longer than the Hubble time.",
01640 timesc.tcmptn );
01641 caunin(chLine);
01642 }
01643 else
01644 {
01645 sprintf( chLine,
01646 " Compton cooling equilibrium timescale (%.2e s) is longer than Hubble time.",
01647 timesc.tcmptn );
01648 notein(chLine);
01649 }
01650 }
01651
01652 if( timesc.time_therm_long > 5e17 )
01653 {
01654 sprintf( chLine,
01655 " C-Thermal equilibrium timescale, %.2e s, longer than Hubble time; this cloud is not time-steady.",
01656 timesc.time_therm_long );
01657 caunin(chLine);
01658 }
01659
01660
01661
01662
01663 if( log10(radius.depth) > colden.rjnmin )
01664 {
01665
01666 aj = pow(10.,colden.ajmmin - log10(SOLAR_MASS));
01667 if( strcmp(dense.chDenseLaw,"CPRE") == 0 )
01668 {
01669 sprintf( chLine,
01670 " C-Cloud thicker than smallest Jeans length=%8.2ecm; stability problems? (smallest Jeans mass=%8.2eMo)",
01671 pow((realnum)10.f,colden.rjnmin), aj );
01672 caunin(chLine);
01673 }
01674 else
01675 {
01676 sprintf( chLine,
01677 " Cloud thicker than smallest Jeans length=%8.2ecm; stability problems? (smallest Jeans mass=%8.2eMo)",
01678 pow((realnum)10.f,colden.rjnmin), aj );
01679 notein(chLine);
01680 }
01681 }
01682
01683
01684 for( size_t nd=0; nd < gv.bin.size(); nd++ )
01685 {
01686 if( gv.bin[nd]->TeGrainMax > gv.bin[nd]->Tsublimat )
01687 {
01688 sprintf( chLine,
01689 " W-Maximum temperature of grain%-12.12s was %.2eK, above its sublimation temperature, %.2eK.",
01690 gv.bin[nd]->chDstLab, gv.bin[nd]->TeGrainMax,
01691 gv.bin[nd]->Tsublimat );
01692 warnin(chLine);
01693 }
01694 else if( gv.bin[nd]->TeGrainMax > gv.bin[nd]->Tsublimat* 0.9 )
01695 {
01696 sprintf( chLine,
01697 " C-Maximum temperature of grain%-12.12s was %.2eK, near its sublimation temperature, %.2eK.",
01698 gv.bin[nd]->chDstLab, gv.bin[nd]->TeGrainMax,
01699 gv.bin[nd]->Tsublimat );
01700 caunin(chLine);
01701 }
01702 }
01703
01704 if( gv.lgNegGrnDrg )
01705 {
01706 sprintf( chLine, " !Grain drag force <0." );
01707 bangin(chLine);
01708 }
01709
01710
01711 if( gv.GrnElecDonateMax > 0.05 )
01712 {
01713 sprintf( chLine,
01714 " !Grains donated %5.1f%% of the total electrons in some regions.",
01715 gv.GrnElecDonateMax*100. );
01716 bangin(chLine);
01717 }
01718 else if( gv.GrnElecDonateMax > 0.005 )
01719 {
01720 sprintf( chLine,
01721 " Grains donated %5.1f%% of the total electrons in some regions.",
01722 gv.GrnElecDonateMax*100. );
01723 notein(chLine);
01724 }
01725
01726
01727 if( gv.GrnElecHoldMax > 0.05 )
01728 {
01729 sprintf( chLine,
01730 " !Grains contained %5.1f%% of the total electrons in some regions.",
01731 gv.GrnElecHoldMax*100. );
01732 bangin(chLine);
01733 }
01734 else if( gv.GrnElecHoldMax > 0.005 )
01735 {
01736 sprintf( chLine,
01737 " Grains contained %5.1f%% of the total electrons in some regions.",
01738 gv.GrnElecHoldMax*100. );
01739 notein(chLine);
01740 }
01741
01742
01743 if( gv.dphmax > 0.5 )
01744 {
01745 sprintf( chLine,
01746 " !Local grain-gas photoelectric heating rate reached %5.1f%% of the total.",
01747 gv.dphmax*100. );
01748 bangin(chLine);
01749 }
01750 else if( gv.dphmax > 0.05 )
01751 {
01752 sprintf( chLine,
01753 " Local grain-gas photoelectric heating rate reached %5.1f%% of the total.",
01754 gv.dphmax*100. );
01755 notein(chLine);
01756 }
01757
01758 if( gv.TotalDustHeat/SDIV(thermal.power) > 0.01 )
01759 {
01760 sprintf( chLine,
01761 " Global grain photoelectric heating of gas was%5.1f%% of the total.",
01762 gv.TotalDustHeat/thermal.power*100. );
01763 notein(chLine);
01764 if( gv.TotalDustHeat/thermal.power > 0.25 )
01765 {
01766 sprintf( chLine,
01767 " !Grain photoelectric heating is VERY important." );
01768 bangin(chLine);
01769 }
01770 }
01771
01772
01773 if( gv.dclmax > 0.05 )
01774 {
01775 sprintf( chLine,
01776 " Local grain-gas cooling of gas rate reached %5.1f%% of the total.",
01777 gv.dclmax*100. );
01778 notein(chLine);
01779 }
01780
01781
01782 if( h2.renorm_max > 1.05 )
01783 {
01784 if( h2.renorm_max > 1.2 )
01785 {
01786 sprintf( chLine,
01787 " !The large H2 molecule - main chemistry network renormalization factor reached %.2f.",
01788 h2.renorm_max);
01789 bangin(chLine);
01790 }
01791 else
01792 {
01793 sprintf( chLine,
01794 " The large H2 molecule - main chemistry network renormalization factor reached %.2f.",
01795 h2.renorm_max);
01796 notein(chLine);
01797 }
01798 }
01799 if( h2.renorm_min < 0.95 )
01800 {
01801 if( h2.renorm_min < 0.8 )
01802 {
01803 sprintf( chLine,
01804 " !The large H2 molecule - main chemistry network renormalization factor reached %.2f.",
01805 h2.renorm_min);
01806 bangin(chLine);
01807 }
01808 else
01809 {
01810 sprintf( chLine,
01811 " The large H2 molecule - main chemistry network renormalization factor reached %.2f.",
01812 h2.renorm_min);
01813 notein(chLine);
01814 }
01815 }
01816
01817
01818 if( hmi.h2pmax > 0.10 )
01819 {
01820 sprintf( chLine,
01821 " !The local H2+ photodissociation heating rate reached %5.1f%% of the total heating.",
01822 hmi.h2pmax*100. );
01823 bangin(chLine);
01824 }
01825
01826 else if( hmi.h2pmax > 0.01 )
01827 {
01828 sprintf( chLine,
01829 " The local H2+ photodissociation heating rate reached %.1f%% of the total heating.",
01830 hmi.h2pmax*100. );
01831 notein(chLine);
01832 }
01833
01834
01835 if( hmi.h2dfrc > 0.1 )
01836 {
01837 sprintf( chLine,
01838 " !The local H2 photodissociation heating rate reached %.1f%% of the total heating.",
01839 hmi.h2dfrc*100. );
01840 bangin(chLine);
01841 }
01842 else if( hmi.h2dfrc > 0.01 )
01843 {
01844 sprintf( chLine,
01845 " The local H2 photodissociation heating rate reached %.1f%% of the total heating.",
01846 hmi.h2dfrc*100. );
01847 notein(chLine);
01848 }
01849
01850
01851 if( hmi.h2line_cool_frac > 0.1 )
01852 {
01853 sprintf( chLine,
01854 " !The local H2 cooling rate reached %.1f%% of the local cooling.",
01855 hmi.h2line_cool_frac*100. );
01856 bangin(chLine);
01857 }
01858 else if( hmi.h2line_cool_frac > 0.01 )
01859 {
01860 sprintf( chLine,
01861 " The local H2 cooling rate reached %.1f%% of the local cooling.",
01862 hmi.h2line_cool_frac*100. );
01863 notein(chLine);
01864 }
01865
01866 if( hmi.h2dtot/SDIV(thermal.power) > 0.01 )
01867 {
01868 sprintf( chLine,
01869 " Global H2 photodissociation heating of gas was %.1f%% of the total heating.",
01870 hmi.h2dtot/thermal.power*100. );
01871 notein(chLine);
01872 if( hmi.h2dtot/thermal.power > 0.25 )
01873 {
01874 sprintf( chLine, " H2 photodissociation heating is VERY important." );
01875 notein(chLine);
01876 }
01877 }
01878
01879
01880 if( co.codfrc > 0.25 )
01881 {
01882 sprintf( chLine,
01883 " !Local CO photodissociation heating rate reached %.1f%% of the total.",
01884 co.codfrc*100. );
01885 bangin(chLine);
01886 }
01887 else if( co.codfrc > 0.05 )
01888 {
01889 sprintf( chLine,
01890 " Local CO photodissociation heating rate reached %.1f%% of the total.",
01891 co.codfrc*100. );
01892 notein(chLine);
01893 }
01894
01895 if( co.codtot/SDIV(thermal.power) > 0.01 )
01896 {
01897 sprintf( chLine,
01898 " Global CO photodissociation heating of gas was %.1f%% of the total.",
01899 co.codtot/thermal.power*100. );
01900 notein(chLine);
01901 if( co.codtot/thermal.power > 0.25 )
01902 {
01903 sprintf( chLine, " CO photodissociation heating is VERY important." );
01904 notein(chLine);
01905 }
01906 }
01907
01908 if( thermal.lgEdnGTcm )
01909 {
01910 sprintf( chLine,
01911 " Energy density of radiation field was greater than the Compton temperature. Is this physical?" );
01912 notein(chLine);
01913 }
01914
01915
01916 if( hydro.cintot/SDIV(thermal.power) > 0.01 )
01917 {
01918 sprintf( chLine, " Induced recombination cooling was %.1f%% of the total.",
01919 hydro.cintot/thermal.power*100. );
01920 notein(chLine);
01921 }
01922
01923
01924 if( thermal.FreeFreeTotHeat/SDIV(thermal.power) > 0.1 )
01925 {
01926 sprintf( chLine, " !Free-free heating was %.1f%% of the total.",
01927 thermal.FreeFreeTotHeat/thermal.power*100. );
01928 bangin(chLine);
01929 }
01930 else if( thermal.FreeFreeTotHeat/SDIV(thermal.power) > 0.01 )
01931 {
01932 sprintf( chLine, " Free-free heating was %.1f%% of the total.",
01933 thermal.FreeFreeTotHeat/thermal.power*100. );
01934 notein(chLine);
01935 }
01936
01937
01938 if( hmi.hmitot/SDIV(thermal.power) > 0.01 )
01939 {
01940 sprintf( chLine, " H- absorption heating was %.1f%% of the total.",
01941 hmi.hmitot/SDIV(thermal.power)*100. );
01942 notein(chLine);
01943 }
01944
01945
01946 if( co.lgH2Ozer )
01947 {
01948 sprintf( chLine, " Water destruction rate zero." );
01949 notein(chLine);
01950 }
01951
01952
01953 if( atoms.nNegOI > 0 )
01954 {
01955 sprintf( chLine, " C-O I negative level populations %ld times.",
01956 atoms.nNegOI );
01957 caunin(chLine);
01958 }
01959
01960
01961
01962 small = 0.;
01963 imas = 0;
01964 isav = 0;
01965 j = 0;
01966 for( nelem=0; nelem<LIMELM; ++nelem )
01967 {
01968 if( dense.lgElmtOn[nelem] )
01969 {
01970
01971 for( ipLo=ipH2p; ipLo < (iso.numLevels_local[ipH_LIKE][nelem] - 1); ipLo++ )
01972 {
01973
01974 for( ipHi=ipLo + 1; ipHi < iso.numLevels_local[ipH_LIKE][nelem]; ipHi++ )
01975 {
01976 if( Transitions[ipH_LIKE][nelem][ipHi][ipLo].Emis->Aul <= iso.SmallA )
01977 continue;
01978
01979 if( Transitions[ipH_LIKE][nelem][ipHi][ipLo].Emis->TauIn < (realnum)small )
01980 {
01981 small = Transitions[ipH_LIKE][nelem][ipHi][ipLo].Emis->TauIn;
01982 imas = ipHi;
01983 j = ipLo;
01984 isav = nelem;
01985 }
01986 }
01987 }
01988 }
01989 }
01990
01991 if( small < -0.05 )
01992 {
01993 sprintf( chLine,
01994 " !Some hydrogenic lines mased, species was %2s%2ld, smallest tau was %.2e, transition %li-%li",
01995 elementnames.chElementSym[isav],
01996 isav+1,small, imas , j );
01997 bangin(chLine);
01998 }
01999
02000
02001 if( opac.lgOpacNeg )
02002 {
02003 sprintf( chLine, " !Some opacities were negative - the SET NEGOPC command will save which ones." );
02004 bangin(chLine);
02005 }
02006
02007
02008 small = 0.;
02009 imas = 0;
02010 isav = 0;
02011 for( nelem=0; nelem<LIMELM; ++nelem )
02012 {
02013 if( dense.lgElmtOn[nelem] )
02014 {
02015
02016 for( i=0; i < iso.numLevels_local[ipH_LIKE][nelem]; i++ )
02017 {
02018 if( opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][nelem][i]-1] < -0.001 )
02019 {
02020 small = MIN2(small,(double)opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][nelem][i]-1]);
02021 imas = i;
02022 isav = nelem;
02023 }
02024 }
02025 }
02026 }
02027
02028 if( small < -0.05 )
02029 {
02030 sprintf( chLine, " !Some hydrogenic (%2s%2ld) continua optical depths were negative; smallest=%.2e level=%3ld",
02031 elementnames.chElementSym[isav],
02032 isav+1,
02033 small, imas );
02034 bangin(chLine);
02035 }
02036
02037
02038 nneg = 0;
02039 tauneg = 0.;
02040 freqn = 0.;
02041 for( i=0; i < rfield.nflux; i++ )
02042 {
02043 if( opac.TauAbsGeo[0][i] < -0.001 )
02044 {
02045 nneg += 1;
02046
02047 if( nneg == 1 )
02048 freqn = rfield.anu[i];
02049 tauneg = MIN2(tauneg,(double)opac.TauAbsGeo[0][i]);
02050 }
02051 }
02052
02053 if( nneg > 0 )
02054 {
02055 sprintf( chLine, " !Some continuous optical depths <0. The lowest freq was %.3e Ryd, and a total of%4ld",
02056 freqn, nneg );
02057 bangin(chLine);
02058 sprintf( chLine, " !The smallest optical depth was %.2e",
02059 tauneg );
02060 bangin(chLine);
02061 }
02062
02063
02064 if( opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1] > 0.05 )
02065 {
02066 sprintf( chLine, " The Balmer continuum optical depth was %.2e.",
02067 opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1] );
02068 notein(chLine);
02069 }
02070
02071
02072 if( opac.stimax[0] > 0.02 && opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][ipH1s]-1] > 0.2 )
02073 {
02074 sprintf( chLine, " The Lyman continuum stimulated emission correction to optical depths reached %.2e.",
02075 opac.stimax[0] );
02076 notein(chLine);
02077 }
02078 else if( opac.stimax[1] > 0.02 && opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][2]-1] > 0.1 )
02079 {
02080 sprintf( chLine, " The Balmer continuum stimulated emission correction to optical depths reached %.2e.",
02081 opac.stimax[1] );
02082 notein(chLine);
02083 }
02084
02085
02086 if( opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][3]-1] > 0.2 )
02087 {
02088 sprintf( chLine,
02089 " The Paschen continuum optical depth was %.2e.",
02090 opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][3]-1] );
02091 notein(chLine);
02092 }
02093
02094
02095 if( opac.TauAbsGeo[0][0] > 1. )
02096 {
02097 sprintf( chLine,
02098 " The continuum optical depth at the lowest energy considered (%.3e Ryd) was %.3e.",
02099 rfield.anu[0], opac.TauAbsGeo[0][0] );
02100 notein(chLine);
02101 }
02102
02103
02104
02105 if( colden.colden[ipCOL_H0]*7e-24 > 0.01 )
02106 {
02107 sprintf( chLine,
02108 " The optical depth to Rayleigh scattering at 1300A is %.2e",
02109 colden.colden[ipCOL_H0]*6.71e-24 );
02110 notein(chLine);
02111 }
02112
02113 if( colden.colden[ipCOL_H2p]*7e-18 > 0.1 )
02114 {
02115 sprintf( chLine,
02116 " !The optical depth to the H2+ molecular ion is %.2e",
02117 colden.colden[ipCOL_H2p]*7e-18 );
02118 bangin(chLine);
02119 }
02120 else if( colden.colden[ipCOL_H2p]*7e-18 > 0.01 )
02121 {
02122 sprintf( chLine,
02123 " The optical depth to the H2+ molecular ion is %.2e",
02124 colden.colden[ipCOL_H2p]*7e-18 );
02125 notein(chLine);
02126 }
02127
02128
02129 if( opac.thmin > 0.1 )
02130 {
02131 sprintf( chLine,
02132 " !Optical depth to negative hydrogen ion is %.2e",
02133 opac.thmin );
02134 bangin(chLine);
02135 }
02136 else if( opac.thmin > 0.01 )
02137 {
02138 sprintf( chLine,
02139 " Optical depth to negative hydrogen ion is %.2e",
02140 opac.thmin );
02141 notein(chLine);
02142 }
02143
02144
02145
02146
02147 if( ionbal.ifail > 0 && ionbal.ifail <= 10 )
02148 {
02149 sprintf( chLine,
02150 " 3 body recombination coefficient outside range %ld", ionbal.ifail );
02151 notein(chLine);
02152 }
02153 else if( ionbal.ifail > 10 )
02154 {
02155 sprintf( chLine,
02156 " C-3 body recombination coefficient outside range %ld", ionbal.ifail );
02157 caunin(chLine);
02158 }
02159
02160
02161 if( phycon.TEnerDen < 2.6 )
02162 {
02163 sprintf( chLine,
02164 " !Incident radiation field energy density is less than 2.7K. Add background with CMB command." );
02165 bangin(chLine);
02166 }
02167
02168
02169 if( !rfield.lgCMB_set )
02170 {
02171 sprintf( chLine,
02172 " !The CMB was not included. This is added with the CMB command." );
02173 bangin(chLine);
02174 }
02175
02176
02177 if( rfield.lgHabing )
02178 {
02179 sprintf( chLine,
02180 " !The intensity of the incident radiation field is less than 10 times the Habing diffuse ISM field. Is this OK?" );
02181 bangin(chLine);
02182 sprintf( chLine,
02183 " ! Consider adding diffuse ISM emission with TABLE ISM command." );
02184 bangin(chLine);
02185 }
02186
02187
02188
02189
02190 if( dense.lgElmtOn[ipOXYGEN] && dense.lgElmtOn[ipCARBON] )
02191 {
02192 if( dense.gas_phase[ipCARBON]/dense.gas_phase[ipOXYGEN] > 1. )
02193 {
02194 sprintf( chLine, " !The C/O abundance ratio, %.1f, is greater than unity. The chemistry will be carbon dominated.",
02195 dense.gas_phase[ipCARBON]/dense.gas_phase[ipOXYGEN] );
02196 bangin(chLine);
02197 }
02198 }
02199
02200
02201 #if 1
02202
02203 if( hmi.BiggestH2 > 0.1 )
02204 {
02205 sprintf( chLine, " !The fraction of %s in %s reached %.1f%% at some point in the cloud.",
02206 "H ",
02207 "H2 ",
02208 hmi.BiggestH2*100. );
02209 bangin(chLine);
02210 }
02211 else if( hmi.BiggestH2>0.01 )
02212 {
02213 sprintf( chLine, " The fraction of %s in %s reached %.2f%% at some point in the cloud.",
02214 "H ",
02215 "H2 ",
02216 hmi.BiggestH2*100. );
02217 notein(chLine);
02218 }
02219 else if( hmi.BiggestH2 > 1e-3 )
02220 {
02221 sprintf( chLine, " The fraction of %s in %s reached %.3f%% at some point in the cloud.",
02222 "H ",
02223 "H2 ",
02224 hmi.BiggestH2*100. );
02225 notein(chLine);
02226 }
02227 #endif
02228
02229 lgLots_of_moles = false;
02230 lgLotsSolids = false;
02231
02232 for( i=0; i<mole.num_comole_calc; ++i )
02233 {
02234 if( COmole[i]->n_nuclei <= 1)
02235 continue;
02236
02237 if( COmole[i]->xMoleFracMax > 0.1 )
02238 {
02239 sprintf( chLine, " !The fraction of %s in %s reached %.1f%% at some point in the cloud.",
02240 elementnames.chElementSym[COmole[i]->nelem_hevmol],
02241 COmole[i]->label,
02242 COmole[i]->xMoleFracMax*100. );
02243 bangin(chLine);
02244 lgLots_of_moles = true;
02245
02246 if( !COmole[i]->lgGas_Phase )
02247 lgLotsSolids = true;
02248 }
02249 else if( COmole[i]->xMoleFracMax>0.01 )
02250 {
02251 sprintf( chLine, " The fraction of %s in %s reached %.2f%% at some point in the cloud.",
02252 elementnames.chElementSym[COmole[i]->nelem_hevmol],
02253 COmole[i]->label,
02254 COmole[i]->xMoleFracMax*100. );
02255 notein(chLine);
02256 lgLots_of_moles = true;
02257
02258 if( !COmole[i]->lgGas_Phase )
02259 lgLotsSolids = true;
02260 }
02261 else if( COmole[i]->xMoleFracMax > 1e-3 )
02262 {
02263 sprintf( chLine, " The fraction of %s in %s reached %.3f%% at some point in the cloud.",
02264 elementnames.chElementSym[COmole[i]->nelem_hevmol],
02265 COmole[i]->label,
02266 COmole[i]->xMoleFracMax*100. );
02267 notein(chLine);
02268
02269 if( !COmole[i]->lgGas_Phase )
02270 lgLotsSolids = true;
02271 }
02272 }
02273
02274
02275 if( lgLots_of_moles )
02276 {
02277
02278 for(nelem=ipHYDROGEN; nelem<LIMELM; ++nelem )
02279 {
02280
02281 if( !dense.lgElmtOn[nelem] && mole.lgElem_in_chemistry[nelem] )
02282 {
02283
02284 sprintf( chLine,
02285 " C-Molecules are important, but %s, part of the chemistry network, is turned off.",
02286 elementnames.chElementName[nelem] );
02287 caunin(chLine);
02288 }
02289 # if 0
02290
02291 for( i=NUM_HEAVY_MOLEC+NUM_ELEMENTS; i<NUM_COMOLE_CALC; ++i )
02292 {
02293 if( nelem==COmole[i].nelem_hevmol )
02294 {
02295
02296 sprintf( chLine,
02297 " C-Molecules are important, but %s, part of the chemistry network, is turned off.",
02298 elementnames.chElementName[nelem] );
02299 caunin(chLine);
02300 }
02301 }
02302 }
02303 # endif
02304 }
02305 }
02306
02307
02308
02309 if( lgLotsSolids )
02310 {
02311 sprintf( chLine, " !A significant amount of molecules condensed onto grain surfaces." );
02312 bangin(chLine);
02313 sprintf( chLine, " !These are the molecular species with \"grn\" above." );
02314 bangin(chLine);
02315 }
02316
02317
02318 if( rfield.EnergyBremsThin > 0.09 )
02319 {
02320 sprintf( chLine, " !The cloud is optically thick at optical wavelengths, extending to %.3e Ryd =%.3eA",
02321 rfield.EnergyBremsThin, RYDLAM/rfield.EnergyBremsThin );
02322 bangin(chLine);
02323 }
02324 else if( rfield.EnergyBremsThin > 0.009 )
02325 {
02326 sprintf( chLine, " The continuum of the computed structure may be optically thick in the near infrared." );
02327 notein(chLine);
02328 }
02329
02330
02331 if( radius.Radius > 1e23 && radius.Radius/radius.rinner > 10. )
02332 {
02333 sprintf( chLine, " Is an outer radius of %.2e reasonable?",
02334 radius.Radius );
02335 notein(chLine);
02336 }
02337
02338
02339 if( rt.lgMaserCapHit )
02340 {
02341 sprintf( chLine, " Laser maser optical depths capped in RT_line_one_tauinc." );
02342 notein(chLine);
02343 }
02344
02345
02346 if( rt.lgMaserSetDR )
02347 {
02348 sprintf( chLine, " !Line maser set zone thickness in some zones." );
02349 bangin(chLine);
02350 }
02351
02352
02353
02354 if( (pressure.lgPradCap && (strcmp(dense.chDenseLaw,"CPRE") == 0)) &&
02355 pressure.lgPres_radiation_ON )
02356 {
02357 sprintf( chLine, " Radiation pressure kept below gas pressure on this iteration." );
02358 notein(chLine);
02359 }
02360
02361 if( pressure.RadBetaMax > 0.25 )
02362 {
02363 if( pressure.ipPradMax_line == 0 )
02364 {
02365 sprintf( chLine,
02366 " !The ratio of radiation to gas pressure reached %.2e at zone %li. Caused by Lyman alpha.",
02367 pressure.RadBetaMax,
02368 pressure.ipPradMax_nzone);
02369 bangin(chLine);
02370 }
02371 else
02372 {
02373 sprintf( chLine,
02374 " !The ratio of radiation to gas pressure reached %.2e at zone %li. "
02375 "Caused by line number %ld, label %s",
02376 pressure.RadBetaMax,
02377 pressure.ipPradMax_nzone,
02378 pressure.ipPradMax_line,
02379 pressure.chLineRadPres );
02380 bangin(chLine);
02381 }
02382 }
02383
02384 else if( pressure.RadBetaMax > 0.025 )
02385 {
02386 if( pressure.ipPradMax_line == 0 )
02387 {
02388 sprintf( chLine,
02389 " The ratio of radiation to gas pressure reached %.2e at zone %li. Caused by Lyman alpha.",
02390 pressure.RadBetaMax,
02391 pressure.ipPradMax_nzone);
02392 notein(chLine);
02393 }
02394 else
02395 {
02396 sprintf( chLine,
02397 " The ratio of radiation to gas pressure reached %.2e at zone %li. "
02398 "Caused by line number %ld, label %s",
02399 pressure.RadBetaMax,
02400 pressure.ipPradMax_nzone,
02401 pressure.ipPradMax_line,
02402 pressure.chLineRadPres );
02403 notein(chLine);
02404 }
02405 }
02406
02407 if( opac.telec >= 5. )
02408 {
02409 sprintf( chLine, " W-The model is optically thick to electron "
02410 "scattering; tau=%.2e Cloudy is NOT intended for this regime.",
02411 opac.telec );
02412 warnin(chLine);
02413 }
02414 else if( opac.telec > 2.0 )
02415 {
02416 sprintf( chLine, " C-The model is moderately optically thick to electron scattering; tau=%.1f",
02417 opac.telec );
02418 caunin(chLine);
02419 }
02420 else if( opac.telec > 0.1 )
02421 {
02422 sprintf( chLine, " !The model has modest optical depth to electron scattering; tau=%.2f",
02423 opac.telec );
02424 bangin(chLine);
02425 }
02426 else if( opac.telec > 0.01 )
02427 {
02428 sprintf( chLine, " The optical depth to electron scattering is %.3f",
02429 opac.telec );
02430 notein(chLine);
02431 }
02432
02433
02434 if( HFLines[0].Emis->TauIn > 0.5 )
02435 {
02436 sprintf( chLine, " !The optical depth in the H I 21 cm line is %.2e",HFLines[0].Emis->TauIn );
02437 bangin(chLine);
02438 }
02439
02440
02441
02442 if( nWindLine==0 )
02443 {
02444
02445 sprintf( chLine, " !The level2 lines are disabled. UV pumping of excited levels within ground terms is not treated." );
02446 bangin(chLine);
02447 }
02448
02449
02450 for( nelem=0; nelem < LIMELM; nelem++ )
02451 {
02452 if( dense.lgElmtOn[nelem] )
02453 {
02454 if( Transitions[ipH_LIKE][nelem][ipH3p][ipH2s].Emis->TauIn > 0.2 )
02455 {
02456 differ = fabs(1.-Transitions[ipH_LIKE][nelem][ipH3p][ipH2s].Emis->TauIn*
02457 rt.DoubleTau/Transitions[ipH_LIKE][nelem][ipH3p][ipH2s].Emis->TauTot)*100.;
02458
02459
02460
02461
02462 if( ((iterations.lgLastIt && Transitions[ipH_LIKE][nelem][ipH3p][ipH2s].Emis->TauIn > 0.8) &&
02463 differ > 20.) && wind.lgStatic() )
02464 {
02465 sprintf( chLine,
02466 " C-This is the last iteration and %2s%2ld Bal(a) optical depth"
02467 " changed by%6.1f%% (was %.2e). Try another iteration.",
02468 elementnames.chElementSym[nelem],
02469 nelem+1, differ,
02470 Transitions[ipH_LIKE][nelem][ipH3p][ipH2s].Emis->TauTot );
02471 caunin(chLine);
02472 iterations.lgIterAgain = true;
02473 }
02474
02475
02476 if( Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->TauIn > 0. )
02477 {
02478 differ = fabs(1.-Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->TauIn*
02479 rt.DoubleTau/Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->TauTot)*100.;
02480
02481
02482
02483
02484 if( ((iterations.lgLastIt && Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->TauIn > 0.8) &&
02485 differ > 25.) && wind.lgStatic() )
02486 {
02487 sprintf( chLine,
02488 " C-This is the last iteration and %2s%2ld Ly(a) optical depth"
02489 " changed by%7.0f%% (was %.2e). Try another iteration.",
02490 elementnames.chElementSym[nelem],
02491 nelem+1,differ, Transitions[ipH_LIKE][nelem][ipH2p][ipH1s].Emis->TauTot );
02492 caunin(chLine);
02493 iterations.lgIterAgain = true;
02494 }
02495 }
02496 }
02497 }
02498 }
02499
02500
02501 if( radius.Radius/radius.rinner > 2. && !geometry.lgSphere )
02502 {
02503 sprintf( chLine, " C-R(out)/R(in)=%.2e and SPHERE was not set.",
02504 radius.Radius/radius.rinner );
02505 caunin(chLine);
02506 }
02507
02508
02509 if( iterations.lgLastIt && !opac.lgCaseB )
02510 {
02511
02512
02513 if( rfield.nflux > iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][ipH1s] )
02514 {
02515 if( opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][ipH1s]-1] < 2. &&
02516 opac.TauAbsGeo[1][iso.ipIsoLevNIonCon[ipH_LIKE][ipHYDROGEN][ipH1s]-1] > 2. )
02517 {
02518 sprintf( chLine, " C-The H Lyman continuum is thin, and I assumed"
02519 " that it was thick. Try another iteration." );
02520 caunin(chLine);
02521 iterations.lgIterAgain = true;
02522 }
02523 }
02524
02525
02526 if( rfield.nflux > iso.ipIsoLevNIonCon[ipH_LIKE][ipHELIUM][ipH1s] && dense.lgElmtOn[ipHELIUM] )
02527 {
02528 if( (opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipH_LIKE][ipHELIUM][ipH1s]-1] < 2. &&
02529 opac.TauAbsGeo[1][iso.ipIsoLevNIonCon[ipH_LIKE][ipHELIUM][ipH1s]-1] > 2.) )
02530 {
02531 sprintf( chLine,
02532 " C-The He II continuum is thin and I assumed that it was thick."
02533 " Try another iteration." );
02534 caunin(chLine);
02535 iterations.lgIterAgain = true;
02536 }
02537 }
02538
02539 if( rfield.nflux > iso.ipIsoLevNIonCon[ipHE_LIKE][ipHELIUM][0] && dense.lgElmtOn[ipHELIUM] )
02540 {
02541 if( (opac.TauAbsGeo[0][iso.ipIsoLevNIonCon[ipHE_LIKE][ipHELIUM][0]-1] < 2. &&
02542 opac.TauAbsGeo[1][iso.ipIsoLevNIonCon[ipHE_LIKE][ipHELIUM][0]-1] > 2.) )
02543 {
02544 sprintf( chLine,
02545 " C-The He I continuum is thin and I assumed that it was thick."
02546 " Try another iteration." );
02547 caunin(chLine);
02548 iterations.lgIterAgain = true;
02549 }
02550 }
02551 }
02552
02553
02554 if( iteration > 1 )
02555 {
02556 if( colden.colden_old[ipCOL_HTOT] <= 0. )
02557 {
02558 fprintf( ioQQQ, " colden_old is insane in PrtComment.\n" );
02559 ShowMe();
02560 cdEXIT(EXIT_FAILURE);
02561 }
02562
02563 differ = fabs(1.-colden.colden[ipCOL_HTOT]/
02564 colden.colden_old[ipCOL_HTOT]);
02565
02566 if( differ > 0.1 && differ <= 0.3 )
02567 {
02568 sprintf( chLine,
02569 " The H column density changed by %.2e%% between this and previous iteration.",
02570 differ*100. );
02571 notein(chLine);
02572 }
02573
02574 else if( differ > 0.3 )
02575 {
02576 if( iterations.lgLastIt )
02577 {
02578 sprintf( chLine,
02579 " C-The H column density changed by %.2e%% and this is the last iteration. What happened?",
02580 differ*100. );
02581 caunin(chLine);
02582 }
02583 else
02584 {
02585 sprintf( chLine,
02586 " !The H column density changed by %.2e%% What happened?",
02587 differ*100. );
02588 bangin(chLine);
02589 }
02590 }
02591
02592
02593 if( (colden.colden[ipCOL_H2g]+colden.colden[ipCOL_H2s])/SDIV(colden.colden[ipCOL_HTOT]) > 1e-5 )
02594 {
02595 differ = fabs(1.-colden.colden[ipCOL_H2g]/
02596 SDIV(colden.colden_old[ipCOL_H2g]));
02597
02598 if( differ > 0.1 && differ <= 0.3 )
02599 {
02600 sprintf( chLine,
02601 " The H2 column density changed by %.2e%% between this and previous iteration.",
02602 differ*100. );
02603 notein(chLine);
02604 }
02605
02606 else if( differ > 0.3 )
02607 {
02608 if( iterations.lgLastIt )
02609 {
02610 sprintf( chLine,
02611 " C-The H2 column density changed by %.2e%% and this is the last iteration. What happened?",
02612 differ*100. );
02613 caunin(chLine);
02614 }
02615 else
02616 {
02617 sprintf( chLine,
02618 " !The H2 column density changed by %.2e%% What happened?",
02619 differ*100. );
02620 bangin(chLine);
02621 }
02622 }
02623 }
02624 }
02625
02626
02627 differ = fabs(1.-Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->TauIn/
02628 SDIV(Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->TauTot))*100.;
02629
02630 if( iterations.lgLastIt && (pressure.RadBetaMax > 0.1) &&
02631 (differ > 50.) && (pressure.ipPradMax_line == 1) && (pressure.lgPres_radiation_ON) &&
02632 wind.lgStatic() )
02633 {
02634 sprintf( chLine, " C-This is the last iteration, radiation pressure was significant, and the L-a optical depth changed by %7.2f%% (was %.2e)",
02635 differ, Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->TauTot );
02636 caunin(chLine);
02637 }
02638
02639
02640
02641 if( iterations.lgLastIt &&
02642 ( Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->TauTot * 1.02 -
02643 Transitions[ipH_LIKE][ipHYDROGEN][ipH2p][ipH1s].Emis->TauIn ) < 0. )
02644 {
02645 sprintf( chLine, " C-The Lya optical depth scale was overrun and this is the last iteration - Tspin(21 cm) is not valid." );
02646 caunin(chLine);
02647 sprintf( chLine, " C-Another iteration is needed for Tspin(21 cm) to be valid." );
02648 caunin(chLine);
02649 }
02650
02651
02652 if( pressure.lgPradDen )
02653 {
02654 sprintf( chLine, " Line radiation pressure capped by thermalization length." );
02655 notein(chLine);
02656 }
02657
02658
02659 nline = MIN2(conv.nTeFail,10);
02660 if( conv.nTeFail != 0 )
02661 {
02662 long int _o;
02663 if( conv.failmx < 0.1 )
02664 {
02665 _o = sprintf( chLine, " There were %ld minor temperature failures. zones:",
02666 conv.nTeFail );
02667
02668
02669 for( i=0; i < nline; i++ )
02670 {
02671 _o += sprintf( chLine+_o, " %ld", conv.ifailz[i] );
02672 }
02673 notein(chLine);
02674 }
02675 else
02676 {
02677 _o = sprintf( chLine,
02678 " !There were %ld temperature failures, and some were large. The largest was %.1f%%. What happened?",
02679 conv.nTeFail, conv.failmx*100. );
02680 bangin(chLine);
02681
02682
02683
02684 _o = sprintf( chLine , " !The zones were" );
02685 for( i=0; i < nline; i++ )
02686 {
02687 _o += sprintf( chLine+_o, " %ld", conv.ifailz[i] );
02688 }
02689 bangin(chLine);
02690
02691 if( struc.testr[0] > 8e4 && phycon.te < 6e5 )
02692 {
02693 sprintf( chLine, " !I think they may have been caused by the change from hot to nebular gas phase. The physics of this is unclear." );
02694 bangin(chLine);
02695 }
02696 }
02697 }
02698
02699
02700 big_ion_jump = 0.;
02701 j = 0;
02702 for( i=1; i < nzone; i++ )
02703 {
02704 big = fabs(1.-struc.testr[i-1]/struc.testr[i]);
02705 if( big > big_ion_jump )
02706 {
02707 j = i;
02708 big_ion_jump = big;
02709 }
02710 }
02711
02712 if( big_ion_jump > 0.2 )
02713 {
02714
02715 if( j < 1 )
02716 {
02717 fprintf( ioQQQ, " j too small big jump check\n" );
02718 ShowMe();
02719 cdEXIT(EXIT_FAILURE);
02720 }
02721
02722 if( big_ion_jump > 0.4 )
02723 {
02724 sprintf( chLine, " C-A temperature discontinuity occurred at zone %ld from %.2eK to %.2eK.",
02725 j, struc.testr[j-1], struc.testr[j] );
02726 caunin(chLine);
02727
02728
02729
02730
02731 if( struc.testr[j]>100. && struc.testr[j] < 1000. )
02732 {
02733 sprintf( chLine, " C-This was probably due to a thermal front." );
02734 caunin(chLine);
02735 }
02736 }
02737 else if( big_ion_jump > 0.2 )
02738 {
02739 sprintf( chLine, " !A temperature discontinuity occurred at zone %ld from %.2eK to %.2eK.",
02740 j, struc.testr[j-1], struc.testr[j] );
02741 bangin(chLine);
02742
02743
02744
02745
02746 if( struc.testr[j]>100. && struc.testr[j] < 1000. )
02747 {
02748 sprintf( chLine, " !This was probably due to a thermal front." );
02749 bangin(chLine);
02750 }
02751 }
02752 }
02753
02754
02755 if( fabs(conv.BigEdenError) > conv.EdenErrorAllowed )
02756 {
02757
02758 if( fabs(conv.BigEdenError) > conv.EdenErrorAllowed*20. && dense.nzEdenBad !=
02759 nzone )
02760 {
02761 sprintf( chLine, " W-The local error in the electron density reached %.1f%% at zone %ld",
02762 conv.BigEdenError*100, dense.nzEdenBad );
02763 warnin(chLine);
02764 }
02765 else if( fabs(conv.BigEdenError) > conv.EdenErrorAllowed*5. )
02766 {
02767 sprintf( chLine, " C-The local error in the electron density reached %.1f%% at zone %ld",
02768 conv.BigEdenError*100, dense.nzEdenBad );
02769 caunin(chLine);
02770 }
02771 else
02772 {
02773 sprintf( chLine, " The local error in the electron density reached %.1f%% at zone %ld",
02774 conv.BigEdenError*100, dense.nzEdenBad );
02775 notein(chLine);
02776 }
02777 }
02778
02779
02780 big_ion_jump = 0.;
02781 j = 0;
02782 for( i=1; i < (nzone - 1); i++ )
02783 {
02784 big = fabs( (struc.testr[i-1] - struc.testr[i])/struc.testr[i] );
02785 bigm = fabs( (struc.testr[i] - struc.testr[i+1])/struc.testr[i] );
02786
02787
02788 rel = ( (struc.testr[i-1] - struc.testr[i])/struc.testr[i])*
02789 ( (struc.testr[i] - struc.testr[i+1])/struc.testr[i] );
02790
02791 if( rel < 0. && MIN2( bigm , big ) > big_ion_jump )
02792 {
02793 j = i;
02794 big_ion_jump = MIN2( bigm , big );
02795 }
02796 }
02797
02798 if( big_ion_jump > 0.1 )
02799 {
02800
02801 if( j < 1 )
02802 {
02803 fprintf( ioQQQ, " j too small bigjump2 check\n" );
02804 ShowMe();
02805 cdEXIT(EXIT_FAILURE);
02806 }
02807
02808 if( big_ion_jump > 0.3 )
02809 {
02810 sprintf( chLine,
02811 " C-A temperature oscillation occurred at zone%4ld by%5.0f%% from %.2e to %.2e to %.2e",
02812 j, big_ion_jump*100., struc.testr[j-1], struc.testr[j], struc.testr[j+1] );
02813 caunin(chLine);
02814 }
02815 else if( big_ion_jump > 0.1 )
02816 {
02817 sprintf( chLine,
02818 " !A temperature oscillation occurred at zone%4ld by%5.0f%% from %.2e to %.2e to %.2e",
02819 j, big_ion_jump*100., struc.testr[j-1], struc.testr[j], struc.testr[j+1] );
02820 bangin(chLine);
02821 }
02822 }
02823
02824
02825 if( strcmp(dense.chDenseLaw,"CDEN") == 0 )
02826 {
02827 j = 0;
02828 big_ion_jump = 0.;
02829 for( i=1; i < (nzone - 1); i++ )
02830 {
02831 big = (struc.ednstr[i-1] - struc.ednstr[i])/struc.ednstr[i];
02832 if( fabs(big) < conv.EdenErrorAllowed )
02833 big = 0.;
02834 bigm = (struc.ednstr[i] - struc.ednstr[i+1])/struc.ednstr[i];
02835 if( fabs(bigm) < conv.EdenErrorAllowed )
02836 bigm = 0.;
02837 if( big*bigm < 0. &&
02838 fabs(struc.ednstr[i-1]-struc.ednstr[i])/struc.ednstr[i] > big_ion_jump )
02839 {
02840 j = i;
02841 big_ion_jump = fabs(struc.ednstr[i-1]-struc.ednstr[i])/
02842 struc.ednstr[i];
02843 }
02844 }
02845
02846
02847
02848 if( big_ion_jump > conv.EdenErrorAllowed*3. )
02849 {
02850 if( j < 1 )
02851 {
02852 fprintf( ioQQQ, " j too small bigjump3 check\n" );
02853 ShowMe();
02854 cdEXIT(EXIT_FAILURE);
02855 }
02856
02857 if( big_ion_jump > conv.EdenErrorAllowed*10. )
02858 {
02859 sprintf( chLine, " C-An electron density oscillation occurred at zone%4ld by%5.0f%% from %.2e to %.2e to %.2e",
02860 j, big_ion_jump*100., struc.ednstr[j-1], struc.ednstr[j],
02861 struc.ednstr[j+1] );
02862 caunin(chLine);
02863 }
02864 else if( big_ion_jump > conv.EdenErrorAllowed*3. )
02865 {
02866 sprintf( chLine, " !An electron density oscillation occurred at zone%4ld by%5.0f%% from %.2e to %.2e to %.2e",
02867 j, big_ion_jump*100., struc.ednstr[j-1], struc.ednstr[j],
02868 struc.ednstr[j+1] );
02869 bangin(chLine);
02870 }
02871 }
02872 }
02873
02874
02875
02876 prt_smooth_predictions();
02877
02878
02879
02880
02881
02882
02883
02884
02885 i = cdLine( "Unit" , 1 , &rate , &absint );
02886 ASSERT( i> 0 );
02887
02888
02889 VolComputed = LineSv[i].SumLine[0] / 1e-10;
02890
02891
02892 if( radius.Radius/radius.rinner > 1.0001 )
02893 {
02894
02895
02896
02897
02898 VolExpected = geometry.covaper*geometry.FillFac*radius.rinner/(geometry.iEmissPower+1)*
02899 ( powi( radius.Radius/radius.rinner,geometry.iEmissPower+1 ) - 1. );
02900 }
02901 else
02902 {
02903
02904
02905 VolExpected = geometry.covaper*geometry.FillFac*(radius.depth-DEPTH_OFFSET);
02906 }
02907
02908
02909 error = fabs(VolComputed - VolExpected)/SDIV(VolExpected);
02910
02911
02912
02913 if( radius.lgCylnOn || geometry.filpow!=0. )
02914 {
02915 error = 0.;
02916 }
02917
02918
02919 if( error > 0.001 && !lgAbort )
02920 {
02921 sprintf( chLine,
02922 " W-PrtComment insanity - Line unit integration did not verify \n");
02923 warnin(chLine);
02924 fprintf( ioQQQ,
02925 " PROBLEM PrtComment insanity - Line unit integration did not verify \n");
02926 fprintf( ioQQQ,
02927 " expected, derived vols were %g %g \n",
02928 VolExpected , VolComputed );
02929 fprintf( ioQQQ,
02930 " relative difference is %g, ratio is %g.\n",error,VolComputed/VolExpected);
02931 TotalInsanity();
02932 }
02933
02934
02935
02936
02937
02938
02939
02940
02941
02942
02943
02944
02945 ConComputed = rfield.ConInterOut[rfield.nflux]/ 1e-10;
02946
02947
02948 ConComputed /= ( (1. + geometry.covrt)/2. );
02949
02950
02951 if( radius.Radius/radius.rinner < 1.001 )
02952 {
02953
02954 ConExpected = (radius.depth-DEPTH_OFFSET)*geometry.FillFac;
02955 }
02956 else
02957 {
02958
02959 ConExpected = radius.rinner*geometry.FillFac * (1. - radius.rinner/radius.Radius );
02960 }
02961
02962 ASSERT( ConExpected > 0. );
02963
02964
02965 error = fabs(ConComputed - ConExpected)/ConExpected;
02966
02967
02968
02969 if( radius.lgCylnOn || geometry.filpow!=0. )
02970 {
02971 error = 0.;
02972 }
02973
02974
02975
02976
02977
02978
02979
02980 if( error > 0.001 && !lgAbort)
02981 {
02982 sprintf( chLine,
02983 " W-PrtComment insanity - Continuum unit integration did not verify \n");
02984 warnin(chLine);
02985 fprintf( ioQQQ," PROBLEM PrtComment insanity - Continuum unit integration did not verify \n");
02986 fprintf( ioQQQ," exact vol= %g, derived vol= %g relative difference is %g \n",
02987 ConExpected , ConComputed ,error);
02988 fprintf( ioQQQ," ConInterOut= %g, \n",
02989 rfield.ConInterOut[rfield.nflux]);
02990 TotalInsanity();
02991 }
02992
02993
02994 cdNwcns(&lgAbort_flag,&nw,&nc,&nn,&ns,&i,&j,&dum1,&dum2);
02995
02996 warnings.lgWarngs = ( nw > 0 );
02997 warnings.lgCautns = ( nc > 0 );
02998
02999 if( called.lgTalk )
03000 {
03001
03002 fprintf( ioQQQ, " %s\n", input.chTitle );
03003
03004 cdReasonGeo(ioQQQ);
03005
03006 cdWarnings(ioQQQ);
03007
03008 cdCautions(ioQQQ);
03009
03010 cdSurprises(ioQQQ);
03011
03012 cdNotes(ioQQQ);
03013 }
03014
03015
03016 if( lgPrnErr )
03017 {
03018 cdWarnings(ioPrnErr);
03019 }
03020
03021 if( trace.lgTrace )
03022 {
03023 fprintf( ioQQQ, " PrtComment returns.\n" );
03024 }
03025 return;
03026 }
03027
03028
03029
03030 STATIC void chkCaHeps(double *totwid)
03031 {
03032 double conca,
03033 conalog ,
03034 conhe;
03035
03036 DEBUG_ENTRY( "chkCaHeps()" );
03037
03038 *totwid = 0.;
03039
03040 if( iso.n_HighestResolved_local[ipH_LIKE][ipHYDROGEN] +
03041 iso.nCollapsed_local[ipH_LIKE][ipHYDROGEN] >= 6 )
03042 {
03043
03044 long ip6p = iso.QuantumNumbers2Index[ipH_LIKE][ipHYDROGEN][6][1][2];
03045
03046 if( TauLines[ipT3969].Emis->TauIn > 0. &&
03047 Transitions[ipH_LIKE][ipHYDROGEN][ip6p][ipH2s].Emis->TauIn > 0. )
03048 {
03049
03050 conca = pow(6.1e-5* TauLines[ipT3969].Emis->TauIn,0.5);
03051 conalog = log((double)TauLines[ipT3969].Emis->TauIn);
03052 conalog = sqrt(MAX2(1., conalog));
03053 conca = MAX2(conalog,conca);
03054
03055 conalog = log((double)Transitions[ipH_LIKE][ipHYDROGEN][ip6p][ipH2s].Emis->TauIn);
03056 conalog = sqrt(MAX2(1.,conalog));
03057 conhe = pow(1.7e-6*Transitions[ipH_LIKE][ipHYDROGEN][ip6p][ipH2s].Emis->TauIn,0.5);
03058 conhe = MAX2(conalog, conhe);
03059
03060 *totwid = 10.*conhe + 1.6*conca;
03061 }
03062 }
03063 return;
03064 }
03065
03066
03067 STATIC void prt_smooth_predictions(void)
03068 {
03069 long int i,
03070 nzone_oscillation,
03071 nzone_ion_jump,
03072 nzone_den_jump,
03073 nelem,
03074 ion;
03075 double BigOscillation ,
03076 big_ion_jump,
03077 big_jump,
03078 rel,
03079 big,
03080 bigm;
03081
03082 char chLine[INPUT_LINE_LENGTH];
03083
03084 DEBUG_ENTRY( "prt_smooth_predictions()" );
03085
03086
03087 nzone_oscillation = 0;
03088 nzone_ion_jump = 0;
03089
03090 for( nelem=ipHYDROGEN; nelem<LIMELM; ++nelem )
03091 {
03092 if( dense.lgElmtOn[nelem] )
03093 {
03094 for( ion=0; ion<=nelem+1; ++ion)
03095 {
03096 BigOscillation = 0.;
03097 big_ion_jump = -15.;
03098
03099 for( i=1; i < (nzone - 1-(int)lgAbort); i++ )
03100 {
03101
03102
03103 if( struc.xIonDense[nelem][ion][i-1]/struc.gas_phase[nelem][i-1]>struc.dr_ionfrac_limit &&
03104 struc.xIonDense[nelem][ion][i ]/struc.gas_phase[nelem][i ]>struc.dr_ionfrac_limit &&
03105 struc.xIonDense[nelem][ion][i+1]/struc.gas_phase[nelem][i+1]>struc.dr_ionfrac_limit )
03106 {
03107
03108
03109 big = fabs( (struc.xIonDense[nelem][ion][i-1] - struc.xIonDense[nelem][ion][i])/struc.xIonDense[nelem][ion][i] );
03110 bigm = fabs( (struc.xIonDense[nelem][ion][i] - struc.xIonDense[nelem][ion][i+1])/struc.xIonDense[nelem][ion][i] );
03111
03112
03113 rel = ( (struc.xIonDense[nelem][ion][i-1] - struc.xIonDense[nelem][ion][i] )/struc.xIonDense[nelem][ion][i])*
03114 ( (struc.xIonDense[nelem][ion][i] - struc.xIonDense[nelem][ion][i+1])/struc.xIonDense[nelem][ion][i] );
03115
03116 if( rel < 0. && MIN2( bigm , big ) > BigOscillation )
03117 {
03118 nzone_oscillation = i;
03119 BigOscillation = MIN2( bigm , big );
03120 }
03121
03122
03123
03124
03125 rel = -log10( (struc.xIonDense[nelem][ion][i]/struc.gas_phase[nelem][i]) /
03126 (struc.xIonDense[nelem][ion][i+1]/struc.gas_phase[nelem][i+1] ) );
03127
03128 if( rel > big_ion_jump )
03129 {
03130 big_ion_jump = rel;
03131 nzone_ion_jump = i;
03132 }
03133 }
03134 }
03135
03136
03137
03138 if( BigOscillation > 0.2 )
03139 {
03140
03141 if( nzone_oscillation < 1 )
03142 {
03143 fprintf( ioQQQ, " nzone_oscillation too small bigjump2 check\n" );
03144 ShowMe();
03145 cdEXIT(EXIT_FAILURE);
03146 }
03147 if( BigOscillation > 3. )
03148 {
03149 sprintf( chLine,
03150 " W-An ionization oscillation occurred at zone %ld, elem %.2s%2li, by %.0f%% from %.2e to %.2e to %.2e",
03151 nzone_oscillation,
03152 elementnames.chElementSym[nelem], ion+1,
03153 BigOscillation*100.,
03154 struc.xIonDense[nelem][ion][nzone_oscillation-1]/struc.gas_phase[nelem][nzone_oscillation-1],
03155 struc.xIonDense[nelem][ion][nzone_oscillation]/struc.gas_phase[nelem][nzone_oscillation],
03156 struc.xIonDense[nelem][ion][nzone_oscillation+1]/struc.gas_phase[nelem][nzone_oscillation+1] );
03157 warnin(chLine);
03158 }
03159
03160 else if( BigOscillation > 0.7 )
03161 {
03162 sprintf( chLine,
03163 " C-An ionization oscillation occurred at zone %ld, elem %.2s%2li, by %.0f%% from %.2e to %.2e to %.2e",
03164 nzone_oscillation,
03165 elementnames.chElementSym[nelem], ion+1,
03166 BigOscillation*100.,
03167 struc.xIonDense[nelem][ion][nzone_oscillation-1]/struc.gas_phase[nelem][nzone_oscillation-1],
03168 struc.xIonDense[nelem][ion][nzone_oscillation]/struc.gas_phase[nelem][nzone_oscillation],
03169 struc.xIonDense[nelem][ion][nzone_oscillation+1]/struc.gas_phase[nelem][nzone_oscillation+1] );
03170 caunin(chLine);
03171 }
03172 else if( BigOscillation > 0.2 )
03173 {
03174 sprintf( chLine,
03175 " !An ionization oscillation occurred at zone %ld, elem %.2s%2li, by %.0f%% from %.2e to %.2e to %.2e",
03176 nzone_oscillation,
03177 elementnames.chElementSym[nelem], ion+1,
03178 BigOscillation*100.,
03179 struc.xIonDense[nelem][ion][nzone_oscillation-1]/struc.gas_phase[nelem][nzone_oscillation-1],
03180 struc.xIonDense[nelem][ion][nzone_oscillation]/struc.gas_phase[nelem][nzone_oscillation],
03181 struc.xIonDense[nelem][ion][nzone_oscillation+1]/struc.gas_phase[nelem][nzone_oscillation+1] );
03182 bangin(chLine);
03183 }
03184 }
03185
03186
03187
03188 big_ion_jump = pow(10., big_ion_jump );
03189 if( big_ion_jump > 1.5 && nzone_ion_jump > 0 )
03190 {
03191 if( big_ion_jump > 10. )
03192 {
03193 sprintf( chLine,
03194 " C-An ionization jump occurred at zone %ld, elem %.2s%2li, by %.0f%% from %.2e to %.2e to %.2e",
03195 nzone_ion_jump,
03196 elementnames.chElementSym[nelem], ion+1,
03197 big_ion_jump*100.,
03198 struc.xIonDense[nelem][ion][nzone_ion_jump-1]/struc.gas_phase[nelem][nzone_ion_jump-1],
03199 struc.xIonDense[nelem][ion][nzone_ion_jump]/struc.gas_phase[nelem][nzone_ion_jump],
03200 struc.xIonDense[nelem][ion][nzone_ion_jump+1]/struc.gas_phase[nelem][nzone_ion_jump+1] );
03201 caunin(chLine);
03202 }
03203 else
03204 {
03205 sprintf( chLine,
03206 " !An ionization jump occurred at zone %ld, elem %.2s%2li, by %.0f%% from %.2e to %.2e to %.2e",
03207 nzone_ion_jump,
03208 elementnames.chElementSym[nelem], ion+1,
03209 big_ion_jump*100.,
03210 struc.xIonDense[nelem][ion][nzone_ion_jump-1]/struc.gas_phase[nelem][nzone_ion_jump-1],
03211 struc.xIonDense[nelem][ion][nzone_ion_jump]/struc.gas_phase[nelem][nzone_ion_jump],
03212 struc.xIonDense[nelem][ion][nzone_ion_jump+1]/struc.gas_phase[nelem][nzone_ion_jump+1] );
03213 bangin(chLine);
03214 }
03215 }
03216 }
03217 }
03218 }
03219
03220 big_jump = -15;
03221 nzone_den_jump = 0;
03222
03223
03224 for( i=1; i < (nzone - 1 - (int)lgAbort); i++ )
03225 {
03226
03227 rel = fabs(log10( struc.gas_phase[ipHYDROGEN][i] /
03228 struc.gas_phase[ipHYDROGEN][i+1] ) );
03229
03230 if( rel > big_jump )
03231 {
03232 big_jump = rel;
03233 nzone_den_jump = i;
03234 }
03235 }
03236
03237
03238 big_jump = pow( 10., big_jump );
03239 if( big_jump > 1.2 )
03240 {
03241 if( big_jump > 3. )
03242 {
03243 sprintf( chLine,
03244 " C-The H density jumped at by %.0f%% at zone %ld, from %.2e to %.2e to %.2e",
03245 big_jump*100.,
03246 nzone_den_jump,
03247 struc.gas_phase[ipHYDROGEN][nzone_den_jump-1],
03248 struc.gas_phase[ipHYDROGEN][nzone_den_jump],
03249 struc.gas_phase[ipHYDROGEN][nzone_den_jump+1] );
03250 caunin(chLine);
03251 }
03252 else
03253 {
03254 sprintf( chLine,
03255 " !An H density jump occurred at zone %ld, by %.0f%% from %.2e to %.2e to %.2e",
03256 nzone_den_jump,
03257 big_jump*100.,
03258 struc.gas_phase[ipHYDROGEN][nzone_den_jump-1],
03259 struc.gas_phase[ipHYDROGEN][nzone_den_jump],
03260 struc.gas_phase[ipHYDROGEN][nzone_den_jump+1] );
03261 bangin(chLine);
03262 }
03263 }
03264
03265
03266 big_jump = -15;
03267 nzone_den_jump = 0;
03268
03269
03270
03271
03272 for( i=3; i < (nzone - 2 - (int)lgAbort); i++ )
03273 {
03274
03275 rel = fabs(log10( SDIV(struc.pres_radiation_lines_curr[i]) /
03276 SDIV(0.5*(struc.pres_radiation_lines_curr[i-1]+struc.pres_radiation_lines_curr[i+1])) ) );
03277
03278 if( rel > big_jump )
03279 {
03280 big_jump = rel;
03281 nzone_den_jump = i;
03282 }
03283 }
03284
03285
03286
03287 if( pressure.RadBetaMax > 0.01 )
03288 {
03289 big_jump = pow( 10., big_jump );
03290 if( big_jump > 1.2 )
03291 {
03292
03293
03294 if( big_jump > 3. && strcmp(dense.chDenseLaw,"CPRE") == 0)
03295 {
03296 sprintf( chLine,
03297 " C-The radiation pressure jumped by %.0f%% at zone %ld, from %.2e to %.2e to %.2e",
03298 big_jump*100.,
03299 nzone_den_jump,
03300 struc.pres_radiation_lines_curr[nzone_den_jump-1],
03301 struc.pres_radiation_lines_curr[nzone_den_jump],
03302 struc.pres_radiation_lines_curr[nzone_den_jump+1] );
03303 caunin(chLine);
03304 }
03305 else
03306 {
03307 sprintf( chLine,
03308 " !The radiation pressure jumped by %.0f%% at zone %ld, from %.2e to %.2e to %.2e",
03309 big_jump*100.,
03310 nzone_den_jump,
03311 struc.pres_radiation_lines_curr[nzone_den_jump-1],
03312 struc.pres_radiation_lines_curr[nzone_den_jump],
03313 struc.pres_radiation_lines_curr[nzone_den_jump+1] );
03314 bangin(chLine);
03315 }
03316 }
03317 }
03318
03319
03320 phycon.BigJumpTe = 0.;
03321 phycon.BigJumpne = 0.;
03322 phycon.BigJumpH2 = 0.;
03323 phycon.BigJumpCO = 0.;
03324
03325 for( i=1; i < (nzone - 1 - (int)lgAbort); i++ )
03326 {
03327
03328 rel = fabs(log10( struc.testr[i] / struc.testr[i+1] ) );
03329 if( rel > phycon.BigJumpTe )
03330 {
03331 phycon.BigJumpTe = (realnum)rel;
03332 }
03333
03334
03335 rel = fabs(log10( struc.ednstr[i] / struc.ednstr[i+1] ) );
03336 if( rel > phycon.BigJumpne )
03337 {
03338 phycon.BigJumpne = (realnum)rel;
03339 }
03340
03341
03342 if( (struc.H2_molec[ipMH2g][i]+struc.H2_molec[ipMH2s][i])>SMALLFLOAT &&
03343 (struc.H2_molec[ipMH2g][i+1]+struc.H2_molec[ipMH2s][i+1]) > SMALLFLOAT
03344
03345 && (struc.H2_molec[ipMH2g][i]+struc.H2_molec[ipMH2s][i])/struc.gas_phase[ipHYDROGEN][i]>1e-3)
03346 {
03347 rel = fabs(log10( (struc.H2_molec[ipMH2g][i]+struc.H2_molec[ipMH2s][i]) /
03348 SDIV(struc.H2_molec[ipMH2g][i+1]+struc.H2_molec[ipMH2s][i+1]) ) );
03349 if( rel > phycon.BigJumpH2 )
03350 {
03351 phycon.BigJumpH2 = (realnum)rel;
03352 }
03353 }
03354
03355 {
03356 int ipCO = findspecies("CO")->index;
03357
03358 if( struc.CO_molec[ipCO][i]>SMALLFLOAT &&
03359 struc.CO_molec[ipCO][i+1]>SMALLFLOAT &&
03360 struc.CO_molec[ipCO][i]/SDIV(struc.gas_phase[ipCARBON][i])>1e-3 )
03361 {
03362 rel = fabs(log10( struc.CO_molec[ipCO][i] / struc.CO_molec[ipCO][i+1] ) );
03363 if( rel > phycon.BigJumpCO )
03364 {
03365 phycon.BigJumpCO = (realnum)rel;
03366 }
03367 }
03368 }
03369 }
03370
03371
03372 if(phycon.BigJumpTe>0. )
03373 phycon.BigJumpTe = (realnum)pow( 10., (double)phycon.BigJumpTe ) -1.f;
03374
03375 if(phycon.BigJumpne>0. )
03376 phycon.BigJumpne = (realnum)pow( 10., (double)phycon.BigJumpne )-1.f;
03377
03378 if(phycon.BigJumpH2>0. )
03379 phycon.BigJumpH2 = (realnum)pow( 10., (double)phycon.BigJumpH2 )-1.f;
03380
03381 if(phycon.BigJumpCO>0. )
03382 phycon.BigJumpCO = (realnum)pow( 10., (double)phycon.BigJumpCO )-1.f;
03383
03384
03385
03386 if( phycon.BigJumpTe > 0.3 )
03387 {
03388 sprintf( chLine,
03389 " C-The temperature varied by %.1f%% between two zones",
03390 phycon.BigJumpTe*100.);
03391 caunin(chLine);
03392 }
03393 else if( phycon.BigJumpTe > 0.1 )
03394 {
03395 sprintf( chLine,
03396 " !The temperature varied by %.1f%% between two zones",
03397 phycon.BigJumpTe*100.);
03398 bangin(chLine);
03399 }
03400
03401 if( phycon.BigJumpne > 0.3 )
03402 {
03403 sprintf( chLine,
03404 " C-The electron density varied by %.1f%% between two zones",
03405 phycon.BigJumpne*100.);
03406 caunin(chLine);
03407 }
03408 else if( phycon.BigJumpne > 0.1 )
03409 {
03410 sprintf( chLine,
03411 " !The electron density varied by %.1f%% between two zones",
03412 phycon.BigJumpne*100.);
03413 bangin(chLine);
03414 }
03415
03416 if( phycon.BigJumpH2 > 0.8 )
03417 {
03418 sprintf( chLine,
03419 " C-The H2 density varied by %.1f%% between two zones",
03420 phycon.BigJumpH2*100.);
03421 caunin(chLine);
03422 }
03423 else if( phycon.BigJumpH2 > 0.1 )
03424 {
03425 sprintf( chLine,
03426 " !The H2 density varied by %.1f%% between two zones",
03427 phycon.BigJumpH2*100.);
03428 bangin(chLine);
03429 }
03430
03431 if( phycon.BigJumpCO > 0.8 )
03432 {
03433 sprintf( chLine,
03434 " C-The CO density varied by %.1f%% between two zones",
03435 phycon.BigJumpCO*100.);
03436 caunin(chLine);
03437 }
03438 else if( phycon.BigJumpCO > 0.2 )
03439 {
03440 sprintf( chLine,
03441 " !The CO density varied by %.1f%% between two zones",
03442 phycon.BigJumpCO*100.);
03443 bangin(chLine);
03444 }
03445 return;
03446 }