cloudy  trunk
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
prt_comment.cpp
Go to the documentation of this file.
1 /* This file is part of Cloudy and is copyright (C)1978-2022 by Gary J. Ferland and
2  * others. For conditions of distribution and use see copyright notice in license.txt */
3 /*PrtComment analyze model, generating comments on its features */
4 /*chkCaHeps check whether CaII K and H epsilon overlap */
5 /*prt_smooth_predictions check whether fluctuations in any predicted quantities occurred */
6 #include "cddefines.h"
7 #include "prt.h"
8 #include "cddrive.h"
9 #include "iso.h"
10 #include "continuum.h"
11 #include "stopcalc.h"
12 #include "hyperfine.h"
13 #include "grainvar.h"
14 #include "version.h"
15 #include "rt.h"
16 #include "he.h"
17 #include "taulines.h"
18 #include "hydrogenic.h"
19 #include "lines.h"
20 #include "trace.h"
21 #include "hcmap.h"
22 #include "hmi.h"
23 #include "save.h"
24 #include "h2.h"
25 #include "conv.h"
26 #include "dynamics.h"
27 #include "opacity.h"
28 #include "geometry.h"
29 #include "elementnames.h"
30 #include "ca.h"
31 #include "pressure.h"
32 #include "co.h"
33 #include "atoms.h"
34 #include "abund.h"
35 #include "colden.h"
36 #include "phycon.h"
37 #include "timesc.h"
38 #include "hextra.h"
39 #include "radius.h"
40 #include "iterations.h"
41 #include "fudgec.h"
42 #include "called.h"
43 #include "magnetic.h"
44 #include "wind.h"
45 #include "secondaries.h"
46 #include "struc.h"
47 #include "oxy.h"
48 #include "input.h"
49 #include "thermal.h"
50 #include "atmdat.h"
51 #include "warnings.h"
52 #include "mole.h"
53 #include "rfield.h"
54 #include "doppvel.h"
55 #include "freebound.h"
56 #include "two_photon.h"
57 #include "dense.h"
58 #include "lines_service.h"
59 
60 /*chkCaHeps check whether CaII K and H epsilon overlap */
61 STATIC void chkCaHeps(double *totwid);
62 
63 /*prt_smooth_predictions check whether fluctuations in any predicted quantities occurred */
64 STATIC void prt_smooth_predictions(void);
65 
66 void PrtComment(void)
67 {
68  char
69  chLine[2*INPUT_LINE_LENGTH];
70 
71  string chLbl;
72 
73  bool lgAbort_flag,
74  lgThick;
75 
76  long int dum1,
77  dum2,
78  i,
79  imas,
80  ipLo ,
81  ipHi ,
82  ipISO,
83  nelem,
84  isav,
85  j,
86  nc,
87  nline,
88  nn,
89  nneg,
90  ns,
91  nw;
92 
93  double big_ion_jump,
94  absint ,
95  aj,
96  alpha,
97  big,
98  bigm,
99  comfrc,
100  differ,
101  error,
102  flur,
103  freqn,
104  rate,
105  ratio,
106  rel,
107  small,
108  tauneg,
109  ts ,
110  HBeta,
111  relfl ,
112  relhm,
113  fedest,
114  GetHeat,
115  SumNeg,
116  c4363,
117  t4363,
118  totwid;
119 
120  double VolComputed , VolExpected , ConComputed , ConExpected;
121 
122  DEBUG_ENTRY( "PrtComment()" );
123 
124  if( 0 && lgAbort )
125  {
126  return;
127  }
128 
129  if( trace.lgTrace )
130  {
131  fprintf( ioQQQ, " PrtComment called.\n" );
132  }
133 
134  /*
135  * enter all comments cautions warnings and surprises into a large
136  * stack of statements
137  * at end of this routine is master call to printing routines
138  */
139  iterations.lgIterAgain = false;
140 
141  /* initialize the line saver */
142  warnings.zero();
143 
144  if( t_version::Inst().nBetaVer > 0 )
145  {
146  sprintf( chLine,
147  " !This is beta test version %ld and is intended for testing only.",
148  t_version::Inst().nBetaVer );
149  bangin(chLine);
150  }
151 
152  vector<module*>& mods=module_list::Inst().m_l;
153  for (vector<module*>::iterator mi = mods.begin(); mi != mods.end(); ++mi)
154  {
155  (*mi)->comment(warnings);
156  }
157 
158  /* say why calculation stopped */
159  if( conv.lgBadStop )
160  {
161  /* this case stop probably was not intended */
162  sprintf( warnings.chRgcln[0], " W-Calculation stopped because %s Iteration%3ld of %ld",
163  StopCalc.chReasonStop, iteration, iterations.itermx + 1 );
164  sprintf( chLine, " W-Calculation stopped because %s",
165  StopCalc.chReasonStop );
166  warnin(chLine);
167  sprintf( chLine, " W-This was not intended." );
168  warnin(chLine);
169  }
170  else
171  {
172  /* for iterate to convergence, print reason why it was not converged on 3rd and higher iterations */
173  if( (conv.lgAutoIt && iteration != iterations.itermx + 1) &&
174  iteration > 2 )
175  {
176  sprintf( warnings.chRgcln[0],
177  " Calculation stopped because %s Iteration %ld of %ld, not converged due to %s",
178  StopCalc.chReasonStop,
179  iteration,
180  iterations.itermx + 1,
181  conv.chNotConverged );
182  }
183  else
184  {
185  sprintf( warnings.chRgcln[0],
186  " Calculation stopped because %s Iteration %ld of %ld",
187  StopCalc.chReasonStop, iteration, iterations.itermx + 1 );
188  }
189  }
190 
191  /* check whether stopped because default number of zones hit,
192  * not intended?? */
193  if( (!geometry.lgZoneSet) && geometry.lgZoneTrp )
194  {
195  conv.lgBadStop = true;
196  sprintf( chLine,
197  " W-Calculation stopped because default number of zones reached. Was this intended???" );
198  warnin(chLine);
199  sprintf( chLine,
200  " W-Default limit can be increased while retaining this check with the SET NEND command." );
201  warnin(chLine);
202  }
203 
204  /* check whether stopped because zones too thin - only happens when glob set
205  * and depth + dr = depth
206  * not intended */
207  if( radius.lgDrMinUsed || radius.lgdR2Small )
208  {
209  conv.lgBadStop = true;
210  sprintf( chLine,
211  " W-Calculation stopped zone thickness became too thin. This was not intended." );
212  warnin(chLine);
213  sprintf( chLine,
214  " W-The most likely reason was an uncontrolled oscillation." );
215  warnin(chLine);
216  ShowMe();
217  }
218 
219  if( radius.lgdR2Small )
220  {
221  sprintf( chLine,
222  " W-This happened because the globule scale became very small relative to the depth." );
223  warnin(chLine);
224  sprintf( chLine,
225  " W-This problem is described in Hazy." );
226  warnin(chLine);
227  }
228 
229  /* possible that last zone does not have stored temp - if so
230  * then make it up - this happens for some bad stops */
231  ASSERT( nzone < struc.nzlim );
232 
233  if( struc.testr[nzone-1] == 0. && nzone > 1)
234  struc.testr[nzone-1] = struc.testr[nzone-2];
235 
236  if( struc.ednstr[nzone-1] == 0. && nzone > 1)
237  struc.ednstr[nzone-1] = struc.ednstr[nzone-2];
238 
239  /* give indication of geometry */
240  rel = radius.depth/radius.rinner;
241  if( rel < 0.1 )
242  {
243  sprintf( warnings.chRgcln[1], " The geometry is plane-parallel." );
244  }
245  else if( rel >= 0.1 && rel < 3. )
246  {
247  sprintf( warnings.chRgcln[1], " The geometry is a thick shell." );
248  }
249  else
250  {
251  sprintf( warnings.chRgcln[1], " The geometry is spherical." );
252  }
253  /* levels of warnings: Warning (possibly major problems)
254  * Caution (not likely to invalidate the results)
255  * [ !] surprise, but not a problem
256  * [nothing] interesting note
257  */
258 
259  /* incorrect electron density detected */
260  if( dense.lgEdenBad )
261  {
262  if( dense.nzEdenBad == nzone )
263  {
264  sprintf( chLine, " C-The assumed electron density was incorrect for the last zone." );
265  caunin(chLine);
266  sprintf( chLine, " C-Did a temperature discontinuity occur??" );
267  caunin(chLine);
268  }
269  else
270  {
271  sprintf( chLine, " W-The assumed electron density was incorrect during the calculation. This is bad." );
272  warnin(chLine);
273  ShowMe();
274  }
275  }
276 
277  if( lgAbort )
278  {
279  sprintf( chLine, " W-The calculation aborted. Something REALLY went wrong!" );
280  warnin(chLine);
281  }
282 
283  /* thermal map was done but results were not ok */
284  if( hcmap.lgMapDone && !hcmap.lgMapOK )
285  {
286  sprintf( chLine, " !The thermal map had changes in slope - check map output." );
287  bangin(chLine);
288  }
289 
290  /* first is greater than zero if fudge factors were entered, second is
291  * true if fudge ever evaluated, even to see if fudge factors are in place */
292  if( fudgec.nfudge > 0 || fudgec.lgFudgeUsed )
293  {
294  sprintf( chLine, " !Fudge factors were used or were checked. Why?" );
295  bangin(chLine);
296  }
297 
298  if( cpu.i().foundMD5Mismatch() )
299  {
300  sprintf( chLine, " !Modified data files were used in this simulation."
301  " This is fine if it was done intentionally." );
302  bangin(chLine);
303  }
304 
305  if( dense.gas_phase[ipHYDROGEN] > 1.1e13 )
306  {
307  if( dense.gas_phase[ipHYDROGEN] > 1e15 )
308  {
309  sprintf( chLine, " C-Density greater than 10**15, heavy elements are very uncertain." );
310  caunin(chLine);
311  }
312  else
313  {
314  sprintf( chLine, " C-Density greater than 10**13" );
315  caunin(chLine);
316  }
317  }
318 
319  /* HBeta is used later in the code to check on line intensities */
320  if( cdLine("Pump",4861.33f,&relfl,&absint)<=0 )
321  {
322  fprintf( ioQQQ, " PROBLEM Did not find Pump H-beta, set to unity\n" );
323  relfl = 1.;
324  absint = 1.;
325  }
326 
327  /* now find total Hbeta */
328  /* >>chng from "totl" Hbeta which was a special entry, to "H 1" Hbeta, which
329  * is the general case */
330  if( cdLine( "H 1",wlAirVac(4861.33),&HBeta,&absint)<=0 )
331  {
332  fprintf( ioQQQ, " NOTE Did not find H 1 H-beta - set intensity to unity, "
333  "will not check on importance of H 1 pumping.\n" );
334  HBeta = 1.;
335  absint = 1.;
336  }
337  else
338  {
339  /* check on continuum pumping of Balmer lines */
340  if( HBeta>SMALLFLOAT )
341  {
342  flur = relfl/HBeta;
343  if( flur > 0.1 )
344  {
345  sprintf( chLine, " !Continuum fluorescent production of H-beta was very important." );
346  bangin(chLine);
347  }
348  else if(flur > 0.01 )
349  {
350  sprintf( chLine, " Continuum fluorescent production of H-beta was significant." );
351  notein(chLine);
352  }
353  }
354  }
355 
356  // iterate to convergence - status of this iteration
357  if( iteration > 1 && conv.lgAutoIt && iteration<iterations.itermx)
358  {
359  sprintf( chLine , " Iteration not converged because %s.",
360  conv.chNotConverged );
361  notein(chLine);
362  }
363 
364  /* advice about "iterate to convergence" no converging at end
365  * lgIterAgain -> not converged, not another iteration */
366  if( conv.lgAutoIt && (iteration == iterations.itermx + 1) && !iterations.lgOpticalDepthonverged )
367  {
368  sprintf( chLine, " C-Iterate to convergence did not converge in %li iterations.",
369  iteration );
370  caunin(chLine);
371 
372  if( conv.lgAutoIt && (phycon.te < StopCalc.TempLoStopZone) )
373  {
374  sprintf( chLine, " C-Iterate to convergence requested but sim stopped due to too-low temperature.");
375  caunin(chLine);
376  sprintf( chLine, " C-This may have convergence problems due to outer radius changing.");
377  caunin(chLine);
378  sprintf( chLine, " C-Consider setting a different stop criterion.");
379  caunin(chLine);
380  }
381  }
382 
383  /* check if there were problems with initial wind velocity */
384  if( wind.lgBallistic() && ((!wind.lgWindOK) || wind.windv < 1e6) )
385  {
386  sprintf( chLine, " C-Wind velocity below sonic point; solution is not valid." );
387  caunin(chLine);
388  }
389 
390  /* now confirm that mass flux here is correct */
391  if( !wind.lgStatic() )
392  {
393  rel = wind.emdot/(wind.windv*dense.gas_phase[ipHYDROGEN])/radius.r1r0sq;
394  if( fabs(1.-rel)> 0.02 )
395  {
396  sprintf( chLine, " C-Wind mass flux error is %g%%",fabs(1.-rel)*100. );
397  caunin(chLine);
398  fprintf(ioQQQ,"DEBUG emdot\t%.3e\t%.3e\t%.3e\t%.3e\n",
399  wind.emdot , wind.windv*dense.gas_phase[ipHYDROGEN],wind.windv,dense.gas_phase[ipHYDROGEN]);
400  }
401  }
402 
403  /* check that we didn't overrun zone scale */
404  if( nzone >= struc.nzlim )
405  {
406  TotalInsanity();
407  }
408 
409  // check on energy conservation
410  if( !lgConserveEnergy() )
411  {
412  ShowMe();
413  lgAbort = true;
414  fprintf(ioQQQ,"\n\n Enable per zone energy conservation check by setting "
415  "CHECK_ENERGY_EVERY_ZONE=true in cloudy.cpp, recompile, then rerun.\n");
416  }
417 
418  /* comments having to do with cosmic rays */
419  /* comment if cosmic rays and magnetic field both present */
420  if( hextra.cryden*magnetic.lgB > 0. )
421  {
422  sprintf( chLine,
423  " !Magnetic field & cosmic rays both present. Their interactions are not treated." );
424  bangin(chLine);
425  }
426 
427  /* comment if cosmic rays are not included and stop temp has been lowered to go into neutral gas */
428  if( hextra.cryden== 0. && StopCalc.TempLoStopZone < phycon.TEMP_STOP_DEFAULT)
429  {
430  sprintf( chLine,
431  " !Background cosmic rays are not included - is this physical? It affects the chemistry." );
432  bangin(chLine);
433  }
434 
435  /* check whether cosmic rays on, but model thick to them */
436  if( hextra.cryden > 0. && (findspecieslocal("H")->column/10. + findspecieslocal("H+")->column) > 1e23 )
437  {
438  sprintf( chLine,
439  " C-Model is thick to cosmic rays, which are on." );
440  caunin(chLine);
441  }
442 
443  /* was ionization rate less than cosmic ray ionization rate in ISM? */
444  if( hextra.cryden == 0. && iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH1s].gamnc < 1e-17 )
445  {
446  sprintf( chLine,
447  " !Ionization rate fell below background cosmic ray ionization rate. Should this be added too?" );
448  bangin(chLine);
449  sprintf( chLine,
450  " ! Use the COSMIC RAY BACKGROUND command." );
451  bangin(chLine);
452  }
453 
454  /* PrtComment if test code is in place */
455  if( lgTestCodeCalled && !t_version::Inst().lgReleaseBranch && !t_version::Inst().lgRelease )
456  {
457  sprintf( chLine, " !Test code is in place." );
458  bangin(chLine);
459  }
460 
461  /* lgComUndr set to .true. if Compton cooling rate underflows to 0 */
462  if( rfield.lgComUndr )
463  {
464  sprintf( chLine,
465  " !Compton cooling rate underflows to zero. Is this important?" );
466  bangin(chLine);
467  }
468 
469  /* make note if input stream contained an underscore, which was converted into a space */
470  if( input.lgUnderscoreFound )
471  {
472  sprintf( chLine,
473  " !Some input lines contained underscores, these were changed to spaces." );
474  bangin(chLine);
475  }
476 
477  /* make note if input stream contained a left or right bracket, which was converted into a space */
478  if( input.lgBracketFound )
479  {
480  sprintf( chLine,
481  " !Some input lines contained [ or ], these were changed to spaces." );
482  bangin(chLine);
483  }
484 
485  /* lgHionRad set to .true. if no hydrogen ionizing radiation */
486  if( rfield.lgHionRad )
487  {
488  sprintf( chLine,
489  " !There is no hydrogen-ionizing radiation. Was this intended?" );
490  bangin(chLine);
491  }
492 
493  /* check whether certain zones were thermally unstable */
494  if( thermal.nUnstable > 0 )
495  {
496  sprintf( chLine,
497  " Derivative of net cooling negative and so possibly thermally unstable in%4ld zones.",
498  thermal.nUnstable );
499  notein(chLine);
500  }
501 
502  /* generate a bang if a large fraction of the zones were unstable */
503  if( nzone > 1 &&
504  (realnum)(thermal.nUnstable)/(realnum)(nzone) > 0.25 )
505  {
506  sprintf( chLine,
507  " !A large fraction of the zones were possibly thermally unstable,%4ld out of%4ld",
508  thermal.nUnstable, nzone );
509  bangin(chLine);
510  }
511 
512  /* comment if negative coolants were ever significant */
513  if( thermal.CoolHeatMax > 0.2 )
514  {
515  sprintf( chLine,
516  " !Negative cooling reached %6.1f%% of the local heating, due to %4.4s %.1f",
517  thermal.CoolHeatMax*100., thermal.chCoolHeatMax, thermal.wlCoolHeatMax );
518  bangin(chLine);
519  }
520  else if( thermal.CoolHeatMax > 0.05 )
521  {
522  sprintf( chLine,
523  " Negative cooling reached %6.1f%% of the local heating, due to %4.4s %.2f",
524  thermal.CoolHeatMax*100., thermal.chCoolHeatMax, thermal.wlCoolHeatMax );
525  notein(chLine);
526  }
527 
528  /* check if advection heating was important */
529  if( dynamics.HeatMax > 0.05 )
530  {
531  sprintf( chLine,
532  " !Advection heating reached %.2f%% of the local heating.",
533  dynamics.HeatMax*100. );
534  bangin(chLine);
535  }
536  else if( dynamics.HeatMax > 0.005 )
537  {
538  sprintf( chLine,
539  " Advection heating reached %.2f%% of the local heating.",
540  dynamics.HeatMax*100. );
541  notein(chLine);
542  }
543 
544  /* check if advection cooling was important */
545  if( dynamics.CoolMax > 0.05 )
546  {
547  sprintf( chLine,
548  " !Advection cooling reached %.2f%% of the local cooling.",
549  dynamics.CoolMax*100. );
550  bangin(chLine);
551  }
552  else if( dynamics.CoolMax > 0.005 )
553  {
554  sprintf( chLine,
555  " Advection cooling reached %.2f%% of the local heating.",
556  dynamics.CoolMax*100. );
557  notein(chLine);
558  }
559 
560  /* >>chng 06 mar 22, add this comment
561  * check if time dependent ionization front being done with too large a U */
562  if( dynamics.lgTimeDependentStatic && dynamics.lgRecom )
563  {
564  if( rfield.uh > 1. )
565  {
566  sprintf( chLine,
567  " W-Time dependent ionization front cannot now handle strong-R cases - the ionization parameter is too large." );
568  warnin(chLine);
569  }
570  else if( rfield.uh > 0.1 )
571  {
572  sprintf( chLine,
573  " C-Time dependent ionization front cannot now handle strong-R cases - the ionization parameter is too large." );
574  caunin(chLine);
575  }
576  }
577 
578  /* check if thermal ionization of ground state of hydrogen was important */
579  if( hydro.HCollIonMax > 0.10 )
580  {
581  sprintf( chLine,
582  " !Thermal collisional ionization of H reached %.2f%% of the local ionization rate.",
583  hydro.HCollIonMax*100. );
584  bangin(chLine);
585  }
586  else if( hydro.HCollIonMax > 0.005 )
587  {
588  sprintf( chLine,
589  " Thermal collisional ionization of H reached %.2f%% of the local ionization rate.",
590  hydro.HCollIonMax*100. );
591  notein(chLine);
592  }
593 
594  /* check if lookup table for Hummer & Storey case B was exceeded */
595  if( !atmdat.lgHCaseBOK[1][ipHYDROGEN] )
596  {
597  sprintf( chLine,
598  " Te-ne bounds of Case B lookup table exceeded, H I Case B line intensities set to zero." );
599  notein(chLine);
600  }
601  if( !atmdat.lgHCaseBOK[1][ipHELIUM] )
602  {
603  sprintf( chLine,
604  " Te-ne bounds of Case B lookup table exceeded, He II Case B line intensities set to zero." );
605  notein(chLine);
606  }
607 
608  if( dense.EdenMax>1e8 )
609  {
610  sprintf( chLine,
611  " !The high electron density makes the Nussbaumer/Storey CNO recombination predictions unreliable." );
612  bangin(chLine);
613  }
614 
615  /* check if secondary ionization of hydrogen was important */
616  if( secondaries.SecHIonMax > 0.10 )
617  {
618  sprintf( chLine,
619  " !Suprathermal collisional ionization of H reached %.2f%% of the local H ionization rate.",
620  secondaries.SecHIonMax*100. );
621  bangin(chLine);
622  }
623  else if( secondaries.SecHIonMax > 0.005 )
624  {
625  sprintf( chLine,
626  " Suprathermal collisional ionization of H reached %.2f%% of the local H ionization rate.",
627  secondaries.SecHIonMax*100. );
628  notein(chLine);
629  }
630 
631  /* check if H2 vib-deexcitation heating was important */
632  if( hmi.HeatH2DexcMax > 0.05 )
633  {
634  sprintf( chLine,
635  " !H2 vib deexec heating reached %.2f%% of the local heating.",
636  hmi.HeatH2DexcMax*100. );
637  bangin(chLine);
638  }
639  else if( hmi.HeatH2DexcMax > 0.005 )
640  {
641  sprintf( chLine,
642  " H2 vib deexec heating reached %.2f%% of the local heating.",
643  hmi.HeatH2DexcMax*100. );
644  notein(chLine);
645  }
646 
647  /* check if H2 vib-deexcitation heating was important */
648  if( hmi.CoolH2DexcMax > 0.05 )
649  {
650  sprintf( chLine,
651  " !H2 deexec cooling reached %.2f%% of the local heating.",
652  hmi.CoolH2DexcMax*100. );
653  bangin(chLine);
654  }
655  else if( hmi.CoolH2DexcMax > 0.005 )
656  {
657  sprintf( chLine,
658  " H2 deexec cooling reached %.2f%% of the local heating.",
659  hmi.CoolH2DexcMax*100. );
660  notein(chLine);
661  }
662 
663  /* check if charge transfer ionization of hydrogen was important */
664  if( atmdat.HIonFracMax > 0.10 )
665  {
666  sprintf( chLine,
667  " !Charge transfer H => H+ reached %.1f%% of the local H ionization rate.",
668  atmdat.HIonFracMax*100. );
669  bangin(chLine);
670  }
671  else if( atmdat.HIonFracMax > 0.005 )
672  {
673  sprintf( chLine,
674  " Charge transfer H => H+ reached %.2f%% of the local H ionization rate.",
675  atmdat.HIonFracMax*100. );
676  notein(chLine);
677  }
678 
679  /* check if charge transfer heating cooling was important */
680  if( atmdat.HCharHeatMax > 0.05 )
681  {
682  sprintf( chLine,
683  " !Charge transfer heating reached %.2f%% of the local heating.",
684  atmdat.HCharHeatMax*100. );
685  bangin(chLine);
686  }
687  else if( atmdat.HCharHeatMax > 0.005 )
688  {
689  sprintf( chLine,
690  " Charge transfer heating reached %.2f%% of the local heating.",
691  atmdat.HCharHeatMax*100. );
692  notein(chLine);
693  }
694 
695  if( atmdat.HCharCoolMax > 0.05 )
696  {
697  sprintf( chLine,
698  " !Charge transfer cooling reached %.2f%% of the local heating.",
699  atmdat.HCharCoolMax*100. );
700  bangin(chLine);
701  }
702  else if( atmdat.HCharCoolMax > 0.005 )
703  {
704  sprintf( chLine,
705  " Charge transfer cooling reached %.2f%% of the local heating.",
706  atmdat.HCharCoolMax*100. );
707  notein(chLine);
708  }
709 
710  /* check whether photo from up level of Mg2 2798 ever important */
711  if( atoms.xMg2Max > 0.1 )
712  {
713  sprintf( chLine,
714  " !Photoionization of upper level of Mg II 2798 reached %.1f%% of the total Mg+ photo rate.",
715  atoms.xMg2Max*100. );
716  bangin(chLine);
717  }
718  else if( atoms.xMg2Max > 0.01 )
719  {
720  sprintf( chLine,
721  " Photoionization of upper level of Mg II 2798 reached %.1f%% of the total Mg+ photo rate.",
722  atoms.xMg2Max*100. );
723  notein(chLine);
724  }
725 
726  /* check whether photo from up level of [O I] 6300 ever important */
727  if( oxy.poimax > 0.1 )
728  {
729  sprintf( chLine,
730  " !Photoionization of upper levels of [O I] reached %.1f%% of the total O destruction rate.",
731  oxy.poimax*100. );
732  bangin(chLine);
733  }
734  else if( oxy.poimax > 0.01 )
735  {
736  sprintf( chLine,
737  " Photoionization of upper levels of [O I] reached %.1f%% of the total O destruction rate.",
738  oxy.poimax*100. );
739  notein(chLine);
740  }
741 
742  /* check whether photo from up level of [O III] 5007 ever important */
743  if( (oxy.poiii2Max + oxy.poiii3Max) > 0.1 )
744  {
745  sprintf( chLine,
746  " !Photoionization of upper levels of [O III] reached %.1f%% of the total O++ photo rate.",
747  (oxy.poiii2Max + oxy.poiii3Max)*100. );
748  bangin(chLine);
749  }
750  else if( (oxy.poiii2Max + oxy.poiii3Max) > 0.01 )
751  {
752  sprintf( chLine,
753  " Photoionization of upper levels of [O III] reached %.1f%% of the total O++ photo rate.",
754  (oxy.poiii2Max + oxy.poiii3Max)*100. );
755  notein(chLine);
756  }
757 
758  /* check whether photoionization of He 2trip S was important */
759  if( he.frac_he0dest_23S > 0.1 )
760  {
761  sprintf( chLine,
762  " !Destruction of He 2TriS reached %.1f%% of the total He0 dest rate"
763  " at zone %li, %.1f%% of that was photoionization.",
764  he.frac_he0dest_23S*100,
765  he.nzone,
766  he.frac_he0dest_23S_photo*100. );
767  bangin(chLine);
768  }
769  else if( he.frac_he0dest_23S > 0.01 )
770  {
771  sprintf( chLine,
772  " Destruction of He 2TriS reached %.1f%% of the total He0 dest rate"
773  " at zone %li, %.1f%% of that was photoionization.",
774  he.frac_he0dest_23S*100,
775  he.nzone,
776  he.frac_he0dest_23S_photo*100. );
777  notein(chLine);
778  }
779 
780  /* check for critical density for l-mixing */
781  for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
782  {
783  if( !iso_ctrl.lgCritDensLMix[ipISO] && dense.lgElmtOn[ipISO] )
784  {
785  sprintf( chLine,
786  " The density is too low to l-mix the lowest %s I collapsed level. "
787  " More resolved levels are needed for accurate line ratios.",
788  elementnames.chElementSym[ipISO]);
789  notein(chLine);
790  }
791  }
792 
793  for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
794  {
795  /* report continuum lowering for xx-like iso-sequence. */
796  for( nelem=ipISO; nelem<LIMELM; ++nelem )
797  {
798  if( iso_sp[ipISO][nelem].lgLevelsLowered && dense.lgElmtOn[nelem] )
799  {
800  sprintf( chLine, " !Continuum was lowered into model %s-like %s due to high density. Highest n is %li",
801  elementnames.chElementSym[ipISO],
802  elementnames.chElementSym[nelem],
803  iso_sp[ipISO][nelem].n_HighestResolved_local+iso_sp[ipISO][nelem].nCollapsed_local);
804  bangin(chLine);
805  }
806  else if( iso_sp[ipISO][nelem].lgLevelsEverLowered && dense.lgElmtOn[nelem] )
807  {
808  sprintf( chLine, " !Continuum was lowered into model %s-like %s due to high density at SOME point but NOT at the last zone.",
809  elementnames.chElementSym[ipISO],
810  elementnames.chElementNameShort[nelem]);
811  bangin(chLine);
812  }
813  }
814 
815  /* report pop rescaling for xx-like iso-sequence. */
816  for( nelem=ipISO; nelem<LIMELM; ++nelem )
817  {
818  if( iso_sp[ipISO][nelem].lgPopsRescaled )
819  {
820  ASSERT( dense.lgSetIoniz[nelem] );
821  sprintf( chLine, " C-Populations were rescaled for %s-like %s due to \"element ionization\" command.",
822  elementnames.chElementSym[ipISO],
823  elementnames.chElementSym[nelem] );
824  caunin(chLine);
825  }
826  }
827  }
828 
829  /* frequency array may not have been defined for all energies */
830  if( !rfield.lgMMok )
831  {
832  sprintf( chLine,
833  " C-Continuum not defined in extreme infrared - Compton scat, grain heating, not treated properly?" );
834  caunin(chLine);
835  }
836 
837  if( !rfield.lgHPhtOK )
838  {
839  sprintf( chLine,
840  " C-Continuum not defined at photon energies which ionize excited states of H, important for H- and ff heating." );
841  caunin(chLine);
842  }
843 
844  if( !rfield.lgXRayOK )
845  {
846  sprintf( chLine,
847  " C-Continuum not defined at X-Ray energies - Compton scattering and Auger ionization wrong?" );
848  caunin(chLine);
849  }
850 
851  if( !rfield.lgGamrOK )
852  {
853  sprintf( chLine,
854  " C-Continuum not defined at gamma-ray energies - pair production and Compton scattering OK?" );
855  caunin(chLine);
856  }
857 
858  if( continuum.lgCon0 )
859  {
860  sprintf( chLine, " C-Continuum zero at some energies." );
861  caunin(chLine);
862  }
863 
864  if( continuum.lgCoStarInterpolationCaution )
865  {
866  sprintf( chLine , " C-CoStarInterpolate interpolated between non-adjoining tracks, this may not be valid." );
867  caunin(chLine);
868  }
869 
870  if( rfield.lgOcc1Hi )
871  {
872  sprintf( chLine,
873  " !The continuum occupation number at 1 Ryd is greater than unity." );
874  bangin(chLine);
875  }
876 
877  /* this flag set true it set dr forced first zone to be too big */
878  if( radius.lgDR2Big )
879  {
880  sprintf( chLine,
881  " C-The thickness of the first zone was set larger than optimal by a SET DR command." );
882  caunin(chLine);
883  /* this is case where did one zone of specified thickness - but it
884  * was too large */
885  if( nzone<=1 )
886  sprintf( chLine,
887  " C-Consider using the STOP THICKNESS command instead." );
888  caunin(chLine);
889  }
890 
891  /* check whether non-col excitation of 4363 was important */
892  if( atmdat.lgdBaseSourceExists[ipOXYGEN][2] )
893  {
894 
895  if( cdLine("Blnd",4363,&t4363,&absint)<=0 )
896  {
897  fprintf( ioQQQ, " PrtComment could not find total O III 4363 with cdLine.\n" );
898  ShowMe();
899  cdEXIT(EXIT_FAILURE);
900  }
901 
902  if( cdLine("O 3",wlAirVac(4363.21),&c4363,&absint)<=0 )
903  {
904  fprintf( ioQQQ, " PrtComment could not find collisional O III 4363 with cdLine, assuming intensity of zero.\n" );
905  c4363 = 0;
906  absint = -37.;
907  }
908 
909  /* only print this comment if 4363 is significant and density low */
910  if( HBeta > 1e-35 )
911  {
912  /* >>chng 02 feb 27, lower ratio from -4 to -5, and raise density from 7 to 8 */
913  if( t4363/HBeta > 1e-5 && dense.gas_phase[ipHYDROGEN] < 1e8 )
914  {
915  ratio = (t4363 - c4363)/t4363;
916  if( ratio > 0.01 )
917  {
918  sprintf( chLine,
919  " !Non-collisional excitation of [O III] 4363 reached %.2f%% of the total.",
920  ratio*100. );
921  bangin(chLine);
922  }
923  else if( ratio > 0.001 )
924  {
925  sprintf( chLine,
926  " Non-collisional excitation of [O III] 4363 reached %.2f%% of the total.",
927  ratio*100. );
928  notein(chLine);
929  }
930  }
931  }
932  }
933 
934  /* check for plasma shielding */
935  if( rfield.lgPlasNu )
936  {
937  sprintf( chLine,
938  " !The largest plasma frequency was %.2e Ryd = %.2e micron The continuum is set to 0 below this.",
939  rfield.plsfrqmax,
940  /* wavelength in microns */
941  RYDLAM/rfield.plsfrqmax/1e4);
942  bangin(chLine);
943  }
944 
945  if( rfield.occmax > 0.1 )
946  {
947  if( rfield.occmnu > 1e-4 )
948  {
949  sprintf( chLine,
950  " !The largest continuum occupation number was %.3e at %.3e Ryd.",
951  rfield.occmax, rfield.occmnu );
952  bangin(chLine);
953  }
954  else
955  {
956  /* not surprising if occupation number bigger than 1 around 1e-5 Ryd,
957  * since this is the case for 3K background */
958  sprintf( chLine,
959  " The largest continuum occupation number was %.3e at %.3e Ryd.",
960  rfield.occmax, rfield.occmnu );
961  notein(chLine);
962  }
963  }
964 
965  if( rfield.occmax > 1e4 && rfield.occ1nu > 0. )
966  {
967  /* occ1nu is energy (ryd) where continuum occupation number falls below 1 */
968  if( rfield.occ1nu < 0.0912 )
969  {
970  sprintf( chLine,
971  " The continuum occupation number fell below 1 at %.3e microns.",
972  0.0912/rfield.occ1nu );
973  notein(chLine);
974  }
975  else if( rfield.occ1nu < 1. )
976  {
977  sprintf( chLine,
978  " The continuum occupation number fell below 1 at %.3e Angstroms.",
979  912./rfield.occ1nu );
980  notein(chLine);
981  }
982  else
983  {
984  sprintf( chLine,
985  " The continuum occupation number fell below 1 at %.3e Ryd.",
986  rfield.occ1nu );
987  notein(chLine);
988  }
989  }
990 
991  if( rfield.tbrmax > 1e3 )
992  {
993  sprintf( chLine,
994  " !The largest continuum brightness temperature was %.3eK at %.3e Ryd.",
995  rfield.tbrmax, rfield.tbrmnu );
996  bangin(chLine);
997  }
998 
999  if( rfield.tbrmax > 1e4 )
1000  {
1001  /* tbr4nu is energy (ryd) where continuum bright temp falls < 1e4 */
1002  if( rfield.tbr4nu < 0.0912 )
1003  {
1004  sprintf( chLine,
1005  " The continuum brightness temperature fell below 10000K at %.3e microns.",
1006  0.0912/rfield.tbr4nu );
1007  notein(chLine);
1008  }
1009  else if( rfield.tbr4nu < 1. )
1010  {
1011  sprintf( chLine,
1012  " The continuum brightness temperature fell below 10000K at %.3e Angstroms.",
1013  912./rfield.tbr4nu );
1014  notein(chLine);
1015  }
1016  else
1017  {
1018  sprintf( chLine,
1019  " The continuum brightness temperature fell below 10000K at %.3e Ryd.",
1020  rfield.tbr4nu );
1021  notein(chLine);
1022  }
1023  }
1024 
1025  /* turbulence AND constant pressure do not make sense */
1026  if( DoppVel.TurbVel > 0. && strcmp(dense.chDenseLaw,"CPRE") == 0 )
1027  {
1028  sprintf( chLine,
1029  " !Both constant pressure and turbulence makes no physical sense?" );
1030  bangin(chLine);
1031  }
1032 
1033  /* filling factor AND constant pressure do not make sense */
1034  if( geometry.FillFac < 1. && strcmp(dense.chDenseLaw,"CPRE") == 0 )
1035  {
1036  sprintf( chLine,
1037  " !Both constant pressure and a filling factor makes no physical sense?" );
1038  bangin(chLine);
1039  }
1040 
1041  /* grains and solar abundances do not make sense */
1042  if( gv.lgDustOn() && abund.lgAbnSolar )
1043  {
1044  sprintf( chLine,
1045  " !Grains are present, but the gas phase abundances were left at the solar default. This is not physical." );
1046  bangin(chLine);
1047  }
1048 
1049  /* check if depletion command set but no grains, another silly thing to do */
1050  if( abund.lgDepln && !gv.lgDustOn() )
1051  {
1052  sprintf( chLine,
1053  " !Grains are not present, but the gas phase abundances were depleted. This is not physical." );
1054  bangin(chLine);
1055  }
1056 
1057  if( gv.lgDustOn() )
1058  {
1059  long nBin=0L, nFail=0L;
1060  for( size_t nd=0; nd < gv.bin.size(); nd++ )
1061  {
1062  if( gv.bin[nd]->QHeatFailures > 0L )
1063  {
1064  ++nBin;
1065  nFail += gv.bin[nd]->QHeatFailures;
1066  }
1067  }
1068  if( nFail > 0 )
1069  {
1070  sprintf( chLine,
1071  " !The grain quantum heating treatment failed to converge %ld time(s) in %ld bin(s).", nFail, nBin );
1072  bangin(chLine);
1073  }
1074  }
1075 
1076 #if 0
1077  /* check if PAHs were present in the ionized region */
1078  /* >>chng 05 jan 01, disabled this code now that PAH's have varying abundances by default, PvH */
1080  if( gv.lgDustOn() )
1081  {
1082  bool lgPAHsPresent_and_constant = false;
1083  for( size_t nd=0; nd < gv.bin.size(); nd++ )
1084  {
1085  lgPAHsPresent_and_constant = lgPAHsPresent_and_constant ||
1086  /* it is ok to have PAHs in the ionized region if the abundances vary */
1087  (gv.bin[nd]->lgPAHsInIonizedRegion /* && !gv.bin[nd]-> lgDustVary */);
1088  }
1089  if( lgPAHsPresent_and_constant )
1090  {
1091  sprintf( chLine,
1092  " C-PAH's were present in the ionized region, this has never been observed in H II regions." );
1093  caunin(chLine);
1094  }
1095  }
1096 #endif
1097 
1098  /* constant temperature greater than continuum energy density temperature */
1099  if( thermal.lgTemperatureConstant && thermal.ConstTemp*1.0001 < phycon.TEnerDen )
1100  {
1101  sprintf( chLine,
1102  " C-The continuum energy density temperature (%g K)"
1103  " is greater than the gas kinetic temperature (%g K).",
1104  phycon.TEnerDen , thermal.ConstTemp);
1105  caunin(chLine);
1106  sprintf( chLine, " C-This is unphysical." );
1107  caunin(chLine);
1108  }
1109 
1110  /* remark that grains not present but energy density was low */
1111  if( !gv.lgDustOn() && phycon.TEnerDen < 800. )
1112  {
1113  sprintf( chLine,
1114  " Grains were not present but might survive in this environment (energy density temperature was %.2eK)",
1115  phycon.TEnerDen );
1116  notein(chLine);
1117  }
1118 
1119  /* call routine that will check age of cloud */
1120  AgeCheck();
1121 
1122  /* check on Ca H and H-epsilon overlapping
1123  * need to do this since need to refer to lines arrays */
1124  chkCaHeps(&totwid);
1125  if( totwid > 121. )
1126  {
1127  sprintf( chLine, " H-eps and Ca H overlap." );
1128  notein(chLine);
1129  }
1130 
1131  /* warning that something was turned off */
1132  if( !phycon.lgPhysOK )
1133  {
1134  sprintf( chLine, " !A physical process has been disabled." );
1135  bangin(chLine);
1136  }
1137 
1138  /* check on lifetimes of [O III] against photoionization, only for low den */
1139  if( dense.gas_phase[ipHYDROGEN] < 1e8 )
1140  {
1141  if( oxy.r5007Max > 0.0263f )
1142  {
1143  sprintf( chLine,
1144  " !Photoionization of upper level of [O III] 5007 reached %.2e%% of the radiative lifetime.",
1145  oxy.r5007Max*100. );
1146  bangin(chLine);
1147  }
1148  else if( oxy.r5007Max > 0.0263f/10.f )
1149  {
1150  sprintf( chLine,
1151  " Photoionization of upper level of [O III] 5007 reached %.2e%% of the radiative lifetime.",
1152  oxy.r5007Max*100. );
1153  notein(chLine);
1154  }
1155  if( oxy.r4363Max > 1.78f )
1156  {
1157  sprintf( chLine,
1158  " !Photoionization of upper level of [O III] 4363 reached %.2e%% of the radiative lifetime.",
1159  oxy.r4363Max*100. );
1160  bangin(chLine);
1161  }
1162  else if( oxy.r4363Max > 1.78f/10.f )
1163  {
1164  sprintf( chLine,
1165  " Photoionization of upper level of [O III] 4363 reached %.2e%% of the radiative lifetime.",
1166  oxy.r4363Max*100. );
1167  notein(chLine);
1168  }
1169  }
1170 
1171  /* check whether total heating and cooling matched
1172  * >>chng 97 mar 28, added GrossHeat, heat in terms normally heat-cool */
1173  error = fabs(thermal.power-thermal.totcol)/SDIV((thermal.power + thermal.totcol)/2.);
1174  if( thermal.lgTemperatureConstant )
1175  {
1176  if( error > 0.05 )
1177  {
1178  sprintf( chLine,
1179  " !Heating - cooling mismatch =%5.1f%%. Caused by constant temperature assumption. ",
1180  error*100. );
1181  bangin(chLine);
1182  }
1183  }
1184 
1185  else
1186  {
1187  if( error > 0.05 && error < 0.2 )
1188  {
1189  sprintf( chLine, " C-Heating - cooling mismatch =%.1f%%. What\'s wrong?",
1190  error*100. );
1191  caunin(chLine);
1192  }
1193  else if( error >= 0.2 )
1194  {
1195  sprintf( chLine, " W-Heating - cooling mismatch =%.2e%%. What\'s wrong????",
1196  error*100. );
1197  warnin(chLine);
1198  }
1199  }
1200 
1201  /* say if Ly-alpha photo of Ca+ excited levels was important */
1202  if( ca.Ca2RmLya > 0.01 )
1203  {
1204  sprintf( chLine,
1205  " Photoionization of Ca+ 2D level by Ly-alpha reached %6.1f%% of the total rate out.",
1206  ca.Ca2RmLya*100. );
1207  notein(chLine);
1208  }
1209 
1210  /* check if Lya alpha ever hotter than gas */
1211  if( hydro.nLyaHot > 0 )
1212  {
1213  if( hydro.TLyaMax/hydro.TeLyaMax > 1.05 )
1214  {
1215  sprintf( chLine,
1216  " !The excitation temp of Lya exceeded the electron temp, largest value was %.2eK (gas temp there was %.2eK, zone %ld)",
1217  hydro.TLyaMax, hydro.TeLyaMax, hydro.nZTLaMax );
1218  bangin(chLine);
1219  }
1220  }
1221 
1222  /* check if line absorption heating was important */
1223 
1224  /* get all negative lines, check if line absorption significant heat source
1225  * this is used in "final" for energy budget print out */
1226  if( cdLine("Line",0,&SumNeg,&absint)<=0 )
1227  {
1228  fprintf( ioQQQ, " did not get sumneg cdLine\n" );
1229  ShowMe();
1230  cdEXIT(EXIT_FAILURE);
1231  }
1232 
1233  /* this is total heating */
1234  if( cdLine("TotH",0,&GetHeat,&absint)<=0 )
1235  {
1236  fprintf( ioQQQ, " did not get GetHeat cdLine\n" );
1237  ShowMe();
1238  cdEXIT(EXIT_FAILURE);
1239  }
1240 
1241  if( GetHeat > 0. )
1242  {
1243  SumNeg /= GetHeat;
1244  if( SumNeg > 0.1 )
1245  {
1246  sprintf( chLine,
1247  " !Line absorption heating reached %.2f%% of the global heating.",
1248  SumNeg*100. );
1249  bangin(chLine);
1250  }
1251  else if( SumNeg > 0.01 )
1252  {
1253  sprintf( chLine,
1254  " Line absorption heating reached %.2f%% of the global heating.",
1255  SumNeg*100. );
1256  notein(chLine);
1257  }
1258  }
1259 
1260  if( input.lgSetNoBuffering )
1261  {
1262  sprintf( chLine,
1263  " !NO BUFFERING command was entered - this increases exec time by LARGE amounts.");
1264  bangin(chLine);
1265  }
1266 
1267  /* this is check of extra lines added with g-bar */
1268  if( thermal.GBarMax > 0.1 )
1269  {
1270  ASSERT( thermal.ipMaxExtra > 0 );
1271  chLbl = chLineLbl(TauLine2[thermal.ipMaxExtra-1]);
1272 
1273  sprintf( chLine,
1274  " !G-bar cooling lines reached %.2f%% of the local cooling. Line=%.10s",
1275  thermal.GBarMax*100., chLbl.c_str() );
1276  bangin(chLine);
1277  }
1278 
1279  else if( thermal.GBarMax > 0.01 )
1280  {
1281  chLbl = chLineLbl(TauLine2[thermal.ipMaxExtra-1]);
1282 
1283  sprintf( chLine,
1284  " G-bar cooling lines reached %.2f%% of the local cooling. Line=%.10s",
1285  thermal.GBarMax*100., chLbl.c_str() );
1286  notein(chLine);
1287  }
1288 
1289  /* this is check of hyperfine structure lines*/
1290  if( hyperfine.cooling_max > 0.1 )
1291  {
1292  sprintf( chLine,
1293  " !Hyperfine structure line cooling reached %.2f%% of the local cooling.",
1294  hyperfine.cooling_max*100.);
1295  bangin(chLine);
1296  }
1297 
1298  else if( hyperfine.cooling_max > 0.01 )
1299  {
1300  sprintf( chLine,
1301  " Hyperfine structure line cooling reached %.2f%% of the local cooling.",
1302  hyperfine.cooling_max*100. );
1303  notein(chLine);
1304  }
1305 
1306  /* line absorption heating reached more than 10% of local heating?
1307  * HeatLineMax is largest heating(1,23)/htot */
1308  if( thermal.HeatLineMax > 0.1 )
1309  {
1310  long level = -1;
1311  TransitionProxy t = FndLineHt(&level);
1312  chLbl = chLineLbl(t);
1313  sprintf( chLine,
1314  " !Line absorption heating reached %.2f%% of the local heating - largest by level%2ld line %.10s",
1315  thermal.HeatLineMax*100., level, chLbl.c_str() );
1316  bangin(chLine);
1317  }
1318 
1319  else if( thermal.HeatLineMax > 0.01 )
1320  {
1321  sprintf( chLine,
1322  " Line absorption heating reached %.2f%% of the local heating.",
1323  thermal.HeatLineMax*100. );
1324  notein(chLine);
1325  }
1326 
1327  /* check whether any lines in the iso sequences mased */
1328  for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
1329  {
1330  for( nelem=ipISO; nelem<LIMELM; ++nelem )
1331  {
1332  if( dense.lgElmtOn[nelem] )
1333  {
1334  /* >>chng 06 aug 17, should go to numLevels_local instead of _max. */
1335  long int nmax = iso_sp[ipISO][nelem].numLevels_local;
1336 
1337  /* minus one here is to exclude highest level */
1338  for( ipHi=1; ipHi < nmax - 1; ++ipHi )
1339  {
1340  for( ipLo=0; ipLo < ipHi; ++ipLo )
1341  {
1342  if( iso_sp[ipISO][nelem].trans(ipHi,ipLo).Emis().Aul() <= iso_ctrl.SmallA )
1343  continue;
1344 
1345  /* did the line mase */
1346  if( iso_sp[ipISO][nelem].trans(ipHi,ipLo).Emis().TauIn() < -0.1 )
1347  {
1348  sprintf( chLine,
1349  " !Some iso-structure lines mased: %s-like %s, line %li-%li had optical depth %.2e",
1350  elementnames.chElementSym[ipISO],
1351  elementnames.chElementNameShort[nelem],
1352  ipHi , ipLo ,
1353  iso_sp[ipISO][nelem].trans(ipHi,ipLo).Emis().TauIn() );
1354  bangin(chLine);
1355  }
1356  }
1357  }
1358  }
1359  }
1360  }
1361 
1362  if( dense.gas_phase[ipHYDROGEN] < 1e7 )
1363  {
1364  /* check on IR fine structure lines - not necessary if dense since will be in LTE */
1365  lgThick = false;
1366  tauneg = 0.;
1367  alpha = 0.;
1368 
1369  /* now print results, were any fine structure lines optically thick? */
1370  if( lgThick )
1371  {
1372  sprintf( chLine,
1373  " !Some infrared fine structure lines are optically thick: largest tau was %.2e",
1374  alpha );
1375  bangin(chLine);
1376  }
1377  /* did any fine structure lines mase? */
1378  if( tauneg < -0.01 )
1379  {
1380  sprintf( chLine,
1381  " !Some fine structure lines mased: line %s had optical depth %.2e",
1382  chLbl.c_str(), tauneg );
1383  bangin(chLine);
1384  }
1385  }
1386 
1387  /* were any other lines masing? */
1388  /* this is check that at least a second iteration was done with sphere static,
1389  * the test is overridden with the (OK) option on the sphere static command,
1390  * which sets geometry.lgStaticNoIt true */
1391  if( geometry.lgStatic && iterations.lgLastIt && (iteration == 1) &&
1392  !geometry.lgStaticNoIt)
1393  {
1394  sprintf( chLine, " C-I must iterate when SPHERE STATIC is set." );
1395  caunin(chLine);
1396  iterations.lgIterAgain = true;
1397  }
1398 
1399  /* caution if continuum is punched but only one iteration performed */
1400  if( save.lgPunContinuum && iteration == 1 && iterations.lgLastIt)
1401  {
1402  sprintf( chLine, " C-I must iterate when save continuum output is done." );
1403  caunin(chLine);
1404  iterations.lgIterAgain = true;
1405  }
1406 
1408  /* how important was induced two photon?? */
1409  {
1410  two_photon& tnu = iso_sp[ipH_LIKE][ipHYDROGEN].TwoNu[0];
1411  if( tnu.induc_dn_max > 1. )
1412  {
1413  sprintf( chLine, " !Rate of induced H 2-photon emission reached %.2e s^-1",
1414  tnu.induc_dn_max );
1415  bangin(chLine);
1416  }
1417  else if( tnu.induc_dn_max > 0.01 )
1418  {
1419  sprintf( chLine, " Rate of induced H 2-photon emission reached %.2e s^-1",
1420  tnu.induc_dn_max );
1421  notein(chLine);
1422  }
1423  }
1424 
1425  /* how important was induced recombination? */
1426  if( hydro.FracInd > 0.01 )
1427  {
1428  sprintf( chLine,
1429  " Induced recombination was %5.1f%% of the total for H level%3ld",
1430  hydro.FracInd*100., hydro.ndclev );
1431  notein(chLine);
1432  }
1433 
1434  if( hydro.fbul > 0.01 )
1435  {
1436  sprintf( chLine,
1437  " Stimulated emission was%6.1f%% of the total for H transition%3ld -%3ld",
1438  hydro.fbul*100., hydro.nbul + 1, hydro.nbul );
1439  notein(chLine);
1440  }
1441 
1442  /* check whether Fe II destruction of La was important - entry into lines stack
1443  * is in prt_lines_hydro.c */
1444  if( cdLine("Fe 2",1215.67,&fedest,&absint)<=0 )
1445  {
1446  fprintf( ioQQQ, " Did not find Fe II Lya\n" );
1447  ShowMe();
1448  cdEXIT(EXIT_FAILURE);
1449  }
1450 
1451  /* find total Lya for comparison */
1452  if( cdLine("H 1",1215.67f,&relhm,&absint)<=0 )
1453  {
1454  fprintf( ioQQQ, " Did not find Lya\n" );
1455  ShowMe();
1456  cdEXIT(EXIT_FAILURE);
1457  }
1458 
1459  if( relhm > 0. )
1460  {
1461  ratio = fedest/(fedest + relhm);
1462  if( ratio > 0.1 )
1463  {
1464  sprintf( chLine, " !Fe II destruction of Ly-a removed %.1f%% of the line.",
1465  ratio *100.);
1466  bangin(chLine);
1467  }
1468  else if( ratio > 0.01 )
1469  {
1470  sprintf( chLine, " Fe II destruction of Ly-a removed %.1f%% of the line.",
1471  ratio );
1472  notein(chLine);
1473  }
1474  }
1475 
1476  if( cdLine("H-CT",6562.81,&relhm,&absint)<=0 )
1477  {
1478  fprintf( ioQQQ, " Comment did not find H-CT H-alpha\n" );
1479  ShowMe();
1480  cdEXIT(EXIT_FAILURE);
1481  }
1482 
1483  if( HBeta > 0. )
1484  {
1485  if( relhm/HBeta > 0.01 )
1486  {
1487  sprintf( chLine,
1488  " !Mutual neutralization production of H-alpha was significant." );
1489  bangin(chLine);
1490  }
1491  }
1492 
1493  /* note about very high population in H n=2 rel to ground, set in hydrogenic */
1494  if( hydro.lgHiPop2 )
1495  {
1496  sprintf( chLine,
1497  " The population of H n=2 reached %.2e relative to the ground state.",
1498  hydro.pop2mx );
1499  notein(chLine);
1500  }
1501 
1502  /* check where diffuse emission error */
1503  for( ipISO=ipH_LIKE; ipISO<=ipHE_LIKE; ++ipISO )
1504  {
1505  for( nelem=ipISO; nelem < LIMELM; nelem++ )
1506  {
1507  if( iso_sp[ipISO][nelem].CaseBCheck > 1.5 )
1508  {
1509  sprintf( chLine,
1510  " Ratio of computed diffuse emission to case B reached %g for iso %li element %li",
1511  iso_sp[ipISO][nelem].CaseBCheck , ipISO , nelem+1 );
1512  notein(chLine);
1513  }
1514  }
1515  }
1516 
1517  /* check whether electrons were relativistic */
1518  if( thermal.thist > 1e9 )
1519  {
1520  /* >>chng 06 feb 19, from 5e9 K for warning to 1e10K. add test case at 1e10K
1521  * and don't want warning in test suite. nothing is wrong at this temp - eeff
1522  * is in correctly for relativistic temps and will eventually dominate cooling */
1523  if( thermal.thist > 1.0001e10 )
1524  {
1525  sprintf( chLine, " W-Electrons were relativistic; High TE=%.2e",
1526  thermal.thist );
1527  warnin(chLine);
1528  }
1529  else
1530  {
1531  sprintf( chLine, " C-Electrons were mildly relativistic; High TE=%.2e",
1532  thermal.thist );
1533  caunin(chLine);
1534  }
1535  }
1536 
1537  /* check on timescale for photoerosion of elements */
1538  rate = timesc.TimeErode*2e-26;
1539  if( rate > 1e-35 )
1540  {
1541  /* 2E-26 is roughly cross section for photoerosion
1542  * see
1543  * >>refer all photoerode Boyd, R., & Ferland, G.J. ApJ, 318, L21. */
1544  ts = (1./rate)/3e7;
1545  if( ts < 1e3 )
1546  {
1547  sprintf( chLine, " !Timescale-photoerosion of Fe=%.2e yr",
1548  ts );
1549  bangin(chLine);
1550  }
1551  else if( ts < 1e9 )
1552  {
1553  sprintf( chLine, " Timescale-photoerosion of Fe=%.2e yr",
1554  ts );
1555  notein(chLine);
1556  }
1557  }
1558 
1559  /* check whether Compton heating was significant */
1560  comfrc = rfield.comtot/SDIV(thermal.power);
1561  if( comfrc > 0.01 )
1562  {
1563  sprintf( chLine, " Compton heating was %5.1f%% of the total.",
1564  comfrc*100. );
1565  notein(chLine);
1566  }
1567 
1568  /* check on relative importance of induced Compton heating */
1569  if( comfrc > 0.01 && rfield.cinrat > 0.05 )
1570  {
1571  sprintf( chLine,
1572  " !Induced Compton heating was %.2e of the total Compton heating.",
1573  rfield.cinrat );
1574  bangin(chLine);
1575  }
1576 
1577  /* check whether equilibrium timescales are short rel to Hubble time */
1578  if( timesc.tcmptn > 5e17 )
1579  {
1580  if( comfrc > 0.05 )
1581  {
1582  sprintf( chLine,
1583  " C-Compton cooling is significant and the equilibrium timescale (%.2e s) is longer than the Hubble time.",
1584  timesc.tcmptn );
1585  caunin(chLine);
1586  }
1587  else
1588  {
1589  sprintf( chLine,
1590  " Compton cooling equilibrium timescale (%.2e s) is longer than Hubble time.",
1591  timesc.tcmptn );
1592  notein(chLine);
1593  }
1594  }
1595 
1596  if( timesc.time_therm_long > 5e17 )
1597  {
1598  sprintf( chLine,
1599  " C-Thermal equilibrium timescale, %.2e s, longer than Hubble time; this cloud is not time-steady.",
1600  timesc.time_therm_long );
1601  caunin(chLine);
1602  }
1603 
1604  /* check whether model large relative to Jeans length
1605  * DMEAN is mean density (gm per cc)
1606  * mean temp is weighted by mass density */
1607  if( log10(radius.depth) > colden.rjnmin )
1608  {
1609  /* AJMIN is minimum Jeans mass, log in grams */
1610  aj = exp10((double)colden.ajmmin - log10(SOLAR_MASS));
1611  if( strcmp(dense.chDenseLaw,"CPRE") == 0 )
1612  {
1613  sprintf( chLine,
1614  " C-Cloud thicker than smallest Jeans length=%8.2ecm; stability problems? (smallest Jeans mass=%8.2eMo)",
1615  exp10(colden.rjnmin), aj );
1616  caunin(chLine);
1617  }
1618  else
1619  {
1620  sprintf( chLine,
1621  " Cloud thicker than smallest Jeans length=%8.2ecm; stability problems? (smallest Jeans mass=%8.2eMo)",
1622  exp10(colden.rjnmin), aj );
1623  notein(chLine);
1624  }
1625  }
1626 
1627  /* check whether grains too hot to survive */
1628  for( size_t nd=0; nd < gv.bin.size(); nd++ )
1629  {
1630  if( gv.bin[nd]->nDustFunc != DF_SUBLIMATION && gv.bin[nd]->TeGrainMax > gv.bin[nd]->Tsublimat )
1631  {
1632  sprintf( chLine,
1633  " W-Maximum temperature of grain %s was %.2eK, above its sublimation temperature, %.2eK.",
1634  gv.bin[nd]->chDstLab, gv.bin[nd]->TeGrainMax,
1635  gv.bin[nd]->Tsublimat );
1636  warnin(chLine);
1637  }
1638  else if( gv.bin[nd]->nDustFunc != DF_SUBLIMATION && gv.bin[nd]->TeGrainMax > gv.bin[nd]->Tsublimat*0.9 )
1639  {
1640  sprintf( chLine,
1641  " C-Maximum temperature of grain %s was %.2eK, near its sublimation temperature, %.2eK.",
1642  gv.bin[nd]->chDstLab, gv.bin[nd]->TeGrainMax,
1643  gv.bin[nd]->Tsublimat );
1644  caunin(chLine);
1645  }
1646  }
1647 
1648  if( gv.lgNegGrnDrg )
1649  {
1650  sprintf( chLine, " !Grain drag force <0." );
1651  bangin(chLine);
1652  }
1653 
1654  /* largest relative number of electrons donated by grains */
1655  if( gv.GrnElecDonateMax > 0.05 )
1656  {
1657  sprintf( chLine,
1658  " !Grains donated %5.1f%% of the total electrons in some regions.",
1659  gv.GrnElecDonateMax*100. );
1660  bangin(chLine);
1661  }
1662  else if( gv.GrnElecDonateMax > 0.005 )
1663  {
1664  sprintf( chLine,
1665  " Grains donated %5.1f%% of the total electrons in some regions.",
1666  gv.GrnElecDonateMax*100. );
1667  notein(chLine);
1668  }
1669 
1670  /* largest relative number of electrons on grain surface */
1671  if( gv.GrnElecHoldMax > 0.05 )
1672  {
1673  sprintf( chLine,
1674  " !Grains contained %5.1f%% of the total electrons in some regions.",
1675  gv.GrnElecHoldMax*100. );
1676  bangin(chLine);
1677  }
1678  else if( gv.GrnElecHoldMax > 0.005 )
1679  {
1680  sprintf( chLine,
1681  " Grains contained %5.1f%% of the total electrons in some regions.",
1682  gv.GrnElecHoldMax*100. );
1683  notein(chLine);
1684  }
1685 
1686  /* is photoelectric heating of gas by photoionization of grains important */
1687  if( gv.dphmax > 0.5 )
1688  {
1689  sprintf( chLine,
1690  " !Local grain-gas photoelectric heating rate reached %5.1f%% of the total.",
1691  gv.dphmax*100. );
1692  bangin(chLine);
1693  }
1694  else if( gv.dphmax > 0.05 )
1695  {
1696  sprintf( chLine,
1697  " Local grain-gas photoelectric heating rate reached %5.1f%% of the total.",
1698  gv.dphmax*100. );
1699  notein(chLine);
1700  }
1701 
1702  if( gv.TotalDustHeat/SDIV(thermal.power) > 0.01 )
1703  {
1704  sprintf( chLine,
1705  " Global grain photoelectric heating of gas was%5.1f%% of the total.",
1706  gv.TotalDustHeat/thermal.power*100. );
1707  notein(chLine);
1708  if( gv.TotalDustHeat/thermal.power > 0.25 )
1709  {
1710  sprintf( chLine,
1711  " !Grain photoelectric heating is VERY important." );
1712  bangin(chLine);
1713  }
1714  }
1715 
1716  /* grain-gas collisional cooling of gas */
1717  if( gv.dclmax > 0.05 )
1718  {
1719  sprintf( chLine,
1720  " Local grain-gas cooling of gas rate reached %5.1f%% of the total.",
1721  gv.dclmax*100. );
1722  notein(chLine);
1723  }
1724 
1725  /* check how H2 chemistry network performed */
1726  if( h2.renorm_max > 1.05 )
1727  {
1728  if( h2.renorm_max > 1.2 )
1729  {
1730  sprintf( chLine,
1731  " !The large H2 molecule - main chemistry network renormalization factor reached %.2f.",
1732  h2.renorm_max);
1733  bangin(chLine);
1734  }
1735  else
1736  {
1737  sprintf( chLine,
1738  " The large H2 molecule - main chemistry network renormalization factor reached %.2f.",
1739  h2.renorm_max);
1740  notein(chLine);
1741  }
1742  }
1743  if( h2.renorm_min < 0.95 )
1744  {
1745  if( h2.renorm_min < 0.8 )
1746  {
1747  sprintf( chLine,
1748  " !The large H2 molecule - main chemistry network renormalization factor reached %.2f.",
1749  h2.renorm_min);
1750  bangin(chLine);
1751  }
1752  else
1753  {
1754  sprintf( chLine,
1755  " The large H2 molecule - main chemistry network renormalization factor reached %.2f.",
1756  h2.renorm_min);
1757  notein(chLine);
1758  }
1759  }
1760 
1761  /* check whether photodissociation of H_2^+ molecular ion was important */
1762  if( hmi.h2pmax > 0.10 )
1763  {
1764  sprintf( chLine,
1765  " !The local H2+ photodissociation heating rate reached %5.1f%% of the total heating.",
1766  hmi.h2pmax*100. );
1767  bangin(chLine);
1768  }
1769 
1770  else if( hmi.h2pmax > 0.01 )
1771  {
1772  sprintf( chLine,
1773  " The local H2+ photodissociation heating rate reached %.1f%% of the total heating.",
1774  hmi.h2pmax*100. );
1775  notein(chLine);
1776  }
1777 
1778  /* check whether photodissociation of molecular hydrogen (H2)was important */
1779  if( hmi.h2dfrc > 0.1 )
1780  {
1781  sprintf( chLine,
1782  " !The local H2 photodissociation heating rate reached %.1f%% of the total heating.",
1783  hmi.h2dfrc*100. );
1784  bangin(chLine);
1785  }
1786  else if( hmi.h2dfrc > 0.01 )
1787  {
1788  sprintf( chLine,
1789  " The local H2 photodissociation heating rate reached %.1f%% of the total heating.",
1790  hmi.h2dfrc*100. );
1791  notein(chLine);
1792  }
1793 
1794  /* check whether cooling by molecular hydrogen (H2) was important */
1795  if( hmi.h2line_cool_frac > 0.1 )
1796  {
1797  sprintf( chLine,
1798  " !The local H2 cooling rate reached %.1f%% of the local cooling.",
1799  hmi.h2line_cool_frac*100. );
1800  bangin(chLine);
1801  }
1802  else if( hmi.h2line_cool_frac > 0.01 )
1803  {
1804  sprintf( chLine,
1805  " The local H2 cooling rate reached %.1f%% of the local cooling.",
1806  hmi.h2line_cool_frac*100. );
1807  notein(chLine);
1808  }
1809 
1810  if( hmi.h2dtot/SDIV(thermal.power) > 0.01 )
1811  {
1812  sprintf( chLine,
1813  " Global H2 photodissociation heating of gas was %.1f%% of the total heating.",
1814  hmi.h2dtot/thermal.power*100. );
1815  notein(chLine);
1816  if( hmi.h2dtot/thermal.power > 0.25 )
1817  {
1818  sprintf( chLine, " H2 photodissociation heating is VERY important." );
1819  notein(chLine);
1820  }
1821  }
1822 
1823  /* check whether photodissociation of carbon monoxide (co) was important */
1824  if( co.codfrc > 0.25 )
1825  {
1826  sprintf( chLine,
1827  " !Local CO photodissociation heating rate reached %.1f%% of the total.",
1828  co.codfrc*100. );
1829  bangin(chLine);
1830  }
1831  else if( co.codfrc > 0.05 )
1832  {
1833  sprintf( chLine,
1834  " Local CO photodissociation heating rate reached %.1f%% of the total.",
1835  co.codfrc*100. );
1836  notein(chLine);
1837  }
1838 
1839  if( co.codtot/SDIV(thermal.power) > 0.01 )
1840  {
1841  sprintf( chLine,
1842  " Global CO photodissociation heating of gas was %.1f%% of the total.",
1843  co.codtot/thermal.power*100. );
1844  notein(chLine);
1845  if( co.codtot/thermal.power > 0.25 )
1846  {
1847  sprintf( chLine, " CO photodissociation heating is VERY important." );
1848  notein(chLine);
1849  }
1850  }
1851 
1852  if( thermal.lgEdnGTcm )
1853  {
1854  sprintf( chLine,
1855  " Energy density of radiation field was greater than the Compton temperature. Is this physical?" );
1856  notein(chLine);
1857  }
1858 
1859  /* was cooling due to induced recombination important? */
1860  if( hydro.cintot/SDIV(thermal.power) > 0.01 )
1861  {
1862  sprintf( chLine, " Induced recombination cooling was %.1f%% of the total.",
1863  hydro.cintot/thermal.power*100. );
1864  notein(chLine);
1865  }
1866 
1867  /* check whether free-free heating was significant */
1868  if( thermal.FreeFreeTotHeat/SDIV(thermal.power) > 0.1 )
1869  {
1870  sprintf( chLine, " !Free-free heating was %.1f%% of the total.",
1871  thermal.FreeFreeTotHeat/thermal.power*100. );
1872  bangin(chLine);
1873  }
1874  else if( thermal.FreeFreeTotHeat/SDIV(thermal.power) > 0.01 )
1875  {
1876  sprintf( chLine, " Free-free heating was %.1f%% of the total.",
1877  thermal.FreeFreeTotHeat/thermal.power*100. );
1878  notein(chLine);
1879  }
1880 
1881  /* was heating due to H- absorption important? */
1882  if( hmi.hmitot/SDIV(thermal.power) > 0.01 )
1883  {
1884  sprintf( chLine, " H- absorption heating was %.1f%% of the total.",
1885  hmi.hmitot/SDIV(thermal.power)*100. );
1886  notein(chLine);
1887  }
1888 
1889  /* water destruction rate was zero */
1890  if( mole_global.lgH2Ozer )
1891  {
1892  sprintf( chLine, " Water destruction rate zero." );
1893  notein(chLine);
1894  }
1895 
1896  /* check for negative optical depths,
1897  * optical depth in excited state helium lines */
1898  small = 0.;
1899  imas = 0;
1900  isav = 0;
1901  j = 0;
1902  for( nelem=0; nelem<LIMELM; ++nelem )
1903  {
1904  if( dense.lgElmtOn[nelem] )
1905  {
1906  /* >>chng 06 aug 28, from numLevels_max to _local. */
1907  for( ipLo=ipH2p; ipLo < (iso_sp[ipH_LIKE][nelem].numLevels_local - 1); ipLo++ )
1908  {
1909  /* >>chng 06 aug 28, from numLevels_max to _local. */
1910  for( ipHi=ipLo + 1; ipHi < iso_sp[ipH_LIKE][nelem].numLevels_local; ipHi++ )
1911  {
1912  if( iso_sp[ipH_LIKE][nelem].trans(ipHi,ipLo).Emis().Aul() <= iso_ctrl.SmallA )
1913  continue;
1914 
1915  if( iso_sp[ipH_LIKE][nelem].trans(ipHi,ipLo).Emis().TauIn() < (realnum)small )
1916  {
1917  small = iso_sp[ipH_LIKE][nelem].trans(ipHi,ipLo).Emis().TauIn();
1918  imas = ipHi;
1919  j = ipLo;
1920  isav = nelem;
1921  }
1922  }
1923  }
1924  }
1925  }
1926 
1927  if( small < -0.05 )
1928  {
1929  sprintf( chLine,
1930  " !Some hydrogenic lines mased, species was %2s%2ld, smallest tau was %.2e, transition %li-%li",
1931  elementnames.chElementSym[isav],
1932  isav+1,small, imas , j );
1933  bangin(chLine);
1934  }
1935 
1936  /* check for negative opacities */
1937  if( opac.lgOpacNeg )
1938  {
1939  sprintf( chLine, " !Some opacities were negative - the SET NEGOPC command will save which ones." );
1940  bangin(chLine);
1941  }
1942 
1943  /* now check continua */
1944  small = 0.;
1945  imas = 0;
1946  isav = 0;
1947  for( nelem=0; nelem<LIMELM; ++nelem )
1948  {
1949  if( dense.lgElmtOn[nelem] )
1950  {
1951  /* >>chng 06 aug 28, from numLevels_max to _local. */
1952  for( i=0; i < iso_sp[ipH_LIKE][nelem].numLevels_local; i++ )
1953  {
1954  if( opac.TauAbsGeo[0][iso_sp[ipH_LIKE][nelem].fb[i].ipIsoLevNIonCon-1] < -0.001 )
1955  {
1956  small = MIN2(small,(double)opac.TauAbsGeo[0][iso_sp[ipH_LIKE][nelem].fb[i].ipIsoLevNIonCon-1]);
1957  imas = i;
1958  isav = nelem;
1959  }
1960  }
1961  }
1962  }
1963 
1964  if( small < -0.05 )
1965  {
1966  sprintf( chLine, " !Some hydrogenic (%2s%2ld) continua optical depths were negative; smallest=%.2e level=%3ld",
1967  elementnames.chElementSym[isav],
1968  isav+1,
1969  small, imas );
1970  bangin(chLine);
1971  }
1972 
1973  /* check whether any continuum optical depths are negative */
1974  nneg = 0;
1975  tauneg = 0.;
1976  freqn = 0.;
1977  for( i=0; i < rfield.nflux; i++ )
1978  {
1979  if( opac.TauAbsGeo[0][i] < -0.001 )
1980  {
1981  nneg += 1;
1982  /* only remember the smallest freq, and most neg optical depth */
1983  if( nneg == 1 )
1984  freqn = rfield.anu(i);
1985  tauneg = MIN2(tauneg,(double)opac.TauAbsGeo[0][i]);
1986  }
1987  }
1988 
1989  if( nneg > 0 )
1990  {
1991  sprintf( chLine, " !Some continuous optical depths <0. The lowest freq was %.3e Ryd, and a total of%4ld",
1992  freqn, nneg );
1993  bangin(chLine);
1994  sprintf( chLine, " !The smallest optical depth was %.2e",
1995  tauneg );
1996  bangin(chLine);
1997  }
1998 
1999  /* say if Balmer continuum optically thick */
2000  if( opac.TauAbsGeo[0][iso_sp[ipH_LIKE][ipHYDROGEN].fb[2].ipIsoLevNIonCon-1] > 0.05 )
2001  {
2002  sprintf( chLine, " The Balmer continuum optical depth was %.2e.",
2003  opac.TauAbsGeo[0][iso_sp[ipH_LIKE][ipHYDROGEN].fb[2].ipIsoLevNIonCon-1] );
2004  notein(chLine);
2005  }
2006 
2007  /* was correction for stimulated emission significant? */
2008  if( opac.stimax[0] > 0.02 && opac.TauAbsGeo[0][iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH1s].ipIsoLevNIonCon-1] > 0.2 )
2009  {
2010  sprintf( chLine, " The Lyman continuum stimulated emission correction to optical depths reached %.2e.",
2011  opac.stimax[0] );
2012  notein(chLine);
2013  }
2014  else if( opac.stimax[1] > 0.02 && opac.TauAbsGeo[0][iso_sp[ipH_LIKE][ipHYDROGEN].fb[2].ipIsoLevNIonCon-1] > 0.1 )
2015  {
2016  sprintf( chLine, " The Balmer continuum stimulated emission correction to optical depths reached %.2e.",
2017  opac.stimax[1] );
2018  notein(chLine);
2019  }
2020 
2021  /* say if Paschen continuum optically thick */
2022  if( opac.TauAbsGeo[0][iso_sp[ipH_LIKE][ipHYDROGEN].fb[3].ipIsoLevNIonCon-1] > 0.2 )
2023  {
2024  sprintf( chLine,
2025  " The Paschen continuum optical depth was %.2e.",
2026  opac.TauAbsGeo[0][iso_sp[ipH_LIKE][ipHYDROGEN].fb[3].ipIsoLevNIonCon-1] );
2027  notein(chLine);
2028  }
2029 
2030  /* some comments about near IR total optical depth */
2031  if( opac.TauAbsGeo[0][0] > 1. )
2032  {
2033  sprintf( chLine,
2034  " The continuum optical depth at the lowest energy considered (%.3e Ryd) was %.3e.",
2035  rfield.anu(0), opac.TauAbsGeo[0][0] );
2036  notein(chLine);
2037  }
2038 
2039  /* comment if optical depth to Rayleigh scattering is big
2040  * cs from VAL 76 */
2041  if( findspecieslocal("H")->column*7e-24 > 0.01 )
2042  {
2043  sprintf( chLine,
2044  " The optical depth to Rayleigh scattering at 1300A is %.2e",
2045  findspecieslocal("H")->column*6.71e-24 );
2046  notein(chLine);
2047  }
2048 
2049  if( findspecieslocal("H2+")->column*7e-18 > 0.1 )
2050  {
2051  sprintf( chLine,
2052  " !The optical depth to the H2+ molecular ion is %.2e",
2053  findspecieslocal("H2+")->column*7e-18 );
2054  bangin(chLine);
2055  }
2056  else if( findspecieslocal("H2+")->column*7e-18 > 0.01 )
2057  {
2058  sprintf( chLine,
2059  " The optical depth to the H2+ molecular ion is %.2e",
2060  findspecieslocal("H2+")->column*7e-18 );
2061  notein(chLine);
2062  }
2063 
2064  /* warn if optically thick to H- absorption */
2065  if( opac.thmin > 0.1 )
2066  {
2067  sprintf( chLine,
2068  " !Optical depth to negative hydrogen ion is %.2e",
2069  opac.thmin );
2070  bangin(chLine);
2071  }
2072  else if( opac.thmin > 0.01 )
2073  {
2074  sprintf( chLine,
2075  " Optical depth to negative hydrogen ion is %.2e",
2076  opac.thmin );
2077  notein(chLine);
2078  }
2079 
2080  /* check whether energy density less than background */
2081  if( phycon.TEnerDen < 2.6 )
2082  {
2083  sprintf( chLine,
2084  " !Incident radiation field energy density is less than 2.7K. Add background with CMB command." );
2085  bangin(chLine);
2086  }
2087 
2088  /* check whether CMB set at all */
2089  if( !rfield.lgCMB_set )
2090  {
2091  sprintf( chLine,
2092  " !The CMB was not included. This is added with the CMB command." );
2093  bangin(chLine);
2094  }
2095 
2096  /* incident radiation field is less than background Habing ISM field */
2097  if( rfield.lgHabing )
2098  {
2099  sprintf( chLine,
2100  " !The intensity of the incident radiation field is less than 10 times the Habing diffuse ISM field. Is this OK?" );
2101  bangin(chLine);
2102  sprintf( chLine,
2103  " ! Consider adding diffuse ISM emission with TABLE ISM command." );
2104  bangin(chLine);
2105  }
2106 
2107  /* some things dealing with molecules, or molecule formation */
2108 
2109  /* if C/O > 1 then chemistry will be carbon dominated rather than oxygen dominated */
2110  if( dense.lgElmtOn[ipOXYGEN] && dense.lgElmtOn[ipCARBON] )
2111  {
2112  if( dense.gas_phase[ipCARBON]/dense.gas_phase[ipOXYGEN] > 1. )
2113  {
2114  sprintf( chLine, " !The C/O abundance ratio, %.1f, is greater than unity. The chemistry will be carbon dominated.",
2115  dense.gas_phase[ipCARBON]/dense.gas_phase[ipOXYGEN] );
2116  bangin(chLine);
2117  }
2118  }
2119 
2120  bool lgLots_of_moles = false;
2121  bool lgLotsSolids = false;
2122  /* largest fraction in any molecule */
2123  for( i=0; i<mole_global.num_calc; ++i )
2124  {
2125  if( mole.species[i].location == NULL && ( mole_global.list[i]->isIsotopicTotalSpecies() || mole_global.list[i]->charge<0 ) )
2126  {
2127  if( mole.species[i].xFracLim > 0.1 )
2128  {
2129  sprintf( chLine, " !The fraction of %s in %s reached %.1f%% at some point in the cloud.",
2130  mole.species[i].atomLim->label().c_str(),
2131  mole_global.list[i]->label.c_str(),
2132  mole.species[i].xFracLim*100. );
2133  bangin(chLine);
2134  lgLots_of_moles = true;
2135  /* check whether molecules are on grains */
2136  if( !mole_global.list[i]->lgGas_Phase )
2137  lgLotsSolids = true;
2138  }
2139  else if( mole.species[i].xFracLim>0.01 )
2140  {
2141  sprintf( chLine, " The fraction of %s in %s reached %.2f%% at some point in the cloud.",
2142  mole.species[i].atomLim->label().c_str(),
2143  mole_global.list[i]->label.c_str(),
2144  mole.species[i].xFracLim*100. );
2145  notein(chLine);
2146  lgLots_of_moles = true;
2147  /* check whether molecules are on grains */
2148  if( !mole_global.list[i]->lgGas_Phase )
2149  lgLotsSolids = true;
2150  }
2151  else if( mole.species[i].xFracLim > 1e-3 )
2152  {
2153  sprintf( chLine, " The fraction of %s in %s reached %.3f%% at some point in the cloud.",
2154  mole.species[i].atomLim->label().c_str(),
2155  mole_global.list[i]->label.c_str(),
2156  mole.species[i].xFracLim*100. );
2157  notein(chLine);
2158  /* check whether molecules are on grains */
2159  if( !mole_global.list[i]->lgGas_Phase )
2160  lgLotsSolids = true;
2161  }
2162  }
2163  }
2164 
2165  /* generate comment if molecular fraction was significant but some heavy elements are turned off */
2166  if( lgLots_of_moles )
2167  {
2168  /* find all elements that are turned off */
2169  for(nelem=ipHYDROGEN; nelem<LIMELM; ++nelem )
2170  {
2171  /* >>chng 05 dec 23, add mole_global.lgElem_in_chemistry */
2172  if( !dense.lgElmtOn[nelem] )
2173  {
2174  /* this triggers if element turned off but it is part of co chem net */
2175  sprintf( chLine,
2176  " C-Molecules are important, but %s, part of the chemistry network, is turned off.",
2177  elementnames.chElementName[nelem] );
2178  caunin(chLine);
2179  }
2180 # if 0
2181  /* this element has been turned off - now check if part of chemistry */
2182  for( i=NUM_HEAVY_MOLEC+NUM_ELEMENTS; i<NUM_COMOLE_CALC; ++i )
2183  {
2184  if( nelem==mole_global.list[i].nelem_den )
2185  {
2186  /* this triggers if element turned off but it is part of co chem net */
2187  sprintf( chLine,
2188  " C-Molecules are important, but %s, part of the chemistry network, is turned off.",
2189  elementnames.chElementName[nelem] );
2190  caunin(chLine);
2191  }
2192  }
2193 # endif
2194  }
2195  }
2196 
2197  /* say if lots of molecules on grains,
2198  * molecules with labels that *GR */
2199  if( lgLotsSolids )
2200  {
2201  sprintf( chLine, " !A significant amount of molecules condensed onto grain surfaces." );
2202  bangin(chLine);
2203  sprintf( chLine, " !These are the molecular species with \"grn\" above." );
2204  bangin(chLine);
2205  }
2206 
2207  /* bremsstrahlung optical depth */
2208  if( rfield.EnergyBremsThin > 0.09 )
2209  {
2210  sprintf( chLine, " !The cloud is optically thick at optical wavelengths, extending to %.3e Ryd =%.3eA",
2211  rfield.EnergyBremsThin, RYDLAM/rfield.EnergyBremsThin );
2212  bangin(chLine);
2213  }
2214  else if( rfield.EnergyBremsThin > 0.009 )
2215  {
2216  sprintf( chLine, " The continuum of the computed structure may be optically thick in the near infrared." );
2217  notein(chLine);
2218  }
2219 
2220  /* did model run away to very large radius? */
2221  if( radius.Radius > 1e23 && radius.Radius/radius.rinner > 10. )
2222  {
2223  sprintf( chLine, " Is an outer radius of %.2e reasonable?",
2224  radius.Radius );
2225  notein(chLine);
2226  }
2227 
2228  /* following set true in RT_line_one_tauinc if maser capped at tau = -1 */
2229  if( rt.lgMaserCapHit )
2230  {
2231  sprintf( chLine, " Laser maser optical depths capped in RT_line_one_tauinc." );
2232  notein(chLine);
2233  }
2234 
2235  /* following set true in adius_next if maser cap set dr */
2236  if( rt.lgMaserSetDR )
2237  {
2238  sprintf( chLine, " !Line maser set zone thickness in some zones." );
2239  bangin(chLine);
2240  }
2241 
2242  /* lgPradCap is true if radiation pressure was capped on first iteration
2243  * also check that this is a constant total pressure model */
2244  if( (pressure.lgPradCap && (strcmp(dense.chDenseLaw,"CPRE") == 0)) &&
2245  pressure.lgPres_radiation_ON )
2246  {
2247  sprintf( chLine, " Radiation pressure kept below gas pressure on this iteration." );
2248  notein(chLine);
2249  }
2250 
2251  if( pressure.RadBetaMax > 0.25 )
2252  {
2253  if( pressure.ipPradMax_line == 0 )
2254  {
2255  sprintf( chLine,
2256  " !The ratio of radiation to gas pressure reached %.2e at zone %li. Caused by Lyman alpha.",
2257  pressure.RadBetaMax,
2258  pressure.ipPradMax_nzone);
2259  bangin(chLine);
2260  }
2261  else
2262  {
2263  sprintf( chLine,
2264  " !The ratio of radiation to gas pressure reached %.2e at zone %li. "
2265  "Caused by line number %ld, label %s",
2266  pressure.RadBetaMax,
2267  pressure.ipPradMax_nzone,
2268  pressure.ipPradMax_line,
2269  pressure.chLineRadPres.c_str() );
2270  bangin(chLine);
2271  }
2272  }
2273 
2274  else if( pressure.RadBetaMax > 0.025 )
2275  {
2276  if( pressure.ipPradMax_line == 0 )
2277  {
2278  sprintf( chLine,
2279  " The ratio of radiation to gas pressure reached %.2e at zone %li. Caused by Lyman alpha.",
2280  pressure.RadBetaMax,
2281  pressure.ipPradMax_nzone);
2282  notein(chLine);
2283  }
2284  else
2285  {
2286  sprintf( chLine,
2287  " The ratio of radiation to gas pressure reached %.2e at zone %li. "
2288  "Caused by line number %ld, label %s",
2289  pressure.RadBetaMax,
2290  pressure.ipPradMax_nzone,
2291  pressure.ipPradMax_line,
2292  pressure.chLineRadPres.c_str() );
2293  notein(chLine);
2294  }
2295  }
2296 
2297  if( opac.telec >= 5. )
2298  {
2299  sprintf( chLine, " W-The model is optically thick to electron "
2300  "scattering; tau=%.2e Cloudy is NOT intended for this regime.",
2301  opac.telec );
2302  warnin(chLine);
2303  }
2304  else if( opac.telec > 2.0 )
2305  {
2306  sprintf( chLine, " C-The model is moderately optically thick to electron scattering; tau=%.1f",
2307  opac.telec );
2308  caunin(chLine);
2309  }
2310  else if( opac.telec > 0.1 )
2311  {
2312  sprintf( chLine, " !The model has modest optical depth to electron scattering; tau=%.2f",
2313  opac.telec );
2314  bangin(chLine);
2315  }
2316  else if( opac.telec > 0.01 )
2317  {
2318  sprintf( chLine, " The optical depth to electron scattering is %.3f",
2319  opac.telec );
2320  notein(chLine);
2321  }
2322 
2323  /* optical depth to 21 cm */
2324  if( HFLines[0].Emis().TauIn() > 0.5 )
2325  {
2326  sprintf( chLine, " !The optical depth in the H I 21 cm line is %.2e",HFLines[0].Emis().TauIn() );
2327  bangin(chLine);
2328  }
2329 
2330  /* comment if level2 lines are off - they are used to pump excited states
2331  * of ground term by UV light */
2332  if( nWindLine==0 )
2333  {
2334  /* generate comment */
2335  sprintf( chLine, " !The level2 lines are disabled. UV pumping of excited levels within ground terms is not treated." );
2336  bangin(chLine);
2337  }
2338 
2339  /* check on optical depth convergence of all hydrogenic lines */
2340  for( nelem=0; nelem < LIMELM; nelem++ )
2341  {
2342  if( dense.lgElmtOn[nelem] && !dynamics.lgTimeDependentStatic )
2343  {
2344  if( iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauIn() > 0.2 )
2345  {
2346  differ = fabs(1.-iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauTot()/
2347  (iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauIn()*rt.DoubleTau))*100.;
2348 
2349  /* check whether H-alpha optical depth changed by much on last iteration
2350  * no tolerance can be finer than autocv, the tolerance on the
2351  * iterate to convergence command. It is 15% */
2352  if( ((iterations.lgLastIt && iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauIn() > 0.8) &&
2353  differ > 20.) && wind.lgStatic() )
2354  {
2355  sprintf( chLine,
2356  " C-This is the last iteration and %2s%2ld Bal(a) optical depth"
2357  " changed by %.1f%% (was %.2e). Use the ITERATE command to do more iterations.",
2358  elementnames.chElementSym[nelem],
2359  nelem+1, differ,
2360  iso_sp[ipH_LIKE][nelem].trans(ipH3p,ipH2s).Emis().TauTot() );
2361  caunin(chLine);
2362  iterations.lgIterAgain = true;
2363  }
2364 
2365  /* only check on Lya convergence if Balmer lines are thick */
2366  if( iso_sp[ipH_LIKE][nelem].trans(ipH2p,ipH1s).Emis().TauIn() > 0. )
2367  {
2368  differ = fabs(1.-iso_sp[ipH_LIKE][nelem].trans(ipH2p,ipH1s).Emis().TauTot()/
2369  (iso_sp[ipH_LIKE][nelem].trans(ipH2p,ipH1s).Emis().TauIn()*rt.DoubleTau))*100.;
2370 
2371  /* check whether Lya optical depth changed on last iteration
2372  * no tolerance can be finer than autocv, the tolerance on the
2373  * iterate to convergence command. It is 15% */
2374  if( ((iterations.lgLastIt && iso_sp[ipH_LIKE][nelem].trans(ipH2p,ipH1s).Emis().TauIn() > 0.8) &&
2375  differ > 25.) && wind.lgStatic() )
2376  {
2377  sprintf( chLine,
2378  " C-This is the last iteration and %2s%2ld Ly(a) optical depth"
2379  " changed by %.1f%% (was %.2e). Use the ITERATE command to do more iterations.",
2380  elementnames.chElementSym[nelem],
2381  nelem+1,differ, iso_sp[ipH_LIKE][nelem].trans(ipH2p,ipH1s).Emis().TauTot() );
2382  caunin(chLine);
2383  iterations.lgIterAgain = true;
2384  }
2385  }
2386  }
2387  }
2388  }
2389 
2390  /* check whether sphere was set if dr/r large */
2391  if( radius.Radius/radius.rinner > 2. && !geometry.lgSphere )
2392  {
2393  sprintf( chLine, " C-R(out)/R(in)=%.2e and SPHERE was not set.",
2394  radius.Radius/radius.rinner );
2395  caunin(chLine);
2396  }
2397 
2398  /* check if thin in hydrogen or helium continua, but assumed to be thick */
2399  if( iterations.lgLastIt && !opac.lgCaseB )
2400  {
2401 
2402  /* check if thin in Lyman continuum, and assumed thick */
2403  if( rfield.nflux > iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH1s].ipIsoLevNIonCon )
2404  {
2405  if( opac.TauAbsGeo[0][iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH1s].ipIsoLevNIonCon-1] < 2. &&
2406  opac.TauAbsGeo[1][iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH1s].ipIsoLevNIonCon-1] > 2. )
2407  {
2408  sprintf( chLine, " C-The H Lyman continuum is thin, and I assumed"
2409  " that it was thick. Use the ITERATE command to do more iterations." );
2410  caunin(chLine);
2411  iterations.lgIterAgain = true;
2412  }
2413  }
2414 
2415  /* check on the He+ ionizing continuum */
2416  if( rfield.nflux > iso_sp[ipH_LIKE][ipHELIUM].fb[ipH1s].ipIsoLevNIonCon && dense.lgElmtOn[ipHELIUM] )
2417  {
2418  if( (opac.TauAbsGeo[0][iso_sp[ipH_LIKE][ipHELIUM].fb[ipH1s].ipIsoLevNIonCon-1] < 2. &&
2419  opac.TauAbsGeo[1][iso_sp[ipH_LIKE][ipHELIUM].fb[ipH1s].ipIsoLevNIonCon-1] > 2.) )
2420  {
2421  sprintf( chLine,
2422  " C-The He II continuum is thin and I assumed that it was thick."
2423  " Use the ITERATE command to do more iterations." );
2424  caunin(chLine);
2425  iterations.lgIterAgain = true;
2426  }
2427  }
2428 
2429  if( rfield.nflux > iso_sp[ipHE_LIKE][ipHELIUM].fb[0].ipIsoLevNIonCon && dense.lgElmtOn[ipHELIUM] )
2430  {
2431  if( (opac.TauAbsGeo[0][iso_sp[ipHE_LIKE][ipHELIUM].fb[0].ipIsoLevNIonCon-1] < 2. &&
2432  opac.TauAbsGeo[1][iso_sp[ipHE_LIKE][ipHELIUM].fb[0].ipIsoLevNIonCon-1] > 2.) )
2433  {
2434  sprintf( chLine,
2435  " C-The He I continuum is thin and I assumed that it was thick."
2436  " Use the ITERATE command to do more iterations." );
2437  caunin(chLine);
2438  iterations.lgIterAgain = true;
2439  }
2440  }
2441  }
2442 
2443  /* check whether column density changed by much on this iteration */
2444  if( iteration > 1 )
2445  {
2446  if( colden.colden_old[ipCOL_HTOT] <= 0. )
2447  {
2448  fprintf( ioQQQ, " colden_old is insane in PrtComment.\n" );
2449  ShowMe();
2450  cdEXIT(EXIT_FAILURE);
2451  }
2452 
2453  differ = fabs(1.-colden.colden[ipCOL_HTOT]/
2454  colden.colden_old[ipCOL_HTOT]);
2455 
2456  if( differ > 0.1 && differ <= 0.3 )
2457  {
2458  sprintf( chLine,
2459  " The H column density changed by %.2e%% between this and previous iteration.",
2460  differ*100. );
2461  notein(chLine);
2462  }
2463 
2464  else if( differ > 0.3 )
2465  {
2466  if( iterations.lgLastIt )
2467  {
2468  sprintf( chLine,
2469  " C-The H column density changed by %.2e%% and this is the last iteration. What happened?",
2470  differ*100. );
2471  caunin(chLine);
2472  }
2473  else
2474  {
2475  sprintf( chLine,
2476  " !The H column density changed by %.2e%% What happened?",
2477  differ*100. );
2478  bangin(chLine);
2479  }
2480  }
2481 
2482  /* check on H2 column density, but only if significant fraction of H is molecular */
2483  if( (findspecieslocal("H2")->column+findspecieslocal("H2*")->column)/SDIV(colden.colden[ipCOL_HTOT]) > 1e-5 )
2484  {
2485  differ = fabs(1.-findspecieslocal("H2")->column/
2486  SDIV(findspecieslocal("H2")->column_old));
2487 
2488  if( differ > 0.1 && differ <= 0.3 )
2489  {
2490  sprintf( chLine,
2491  " The H2 column density changed by %.2e%% between this and previous iteration.",
2492  differ*100. );
2493  notein(chLine);
2494  }
2495 
2496  else if( differ > 0.3 )
2497  {
2498  if( iterations.lgLastIt )
2499  {
2500  sprintf( chLine,
2501  " C-The H2 column density changed by %.2e%% and this is the last iteration. What happened?",
2502  differ*100. );
2503  caunin(chLine);
2504  }
2505  else
2506  {
2507  sprintf( chLine,
2508  " !The H2 column density changed by %.2e%% What happened?",
2509  differ*100. );
2510  bangin(chLine);
2511  }
2512  }
2513  }
2514  }
2515 
2516  /* say if rad pressure caused by la and la optical depth changed too much */
2517  differ = fabs(1.-iso_sp[ipH_LIKE][ipHYDROGEN].trans(ipH2p,ipH1s).Emis().TauIn()/
2518  SDIV(iso_sp[ipH_LIKE][ipHYDROGEN].trans(ipH2p,ipH1s).Emis().TauTot()))*100.;
2519 
2520  if( iterations.lgLastIt && (pressure.RadBetaMax > 0.1) &&
2521  (differ > 50.) && (pressure.ipPradMax_line == 1) && (pressure.lgPres_radiation_ON) &&
2522  wind.lgStatic() )
2523  {
2524  sprintf( chLine, " C-This is the last iteration, radiation pressure was significant, and the L-a optical depth changed by %7.2f%% (was %.2e)",
2525  differ, iso_sp[ipH_LIKE][ipHYDROGEN].trans(ipH2p,ipH1s).Emis().TauTot() );
2526  caunin(chLine);
2527  }
2528 
2529  /* caution that 21 cm spin temperature is incorrect when Lya optical depth
2530  * scale is overrun */
2531  if( iterations.lgLastIt &&
2532  ( iso_sp[ipH_LIKE][ipHYDROGEN].trans(ipH2p,ipH1s).Emis().TauTot() * 1.02 -
2533  iso_sp[ipH_LIKE][ipHYDROGEN].trans(ipH2p,ipH1s).Emis().TauIn() ) < 0. )
2534  {
2535  sprintf( chLine, " C-The Lya optical depth scale was overrun and this is the last iteration - Tspin(21 cm) is not valid." );
2536  caunin(chLine);
2537  sprintf( chLine, " C-Another iteration is needed for Tspin(21 cm) to be valid. Use the ITERATE command." );
2538  caunin(chLine);
2539  }
2540 
2541  /* say if la rad pressure capped by thermalization length */
2542  if( pressure.lgPradDen )
2543  {
2544  sprintf( chLine, " Line radiation pressure capped by thermalization length." );
2545  notein(chLine);
2546  }
2547 
2548  /* print te failures */
2549  nline = MIN2(conv.nTeFail,10);
2550  if( conv.nTeFail != 0 )
2551  {
2552  if( conv.failmx < 0.1 )
2553  {
2554  long _o = sprintf( chLine, " There were %ld minor temperature failures. zones:",
2555  conv.nTeFail );
2556  /* don't know how many zones we will save, there are nline,
2557  * hence this use of pointer arith */
2558  for( i=0; i < nline; i++ )
2559  {
2560  _o += sprintf( chLine+_o, " %ld", conv.ifailz[i] );
2561  }
2562  notein(chLine);
2563  }
2564  else
2565  {
2566  sprintf( chLine,
2567  " !There were %ld temperature failures, and some were large. The largest was %.1f%%. What happened?",
2568  conv.nTeFail, conv.failmx*100. );
2569  bangin(chLine);
2570 
2571  /* don't know how many zones we will save, there are nline,
2572  * hence this use of pointer arith */
2573  long _o = sprintf( chLine , " !The zones were" );
2574  for( i=0; i < nline; i++ )
2575  {
2576  _o += sprintf( chLine+_o, " %ld", conv.ifailz[i] );
2577  }
2578  bangin(chLine);
2579 
2580  if( struc.testr[0] > 8e4 && phycon.te < 6e5 )
2581  {
2582  sprintf( chLine, " !I think they may have been caused by the change from hot to nebular gas phase. The physics of this is unclear." );
2583  bangin(chLine);
2584  }
2585  }
2586  }
2587 
2588  /* check for temperature jumps */
2589  big_ion_jump = 0.;
2590  j = 0;
2591  for( i=1; i < nzone; i++ )
2592  {
2593  big = fabs(1.-struc.testr[i-1]/struc.testr[i]);
2594  if( big > big_ion_jump )
2595  {
2596  j = i;
2597  big_ion_jump = big;
2598  }
2599  }
2600 
2601  if( big_ion_jump > 0.2 )
2602  {
2603  /* this is a sanity check, but only do it if jump detected */
2604  if( j < 1 )
2605  {
2606  fprintf( ioQQQ, " j too small big jump check\n" );
2607  ShowMe();
2608  cdEXIT(EXIT_FAILURE);
2609  }
2610 
2611  if( big_ion_jump > 0.4 )
2612  {
2613  sprintf( chLine, " C-A temperature discontinuity occurred"
2614  " from %.2eK (zone %ld) to %.2eK (zone %ld).",
2615  struc.testr[j-1], j, struc.testr[j], j+1 );
2616  caunin(chLine);
2617  /* check if the second temperature is between 100 and 1000K */
2618  /* >>chng 05 nov 07, test second not first temperature since second
2619  * will be lower of the two */
2620  /*if( struc.testr[j-1] < 1000. && struc.testr[j-1]>100. )*/
2621  if( struc.testr[j]>100. && struc.testr[j] < 1000. )
2622  {
2623  sprintf( chLine, " C-This was probably due to a thermal front." );
2624  caunin(chLine);
2625  }
2626  }
2627  else if( big_ion_jump > 0.2 )
2628  {
2629  sprintf( chLine, " !A temperature discontinuity occurred"
2630  " from %.2eK (zone %ld) to %.2eK (zone %ld).",
2631  struc.testr[j-1], j, struc.testr[j], j+1 );
2632  bangin(chLine);
2633  /* check if the second temperature is between 100 and 1000K */
2634  /* >>chng 05 nov 07, test second not first temperature since second
2635  * will be lower of the two */
2636  /*if( struc.testr[j-1] < 1000. && struc.testr[j-1]>100. )*/
2637  if( struc.testr[j]>100. && struc.testr[j] < 1000. )
2638  {
2639  sprintf( chLine, " !This was probably due to a thermal front." );
2640  bangin(chLine);
2641  }
2642  }
2643  }
2644 
2645  /* check for largest error in local electron density */
2646  if( fabs(conv.BigEdenError) > conv.EdenErrorAllowed )
2647  {
2648  /* this only produces a warning if not the very last zone */
2649  if( fabs(conv.BigEdenError) > conv.EdenErrorAllowed*20. && dense.nzEdenBad !=
2650  nzone )
2651  {
2652  sprintf( chLine, " W-The local error in the electron density reached %.1f%% at zone %ld",
2653  conv.BigEdenError*100, dense.nzEdenBad );
2654  warnin(chLine);
2655  }
2656  else if( fabs(conv.BigEdenError) > conv.EdenErrorAllowed*5. )
2657  {
2658  sprintf( chLine, " C-The local error in the electron density reached %.1f%% at zone %ld",
2659  conv.BigEdenError*100, dense.nzEdenBad );
2660  caunin(chLine);
2661  }
2662  else
2663  {
2664  sprintf( chLine, " The local error in the electron density reached %.1f%% at zone %ld",
2665  conv.BigEdenError*100, dense.nzEdenBad );
2666  notein(chLine);
2667  }
2668  }
2669 
2670  /* check for temperature oscillations or fluctuations*/
2671  big_ion_jump = 0.;
2672  j = 0;
2673  for( i=1; i < (nzone - 1); i++ )
2674  {
2675  big = fabs( (struc.testr[i-1] - struc.testr[i])/struc.testr[i] );
2676  bigm = fabs( (struc.testr[i] - struc.testr[i+1])/struc.testr[i] );
2677 
2678  /* this is sign of change in temperature, we are looking for change in sign */
2679  rel = ( (struc.testr[i-1] - struc.testr[i])/struc.testr[i])*
2680  ( (struc.testr[i] - struc.testr[i+1])/struc.testr[i] );
2681 
2682  if( rel < 0. && MIN2( bigm , big ) > big_ion_jump )
2683  {
2684  j = i;
2685  big_ion_jump = MIN2( bigm , big );
2686  }
2687  }
2688 
2689  if( big_ion_jump > 0.1 )
2690  {
2691  /* only do sanity check if jump detected */
2692  if( j < 1 )
2693  {
2694  fprintf( ioQQQ, " j too small bigjump2 check\n" );
2695  ShowMe();
2696  cdEXIT(EXIT_FAILURE);
2697  }
2698 
2699  if( big_ion_jump > 0.3 )
2700  {
2701  sprintf( chLine,
2702  " C-A temperature oscillation occurred by %.0f%%"
2703  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
2704  big_ion_jump*100., struc.testr[j-1], j, struc.testr[j], j+1,
2705  struc.testr[j+1], j+2 );
2706  caunin(chLine);
2707  }
2708  else if( big_ion_jump > 0.1 )
2709  {
2710  sprintf( chLine,
2711  " !A temperature oscillation occurred by %.0f%%"
2712  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
2713  big_ion_jump*100., struc.testr[j-1], j, struc.testr[j], j+1,
2714  struc.testr[j+1], j+2 );
2715  bangin(chLine);
2716  }
2717  }
2718 
2719  /* check for eden oscillations */
2720  if( strcmp(dense.chDenseLaw,"CDEN") == 0 )
2721  {
2722  j = 0;
2723  big_ion_jump = 0.;
2724  for( i=1; i < (nzone - 1); i++ )
2725  {
2726  big = (struc.ednstr[i-1] - struc.ednstr[i])/struc.ednstr[i];
2727  if( fabs(big) < conv.EdenErrorAllowed )
2728  big = 0.;
2729  bigm = (struc.ednstr[i] - struc.ednstr[i+1])/struc.ednstr[i];
2730  if( fabs(bigm) < conv.EdenErrorAllowed )
2731  bigm = 0.;
2732  if( big*bigm < 0. &&
2733  fabs(struc.ednstr[i-1]-struc.ednstr[i])/struc.ednstr[i] > big_ion_jump )
2734  {
2735  j = i;
2736  big_ion_jump = fabs(struc.ednstr[i-1]-struc.ednstr[i])/
2737  struc.ednstr[i];
2738  }
2739  }
2740 
2741  /* only check on j if there was a big jump detected, number must be
2742  * smallest jump */
2743  if( big_ion_jump > conv.EdenErrorAllowed*3. )
2744  {
2745  if( j < 1 )
2746  {
2747  fprintf( ioQQQ, " j too small bigjump3 check\n" );
2748  ShowMe();
2749  cdEXIT(EXIT_FAILURE);
2750  }
2751 
2752  if( big_ion_jump > conv.EdenErrorAllowed*10. )
2753  {
2754  sprintf( chLine, " C-An electron density oscillation occurred by %.0f%%"
2755  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
2756  big_ion_jump*100., struc.ednstr[j-1], j, struc.ednstr[j], j+1,
2757  struc.ednstr[j+1], j+2 );
2758  caunin(chLine);
2759  }
2760  else if( big_ion_jump > conv.EdenErrorAllowed*3. )
2761  {
2762  sprintf( chLine, " !An electron density oscillation occurred by %.0f%%"
2763  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
2764  big_ion_jump*100., struc.ednstr[j-1], j, struc.ednstr[j], j+1,
2765  struc.ednstr[j+1], j+2 );
2766  bangin(chLine);
2767  }
2768  }
2769  }
2770 
2771  /*prt_smooth_predictions check whether fluctuations in any predicted quantities occurred */
2772  /* >>chng 03 dec 05, add this test */
2773  prt_smooth_predictions();
2774 
2775  /**********************************************************
2776  * check that the volume integrates out ok *
2777  **********************************************************/
2778 
2779  /* this was the number 1 fed through the line integrators,
2780  * the number 1e-10 is sent to linadd in lineset1 as follows:*/
2781  /*linadd( 1.e-10 , 1 , "Unit" , 'i' );*/
2782  i = cdLine( "Unit" , 1 , &rate , &absint );
2783  ASSERT( i> 0 );
2784 
2785  /* this is now the linear vol, rel to inner radius */
2786  VolComputed = LineSave.lines[i].SumLine(0) / 1e-10;
2787 
2788  /* spherical or plane parallel case? */
2789  if( radius.Radius/radius.rinner > 1.0001 )
2790  {
2791  /* spherical case,
2792  * geometry.iEmissPower is usually 2,
2793  * and can be reset to 1 (long slit) or 0 (beam) with
2794  * slit and beam options on aperture */
2795  VolExpected = geometry.covaper*geometry.FillFac*radius.rinner/(geometry.iEmissPower+1)*
2796  ( powi( radius.Radius/radius.rinner,geometry.iEmissPower+1 ) - 1. );
2797  }
2798  else
2799  {
2800  /* plane parallel case */
2801  /* next can be zero for very thin model, depth is always positive */
2802  VolExpected = geometry.covaper*geometry.FillFac*(radius.depth-DEPTH_OFFSET);
2803  }
2804 
2805  /* now get the relative difference between computed and expected volumes */
2806  error = fabs(VolComputed - VolExpected)/SDIV(VolExpected);
2807 
2808  /* we need to ignore this test if filling factor changes with radius, or
2809  * cylinder geometry in place */
2810  if( radius.lgCylnOn || geometry.filpow!=0. )
2811  {
2812  error = 0.;
2813  }
2814 
2815  /* how large is relative error? */
2816  if( error > 0.001 && !lgAbort )
2817  {
2818  sprintf( chLine,
2819  " W-PrtComment insanity - Line unit integration did not verify \n");
2820  warnin(chLine);
2821  fprintf( ioQQQ,
2822  " PROBLEM PrtComment insanity - Line unit integration did not verify \n");
2823  fprintf( ioQQQ,
2824  " expected, derived vols were %g %g \n",
2825  VolExpected , VolComputed );
2826  fprintf( ioQQQ,
2827  " relative difference is %g, ratio is %g.\n",error,VolComputed/VolExpected);
2828  TotalInsanity();
2829  }
2830 
2831  /* next do same thing for fake continuum point propagated in highest energy cell, plus 1
2832  * =
2833  * the variable rfield.ConEmitLocal[rfield.nflux]
2834  * are set to
2835  * the number 1.e-10f * Dilution in RT_diffuse. this is the outward
2836  * local emissivity, per unit vol. It is then added to the outward beams
2837  * by the rest of the code, and then checked here.
2838  *
2839  * insanity will be detected if diffuse emission is thrown into the outward beam
2840  * in MadeDiffuse. this happens if the range of ionization encompasses the full
2841  * continuum array, up to nflux. */
2842  ConComputed = rfield.ConInterOut[rfield.nflux]/ 1e-10;
2843  /* correct for fraction that went out, as set in ZoneStart,
2844  * this should now be the volume of the emitting region */
2845  ConComputed /= ( (1. + geometry.covrt)/2. );
2846 
2847  /* we expect this to add up to the integral of unity over r^-2 */
2848  if( radius.Radius/radius.rinner < 1.001 )
2849  {
2850  /* plane parallel case, no dilution, use thickness */
2851  ConExpected = (radius.depth-DEPTH_OFFSET)*geometry.FillFac;
2852  }
2853  else
2854  {
2855  /* spherical case */
2856  ConExpected = radius.rinner*geometry.FillFac * (1. - radius.rinner/radius.Radius );
2857  }
2858  /* this is impossible */
2859  ASSERT( ConExpected > 0. );
2860 
2861  /* now get the relative difference between computed and expected volumes */
2862  error = fabs(ConComputed - ConExpected)/ConExpected;
2863 
2864  /* we need to ignore this test if filling factor changes with radius, or
2865  * cylinder geometry in place */
2866  if( radius.lgCylnOn || geometry.filpow!=0. )
2867  {
2868  error = 0.;
2869  }
2870 
2871  /* \todo 2 - These "volumes" seem to be too small by a factor of two.
2872  * rfield.ConInterOut[rfield.nflux] (hence ConComputed) and ConExpected
2873  * should be greater by a factor of 2 if comparison is really of "volume"
2874  * of 1/cc pencil. */
2875 
2876  /* how large is relative error? */
2877  if( error > 0.001 && !lgAbort)
2878  {
2879  sprintf( chLine,
2880  " W-PrtComment insanity - Continuum unit integration did not verify \n");
2881  warnin(chLine);
2882  fprintf( ioQQQ," PROBLEM PrtComment insanity - Continuum unit integration did not verify \n");
2883  fprintf( ioQQQ," exact vol= %g, derived vol= %g relative difference is %g \n",
2884  ConExpected , ConComputed ,error);
2885  fprintf( ioQQQ," ConInterOut= %g, \n",
2886  rfield.ConInterOut[rfield.nflux]);
2887  TotalInsanity();
2888  }
2889 
2890  /* final printout of warnings, cautions, notes, */
2891  cdNwcns(&lgAbort_flag,&nw,&nc,&nn,&ns,&i,&j,&dum1,&dum2);
2892 
2893  warnings.lgWarngs = ( nw > 0 );
2894  warnings.lgCautns = ( nc > 0 );
2895 
2896  if( called.lgTalk )
2897  {
2898  /* print the title of the calculation */
2899  fprintf( ioQQQ, " %s\n", input.chTitle );
2900  /* say why the calculation stopped, and indicate the geometry*/
2901  cdReasonGeo(ioQQQ);
2902  /* print all warnings */
2903  cdWarnings(ioQQQ);
2904  /* all cautions */
2905  cdCautions(ioQQQ);
2906  /* surprises, beginning with a ! */
2907  cdSurprises(ioQQQ);
2908  /* notes about the calculations */
2909  cdNotes(ioQQQ);
2910  }
2911 
2912  /* option to print warnings on special io */
2913  if( lgPrnErr )
2914  {
2915  cdWarnings(ioPrnErr);
2916  }
2917 
2918  if( trace.lgTrace )
2919  {
2920  fprintf( ioQQQ, " PrtComment returns.\n" );
2921  }
2922  return;
2923 }
2924 
2925 
2926 /*chkCaHeps check whether CaII K and H epsilon overlap */
2927 STATIC void chkCaHeps(double *totwid)
2928 {
2929  double conca,
2930  conalog ,
2931  conhe;
2932 
2933  DEBUG_ENTRY( "chkCaHeps()" );
2934 
2935  *totwid = 0.;
2936 
2937  if( !atmdat.lgdBaseSourceExists[ipCALCIUM][1] )
2938  {
2939  return;
2940  }
2941 
2942  /* pumping of CaH overlapping with Hy epsilon, 6-2 of H */
2943  if( iso_sp[ipH_LIKE][ipHYDROGEN].n_HighestResolved_local +
2944  iso_sp[ipH_LIKE][ipHYDROGEN].nCollapsed_local >= 6 )
2945  {
2946  /* this is 6P */
2947  long ip6p = iso_sp[ipH_LIKE][ipHYDROGEN].QuantumNumbers2Index[6][1][2];
2948 
2949  long id_Ca2 = -1;
2950  if( (id_Ca2 = atmdat.ipSpecIon[ipCALCIUM][1]) < 0 )
2951  {
2952  fprintf(ioQQQ,"PROBLEM: Ca II, the species defined by nelem = %i and ion = %i could not be found.\n",ipCALCIUM,2);
2953  cdEXIT(EXIT_FAILURE);
2954  }
2955 
2956  static TransitionList::iterator it;
2957  static bool lgRunOnce = true;
2958  if( lgRunOnce )
2959  {
2960  for( it = dBaseTrans[id_Ca2].begin(); it != dBaseTrans[id_Ca2].end(); ++it)
2961  {
2962  if( it->ipLo()+1 == 1 && it->ipHi()+1 == 5)
2963  {
2964  lgRunOnce = false;
2965  break;
2966  }
2967  }
2968  }
2969 
2970  realnum ca2_3969_TauIn = it->Emis().TauIn();
2971 
2972  if( ca2_3969_TauIn > 0. &&
2973  iso_sp[ipH_LIKE][ipHYDROGEN].trans(ip6p,ipH2s).Emis().TauIn() > 0. )
2974  {
2975  /* casts to double here are to prevent FPE */
2976  conca = sqrt(6.1e-5*ca2_3969_TauIn);
2977  conalog = log((double)ca2_3969_TauIn);
2978  conalog = sqrt(MAX2(1., conalog));
2979  conca = MAX2(conalog,conca);
2980 
2981  conalog = log((double)iso_sp[ipH_LIKE][ipHYDROGEN].trans(ip6p,ipH2s).Emis().TauIn());
2982  conalog = sqrt(MAX2(1.,conalog));
2983  conhe = sqrt(1.7e-6*iso_sp[ipH_LIKE][ipHYDROGEN].trans(ip6p,ipH2s).Emis().TauIn());
2984  conhe = MAX2(conalog, conhe);
2985 
2986  *totwid = 10.*conhe + 1.6*conca;
2987  }
2988  }
2989  return;
2990 }
2991 
2992 /*prt_smooth_predictions check whether fluctuations in any predicted quantities occurred */
2993 STATIC void prt_smooth_predictions(void)
2994 {
2995  long int i,
2996  nzone_oscillation,
2997  nzone_ion_jump,
2998  nzone_den_jump,
2999  nelem,
3000  ion;
3001  double BigOscillation ,
3002  big_ion_jump,
3003  big_jump,
3004  rel,
3005  big,
3006  bigm;
3007 
3008  char chLine[INPUT_LINE_LENGTH];
3009 
3010  DEBUG_ENTRY( "prt_smooth_predictions()" );
3011 
3012  /* check for ionization oscillations or fluctuations and or jumps */
3013  nzone_oscillation = 0;
3014  nzone_ion_jump = 0;
3015 
3016  for( nelem=ipHYDROGEN; nelem<LIMELM; ++nelem )
3017  {
3018  if( dense.lgElmtOn[nelem] )
3019  {
3020  for( ion=0; ion<=nelem+1; ++ion)
3021  {
3022  BigOscillation = 0.;
3023  big_ion_jump = -15.;
3024  /* >>chng 05 mar 12, add -lgAbort, since in some bad aborts current zone never evaluated */
3025  for( i=1; i < (nzone - 1-(int)lgAbort); i++ )
3026  {
3027 
3028  /* only do check if all ions are positive */
3029  if( struc.xIonDense[i-1][nelem][ion]/struc.gas_phase[i-1][nelem]>struc.dr_ionfrac_limit &&
3030  struc.xIonDense[i ][nelem][ion]/struc.gas_phase[i ][nelem]>struc.dr_ionfrac_limit &&
3031  struc.xIonDense[i+1][nelem][ion]/struc.gas_phase[i+1][nelem]>struc.dr_ionfrac_limit )
3032  {
3033 
3034  /* this is check for oscillations */
3035  big = fabs( (struc.xIonDense[i-1][nelem][ion] - struc.xIonDense[i][nelem][ion])/struc.xIonDense[i][nelem][ion] );
3036  bigm = fabs( (struc.xIonDense[i][nelem][ion] - struc.xIonDense[i+1][nelem][ion])/struc.xIonDense[i][nelem][ion] );
3037 
3038  /* this is sign of change in ionization, we are looking for change in sign */
3039  rel = ( (struc.xIonDense[i-1][nelem][ion] - struc.xIonDense[i][nelem][ion] )/struc.xIonDense[i][nelem][ion])*
3040  ( (struc.xIonDense[i][nelem][ion] - struc.xIonDense[i+1][nelem][ion])/struc.xIonDense[i][nelem][ion] );
3041 
3042  if( rel < 0. && MIN2( bigm , big ) > BigOscillation )
3043  {
3044  nzone_oscillation = i;
3045  BigOscillation = MIN2( bigm , big );
3046  }
3047 
3048  /* check whether we tripped over an ionization front - a major source
3049  * of instability in a complete linearization code like this one */
3050  /* neg sign picks up only increases in ionization */
3051  rel = -log10( (struc.xIonDense[i][nelem][ion]/struc.gas_phase[i][nelem]) /
3052  (struc.xIonDense[i+1][nelem][ion]/struc.gas_phase[i+1][nelem] ) );
3053  /* only do significant stages of ionization */
3054  if( rel > big_ion_jump )
3055  {
3056  big_ion_jump = rel;
3057  nzone_ion_jump = i;
3058  }
3059  }
3060  }
3061  /* end loop over zones,
3062  * check whether this ion and element underwent fluctuations or jump */
3063 
3064  if( BigOscillation > 0.2 )
3065  {
3066  /* only do sanity check if jump detected */
3067  if( nzone_oscillation < 1 )
3068  {
3069  fprintf( ioQQQ, " nzone_oscillation too small bigjump2 check\n" );
3070  ShowMe();
3071  cdEXIT(EXIT_FAILURE);
3072  }
3073 
3074  if( BigOscillation > 3. )
3075  {
3076  sprintf( chLine,
3077  " C-An ionization oscillation occurred, elem %.2s%2li, by %.0f%%"
3078  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
3079  elementnames.chElementSym[nelem], ion+1,
3080  BigOscillation*100.,
3081  struc.xIonDense[nzone_oscillation-1][nelem][ion]/struc.gas_phase[nzone_oscillation-1][nelem],
3082  nzone_oscillation,
3083  struc.xIonDense[nzone_oscillation][nelem][ion]/struc.gas_phase[nzone_oscillation][nelem],
3084  nzone_oscillation+1,
3085  struc.xIonDense[nzone_oscillation+1][nelem][ion]/struc.gas_phase[nzone_oscillation+1][nelem],
3086  nzone_oscillation+2 );
3087  caunin(chLine);
3088  }
3089  else if( BigOscillation > 0.7 )
3090  {
3091  sprintf( chLine,
3092  " !An ionization oscillation occurred, elem %.2s%2li, by %.0f%%"
3093  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
3094  elementnames.chElementSym[nelem], ion+1,
3095  BigOscillation*100.,
3096  struc.xIonDense[nzone_oscillation-1][nelem][ion]/struc.gas_phase[nzone_oscillation-1][nelem],
3097  nzone_oscillation,
3098  struc.xIonDense[nzone_oscillation][nelem][ion]/struc.gas_phase[nzone_oscillation][nelem],
3099  nzone_oscillation+1,
3100  struc.xIonDense[nzone_oscillation+1][nelem][ion]/struc.gas_phase[nzone_oscillation+1][nelem],
3101  nzone_oscillation+2 );
3102  bangin(chLine);
3103  }
3104  }
3105 
3106  /* big_ion_jump was a log above, convert to linear quantity */
3107  /* if no jump occurred then big_ion_jump is small and nzone_ion_jump is 0 */
3108  big_ion_jump = exp10( big_ion_jump );
3109  if( big_ion_jump > 1.5 && nzone_ion_jump > 0 )
3110  {
3111  if( big_ion_jump > 10. )
3112  {
3113  sprintf( chLine,
3114  " C-An ionization jump occurred, elem %.2s%2li, by %.0f%%"
3115  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
3116  elementnames.chElementSym[nelem], ion+1,
3117  big_ion_jump*100.,
3118  struc.xIonDense[nzone_ion_jump-1][nelem][ion]/struc.gas_phase[nzone_ion_jump-1][nelem],
3119  nzone_ion_jump,
3120  struc.xIonDense[nzone_ion_jump][nelem][ion]/struc.gas_phase[nzone_ion_jump][nelem],
3121  nzone_ion_jump+1,
3122  struc.xIonDense[nzone_ion_jump+1][nelem][ion]/struc.gas_phase[nzone_ion_jump+1][nelem],
3123  nzone_ion_jump+2 );
3124  caunin(chLine);
3125  }
3126  else
3127  {
3128  sprintf( chLine,
3129  " !An ionization jump occurred elem %.2s%2li, by %.0f%%"
3130  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
3131  elementnames.chElementSym[nelem], ion+1,
3132  big_ion_jump*100.,
3133  struc.xIonDense[nzone_ion_jump-1][nelem][ion]/struc.gas_phase[nzone_ion_jump-1][nelem],
3134  nzone_ion_jump,
3135  struc.xIonDense[nzone_ion_jump][nelem][ion]/struc.gas_phase[nzone_ion_jump][nelem],
3136  nzone_ion_jump+1,
3137  struc.xIonDense[nzone_ion_jump+1][nelem][ion]/struc.gas_phase[nzone_ion_jump+1][nelem],
3138  nzone_ion_jump+2 );
3139  bangin(chLine);
3140  }
3141  }
3142  }
3143  }
3144  }
3145 
3146  big_jump = -15;
3147  nzone_den_jump = 0;
3148 
3149  /* >>chng 05 mar 12, add -lgAbort, since in some bad aborts current zone never evaluated */
3150  for( i=1; i < (nzone - 1 - (int)lgAbort); i++ )
3151  {
3152  /* this first check is on how the total hydrogen density has changed */
3153  rel = fabs(log10( struc.gas_phase[i][ipHYDROGEN] /
3154  struc.gas_phase[i+1][ipHYDROGEN] ) );
3155  /* only do significant stages of ionization */
3156  if( rel > big_jump )
3157  {
3158  big_jump = rel;
3159  nzone_den_jump = i;
3160  }
3161  }
3162 
3163  /* check how stable density was */
3164  big_jump = exp10( big_jump );
3165  if( big_jump > 1.2 )
3166  {
3167  if( big_jump > 3. )
3168  {
3169  sprintf( chLine,
3170  " C-The H density jumped at by %.0f%%,"
3171  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
3172  big_jump*100.,
3173  struc.gas_phase[nzone_den_jump-1][ipHYDROGEN],
3174  nzone_den_jump,
3175  struc.gas_phase[nzone_den_jump][ipHYDROGEN],
3176  nzone_den_jump+1,
3177  struc.gas_phase[nzone_den_jump+1][ipHYDROGEN],
3178  nzone_den_jump+2 );
3179  caunin(chLine);
3180  }
3181  else
3182  {
3183  sprintf( chLine,
3184  " !An H density jump occurred by %.0f%%"
3185  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
3186  big_jump*100.,
3187  struc.gas_phase[nzone_den_jump-1][ipHYDROGEN],
3188  nzone_den_jump,
3189  struc.gas_phase[nzone_den_jump][ipHYDROGEN],
3190  nzone_den_jump+1,
3191  struc.gas_phase[nzone_den_jump+1][ipHYDROGEN],
3192  nzone_den_jump+2 );
3193  bangin(chLine);
3194  }
3195  }
3196 
3197  /* now do check on smoothness of radiation pressure */
3198  big_jump = -15;
3199  nzone_den_jump = 0;
3200 
3201  /* loop starts on zone 3 since dramatic fall in radiation pressure across first
3202  * few zones is normal behavior */
3203  /* >>chng 05 mar 12, add -lgAbort, since in some bad aborts current zone never evaluated */
3204  for( i=3; i < (nzone - 2 - (int)lgAbort); i++ )
3205  {
3206  /* this first check is on how the total hydrogen density has changed */
3207  rel = fabs(log10( SDIV(struc.pres_radiation_lines_curr[i]) /
3208  SDIV(0.5*(struc.pres_radiation_lines_curr[i-1]+struc.pres_radiation_lines_curr[i+1])) ) );
3209  /* only do significant stages of ionization */
3210  if( rel > big_jump )
3211  {
3212  big_jump = rel;
3213  nzone_den_jump = i;
3214  }
3215  }
3216  /* note that changing log big_jump to linear takes place in next branch */
3217 
3218  /* check how stable radiation pressure was, but only if significant */
3219  if( pressure.RadBetaMax > 0.01 )
3220  {
3221  big_jump = exp10( big_jump );
3222  if( big_jump > 1.2 )
3223  {
3224  /* only make it a caution is pressure jumped, and we were trying
3225  * to do a constant pressure model */
3226  if( big_jump > 3. && strcmp(dense.chDenseLaw,"CPRE") == 0)
3227  {
3228  sprintf( chLine,
3229  " C-The radiation pressure jumped by %.0f%%"
3230  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
3231  big_jump*100.,
3232  struc.pres_radiation_lines_curr[nzone_den_jump-1],
3233  nzone_den_jump,
3234  struc.pres_radiation_lines_curr[nzone_den_jump],
3235  nzone_den_jump+1,
3236  struc.pres_radiation_lines_curr[nzone_den_jump+1],
3237  nzone_den_jump+2 );
3238  caunin(chLine);
3239  }
3240  else
3241  {
3242  sprintf( chLine,
3243  " !The radiation pressure jumped by %.0f%%"
3244  " from %.2e (zone %ld) to %.2e (zone %ld) to %.2e (zone %ld)",
3245  big_jump*100.,
3246  struc.pres_radiation_lines_curr[nzone_den_jump-1],
3247  nzone_den_jump,
3248  struc.pres_radiation_lines_curr[nzone_den_jump],
3249  nzone_den_jump+1,
3250  struc.pres_radiation_lines_curr[nzone_den_jump+1],
3251  nzone_den_jump+2 );
3252  bangin(chLine);
3253  }
3254  }
3255  }
3256 
3257  /* these will be used to check on continuity */
3258  phycon.BigJumpTe = 0.;
3259  phycon.BigJumpne = 0.;
3260  phycon.BigJumpH2 = 0.;
3261  phycon.BigJumpCO = 0.;
3262 
3263  for( i=1; i < (nzone - 1 - (int)lgAbort); i++ )
3264  {
3265  /* check on how much temperature has changed */
3266  rel = fabs(log10( struc.testr[i] / struc.testr[i+1] ) );
3267  if( rel > phycon.BigJumpTe )
3268  {
3269  phycon.BigJumpTe = (realnum)rel;
3270  }
3271 
3272  /* check on how much electron density has changed */
3273  rel = fabs(log10( struc.ednstr[i] / struc.ednstr[i+1] ) );
3274  if( rel > phycon.BigJumpne )
3275  {
3276  phycon.BigJumpne = (realnum)rel;
3277  }
3278 
3279  /* check on how much H2 density has changed */
3280  if( (struc.H2_abund[i])>SMALLFLOAT && (struc.H2_abund[i+1]) > SMALLFLOAT
3281  /* only do this test if H2 abund is significant */
3282  && (struc.H2_abund[i])/struc.gas_phase[i][ipHYDROGEN]>1e-3)
3283  {
3284  rel = fabs(log10( (struc.H2_abund[i]) / SDIV(struc.H2_abund[i+1]) ) );
3285  if( rel > phycon.BigJumpH2 )
3286  {
3287  phycon.BigJumpH2 = (realnum)rel;
3288  }
3289  }
3290 
3291  int ipCO = findspecies("CO")->index;
3292  //fprintf(ioQQQ,"PRTCO %d %ld %d\n",ipCO,i,mole_global.num_calc);
3293  /* check on how much CO density has changed */
3294  if( ipCO != -1 &&
3295  struc.molecules[i][ipCO]>SMALLFLOAT &&
3296  struc.molecules[i+1][ipCO]>SMALLFLOAT &&
3297  struc.molecules[i][ipCO]/SDIV(struc.gas_phase[i][ipCARBON])>1e-3 )
3298  {
3299  rel = fabs(log10( struc.molecules[i][ipCO] / struc.molecules[i+1][ipCO] ) );
3300  if( rel > phycon.BigJumpCO )
3301  {
3302  phycon.BigJumpCO = (realnum)rel;
3303  }
3304  }
3305  }
3306 
3307  /* convert to linear change - subtract 1 to make it the residual difference */
3308  if( phycon.BigJumpTe > 0. )
3309  phycon.BigJumpTe = exp10( phycon.BigJumpTe ) - 1.f;
3310 
3311  if( phycon.BigJumpne > 0. )
3312  phycon.BigJumpne = exp10( phycon.BigJumpne ) - 1.f;
3313 
3314  if( phycon.BigJumpH2 > 0. )
3315  phycon.BigJumpH2 = exp10( phycon.BigJumpH2 ) - 1.f;
3316 
3317  if( phycon.BigJumpCO > 0. )
3318  phycon.BigJumpCO = exp10( phycon.BigJumpCO ) - 1.f;
3319  /*fprintf(ioQQQ,"DEBUG continuity large change %.2e %.2e %.2e %.2e \n",
3320  phycon.BigJumpTe , phycon.BigJumpne , phycon.BigJumpH2 , phycon.BigJumpCO );*/
3321 
3322  if( phycon.BigJumpTe > 0.3 )
3323  {
3324  sprintf( chLine,
3325  " C-The temperature varied by %.1f%% between two zones",
3326  phycon.BigJumpTe*100.);
3327  caunin(chLine);
3328  }
3329  else if( phycon.BigJumpTe > 0.1 )
3330  {
3331  sprintf( chLine,
3332  " !The temperature varied by %.1f%% between two zones",
3333  phycon.BigJumpTe*100.);
3334  bangin(chLine);
3335  }
3336 
3337  if( phycon.BigJumpne > 0.3 )
3338  {
3339  sprintf( chLine,
3340  " C-The electron density varied by %.1f%% between two zones",
3341  phycon.BigJumpne*100.);
3342  caunin(chLine);
3343  }
3344  else if( phycon.BigJumpne > 0.1 )
3345  {
3346  sprintf( chLine,
3347  " !The electron density varied by %.1f%% between two zones",
3348  phycon.BigJumpne*100.);
3349  bangin(chLine);
3350  }
3351 
3352  if( phycon.BigJumpH2 > 0.8 )
3353  {
3354  sprintf( chLine,
3355  " C-The H2 density varied by %.1f%% between two zones",
3356  phycon.BigJumpH2*100.);
3357  caunin(chLine);
3358  }
3359  else if( phycon.BigJumpH2 > 0.1 )
3360  {
3361  sprintf( chLine,
3362  " !The H2 density varied by %.1f%% between two zones",
3363  phycon.BigJumpH2*100.);
3364  bangin(chLine);
3365  }
3366 
3367  if( phycon.BigJumpCO > 0.8 )
3368  {
3369  sprintf( chLine,
3370  " C-The CO density varied by %.1f%% between two zones",
3371  phycon.BigJumpCO*100.);
3372  caunin(chLine);
3373  }
3374  else if( phycon.BigJumpCO > 0.2 )
3375  {
3376  sprintf( chLine,
3377  " !The CO density varied by %.1f%% between two zones",
3378  phycon.BigJumpCO*100.);
3379  bangin(chLine);
3380  }
3381 
3382  if( LineSave.lgIsoContSubSignif )
3383  {
3384  sprintf( chLine,
3385  " !Isotropic continuum subtraction significantly"
3386  " affects line intensities" );
3387  bangin(chLine);
3388  }
3389  return;
3390 }
t_phycon::TEnerDen
double TEnerDen
Definition: phycon.h:108
bangin
void bangin(const char *chLine)
Definition: warnings.h:90
t_conv::nTeFail
long int nTeFail
Definition: conv.h:201
t_iterations::lgIterAgain
bool lgIterAgain
Definition: iterations.h:53
t_hmi::h2line_cool_frac
realnum h2line_cool_frac
Definition: hmi.h:57
fudgec
t_fudgec fudgec
Definition: fudgec.cpp:5
mole_global
t_mole_global mole_global
Definition: mole.cpp:7
t_phycon::BigJumpTe
realnum BigJumpTe
Definition: phycon.h:116
t_hmi::h2dtot
realnum h2dtot
Definition: hmi.h:57
cdNotes
void cdNotes(FILE *ioOUT)
Definition: cddrive.cpp:310
t_thermal::ipMaxExtra
long int ipMaxExtra
Definition: thermal.h:163
t_radius::Radius
double Radius
Definition: radius.h:31
t_hydro::cintot
double cintot
Definition: hydrogenic.h:129
t_radius::depth
double depth
Definition: radius.h:31
atmdat
t_atmdat atmdat
Definition: atmdat.cpp:6
t_atmdat::HCharHeatMax
double HCharHeatMax
Definition: atmdat.h:300
t_save::lgPunContinuum
bool lgPunContinuum
Definition: save.h:369
co
t_co co
Definition: co.cpp:5
t_hcmap::lgMapOK
bool lgMapOK
Definition: hcmap.h:30
t_dense::nzEdenBad
long int nzEdenBad
Definition: dense.h:211
thermal
t_thermal thermal
Definition: thermal.cpp:6
diatomics::renorm_min
double renorm_min
Definition: h2_priv.h:343
t_rfield::lgGamrOK
bool lgGamrOK
Definition: rfield.h:442
t_thermal::power
double power
Definition: thermal.h:169
colden
t_colden colden
Definition: colden.cpp:5
t_timesc::TimeErode
realnum TimeErode
Definition: timesc.h:61
two_photon
Definition: two_photon.h:9
t_dense::EdenMax
double EdenMax
Definition: dense.h:204
t_conv::EdenErrorAllowed
double EdenErrorAllowed
Definition: conv.h:263
exp10
double exp10(double x)
Definition: cddefines.h:1368
t_struc::dr_ionfrac_limit
realnum dr_ionfrac_limit
Definition: struc.h:84
ipHE_LIKE
const int ipHE_LIKE
Definition: iso.h:65
t_rfield::comtot
double comtot
Definition: rfield.h:275
chLineLbl
string chLineLbl(const TransitionProxy &t)
Definition: transition.h:599
t_hmi::hmitot
double hmitot
Definition: hmi.h:33
GrainVar::dclmax
realnum dclmax
Definition: grainvar.h:559
TotalInsanity
NORETURN void TotalInsanity(void)
Definition: service.cpp:971
t_hydro::pop2mx
realnum pop2mx
Definition: hydrogenic.h:89
input
t_input input
Definition: input.cpp:12
Wind::lgWindOK
bool lgWindOK
Definition: wind.h:42
t_thermal::chCoolHeatMax
char chCoolHeatMax[NCOLNT_LAB_LEN+1]
Definition: thermal.h:127
opac
t_opac opac
Definition: opacity.cpp:5
t_thermal::GBarMax
realnum GBarMax
Definition: thermal.h:162
two_photon.h
t_mole_global::num_calc
int num_calc
Definition: mole.h:362
struc
t_struc struc
Definition: struc.cpp:6
hyperfine
t_hyperfine hyperfine
Definition: hyperfine.cpp:5
t_rfield::occ1nu
realnum occ1nu
Definition: rfield.h:455
t_geometry::lgZoneTrp
bool lgZoneTrp
Definition: geometry.h:97
t_thermal::FreeFreeTotHeat
double FreeFreeTotHeat
Definition: thermal.h:178
SMALLFLOAT
const realnum SMALLFLOAT
Definition: cpu.h:246
iso_ctrl
t_isoCTRL iso_ctrl
Definition: iso.cpp:9
t_cpu::i
t_cpu_i & i()
Definition: cpu.h:419
t_rfield::EnergyBremsThin
realnum EnergyBremsThin
Definition: rfield.h:227
NISO
const int NISO
Definition: cddefines.h:311
t_co::codtot
realnum codtot
Definition: co.h:22
t_pressure::ipPradMax_nzone
long int ipPradMax_nzone
Definition: pressure.h:106
t_radius::lgCylnOn
bool lgCylnOn
Definition: radius.h:127
t_conv::chNotConverged
char chNotConverged[INPUT_LINE_LENGTH]
Definition: conv.h:128
t_fudgec::lgFudgeUsed
bool lgFudgeUsed
Definition: fudgec.h:19
chkCaHeps
STATIC void chkCaHeps(double *totwid)
Definition: prt_comment.cpp:2927
cdCautions
void cdCautions(FILE *ioOUT)
Definition: cddrive.cpp:219
ipOXYGEN
const int ipOXYGEN
Definition: cddefines.h:356
t_rfield::lgOcc1Hi
bool lgOcc1Hi
Definition: rfield.h:463
t_timesc::tcmptn
double tcmptn
Definition: timesc.h:26
magnetic
t_magnetic magnetic
Definition: magnetic.cpp:17
t_geometry::iEmissPower
long int iEmissPower
Definition: geometry.h:71
two_photon::induc_dn_max
double induc_dn_max
Definition: two_photon.h:55
t_dynamics::lgTimeDependentStatic
bool lgTimeDependentStatic
Definition: dynamics.h:102
t_phycon::BigJumpH2
realnum BigJumpH2
Definition: phycon.h:116
EmissionProxy::TauTot
realnum & TauTot() const
Definition: emission.h:490
warnings
t_warnings warnings
Definition: warnings.cpp:11
cdWarnings
void cdWarnings(FILE *ioPNT)
Definition: cddrive.cpp:192
StopCalc
t_StopCalc StopCalc
Definition: stopcalc.cpp:7
conv
t_conv conv
Definition: conv.cpp:5
t_struc::ednstr
realnum * ednstr
Definition: struc.h:25
t_mole_global::lgH2Ozer
bool lgH2Ozer
Definition: mole.h:331
hextra
t_hextra hextra
Definition: hextra.cpp:5
t_hydro::HCollIonMax
realnum HCollIonMax
Definition: hydrogenic.h:123
HFLines
TransitionList HFLines("HFLines",&AnonStates)
t_atmdat::HCharCoolMax
double HCharCoolMax
Definition: atmdat.h:300
t_struc::molecules
realnum ** molecules
Definition: struc.h:71
phycon
t_phycon phycon
Definition: phycon.cpp:6
LineSave
t_LineSave LineSave
Definition: lines.cpp:10
t_hydro::nLyaHot
long int nLyaHot
Definition: hydrogenic.h:106
module_list::m_l
vector< module * > m_l
Definition: module.h:14
t_DoppVel::TurbVel
realnum TurbVel
Definition: doppvel.h:21
t_opac::stimax
realnum stimax[2]
Definition: opacity.h:317
t_colden::colden_old
realnum colden_old[NCOLD]
Definition: colden.h:32
t_pressure::lgPradDen
bool lgPradDen
Definition: pressure.h:113
t_oxy::poiii2Max
realnum poiii2Max
Definition: oxy.h:19
t_colden::ajmmin
realnum ajmmin
Definition: colden.h:59
t_thermal::lgTemperatureConstant
bool lgTemperatureConstant
Definition: thermal.h:44
t_atmdat::lgHCaseBOK
bool lgHCaseBOK[2][HS_NZ]
Definition: atmdat.h:342
t_dense::lgSetIoniz
bool lgSetIoniz[LIMELM]
Definition: dense.h:167
t_geometry::lgZoneSet
bool lgZoneSet
Definition: geometry.h:94
ioQQQ
FILE * ioQQQ
Definition: cddefines.cpp:7
findspecieslocal
molezone * findspecieslocal(const char buf[])
Definition: mole_species.cpp:831
t_geometry::FillFac
realnum FillFac
Definition: geometry.h:29
t_input::chTitle
char chTitle[INPUT_LINE_LENGTH]
Definition: input.h:80
t_oxy::r4363Max
realnum r4363Max
Definition: oxy.h:19
nzone
long int nzone
Definition: cddefines.cpp:14
t_called::lgTalk
bool lgTalk
Definition: called.h:12
t_warnings::chRgcln
char chRgcln[2][INPUT_LINE_LENGTH]
Definition: warnings.h:34
DoppVel
t_DoppVel DoppVel
Definition: doppvel.cpp:5
t_colden::rjnmin
realnum rjnmin
Definition: colden.h:59
dynamics
t_dynamics dynamics
Definition: dynamics.cpp:42
t_atmdat::HIonFracMax
double HIonFracMax
Definition: atmdat.h:317
t_phycon::BigJumpne
realnum BigJumpne
Definition: phycon.h:116
t_iso_sp::fb
vector< freeBound > fb
Definition: iso.h:481
TauLine2
TransitionList TauLine2("TauLine2",&AnonStates)
MIN2
#define MIN2(a, b)
Definition: cddefines.h:803
t_mesh::anu
double anu(size_t i) const
Definition: mesh.h:120
t_LineSave::lines
vector< LinSv > lines
Definition: lines.h:132
t_radius::lgDrMinUsed
bool lgDrMinUsed
Definition: radius.h:186
t_secondaries::SecHIonMax
realnum SecHIonMax
Definition: secondaries.h:42
ca
t_ca ca
Definition: ca.cpp:5
t_phycon::lgPhysOK
bool lgPhysOK
Definition: phycon.h:111
dense
t_dense dense
Definition: global.cpp:15
t_pressure::ipPradMax_line
long int ipPradMax_line
Definition: pressure.h:103
t_opac::lgOpacNeg
bool lgOpacNeg
Definition: opacity.h:192
GrainVar::GrnElecDonateMax
realnum GrnElecDonateMax
Definition: grainvar.h:531
Singleton< t_version >::Inst
static t_version & Inst()
Definition: cddefines.h:209
t_rfield::lgHabing
bool lgHabing
Definition: rfield.h:361
elementnames
t_elementnames elementnames
Definition: elementnames.cpp:5
t_conv::lgBadStop
bool lgBadStop
Definition: conv.h:246
t_LineSave::lgIsoContSubSignif
bool lgIsoContSubSignif
Definition: lines.h:150
iso_sp
t_iso_sp iso_sp[NISO][LIMELM]
Definition: iso.cpp:11
warnin
void warnin(const char *chLine)
Definition: warnings.h:74
t_ca::Ca2RmLya
realnum Ca2RmLya
Definition: ca.h:20
t_hydro::FracInd
realnum FracInd
Definition: hydrogenic.h:138
GrainVar::lgNegGrnDrg
bool lgNegGrnDrg
Definition: grainvar.h:486
t_isoCTRL::SmallA
realnum SmallA
Definition: iso.h:391
t_continuum::lgCoStarInterpolationCaution
bool lgCoStarInterpolationCaution
Definition: continuum.h:67
AgeCheck
void AgeCheck(void)
Definition: age_check.cpp:13
t_rt::lgMaserSetDR
bool lgMaserSetDR
Definition: rt.h:201
t_radius::lgDR2Big
bool lgDR2Big
Definition: radius.h:174
wind
Wind wind
Definition: wind.cpp:5
t_geometry::lgSphere
bool lgSphere
Definition: geometry.h:34
t_iso_sp::n_HighestResolved_local
long int n_HighestResolved_local
Definition: iso.h:538
iteration
long int iteration
Definition: cddefines.cpp:16
t_abund::lgAbnSolar
bool lgAbnSolar
Definition: abund.h:202
trace
t_trace trace
Definition: trace.cpp:5
lines_service.h
t_thermal::nUnstable
long int nUnstable
Definition: thermal.h:64
Wind::lgBallistic
bool lgBallistic(void) const
Definition: wind.h:31
t_radius::rinner
double rinner
Definition: radius.h:31
abund
t_abund abund
Definition: abund.cpp:5
t_iterations::lgOpticalDepthonverged
bool lgOpticalDepthonverged
Definition: iterations.h:57
geometry
t_geometry geometry
Definition: geometry.cpp:5
t_atmdat::lgdBaseSourceExists
bool lgdBaseSourceExists[LIMELM][LIMELM+1]
Definition: atmdat.h:452
t_phycon::TEMP_STOP_DEFAULT
const double TEMP_STOP_DEFAULT
Definition: phycon.h:119
t_iso_sp::TwoNu
vector< two_photon > TwoNu
Definition: iso.h:598
t_hydro::ndclev
long int ndclev
Definition: hydrogenic.h:139
t_he::nzone
long nzone
Definition: he.h:29
t_magnetic::lgB
bool lgB
Definition: magnetic.h:35
t_hcmap::lgMapDone
bool lgMapDone
Definition: hcmap.h:36
t_oxy::r5007Max
realnum r5007Max
Definition: oxy.h:19
t_rfield::tbrmax
realnum tbrmax
Definition: rfield.h:455
t_rfield::lgHPhtOK
bool lgHPhtOK
Definition: rfield.h:442
ipH1s
const int ipH1s
Definition: iso.h:29
PrtComment
void PrtComment(void)
Definition: prt_comment.cpp:66
STATIC
#define STATIC
Definition: cddefines.h:118
t_trace::lgTrace
bool lgTrace
Definition: trace.h:12
t_he::frac_he0dest_23S_photo
double frac_he0dest_23S_photo
Definition: he.h:25
lgConserveEnergy
bool lgConserveEnergy()
Definition: energy.cpp:311
t_rfield::lgCMB_set
bool lgCMB_set
Definition: rfield.h:106
TransitionProxy::Emis
EmissionList::reference Emis() const
Definition: transition.h:447
t_conv::ifailz
long int ifailz[12]
Definition: conv.h:234
continuum
t_continuum continuum
Definition: continuum.cpp:6
mole
t_mole_local mole
Definition: mole.cpp:8
t_phycon::BigJumpCO
realnum BigJumpCO
Definition: phycon.h:116
findspecies
molecule * findspecies(const char buf[])
Definition: mole_species.cpp:789
pressure
t_pressure pressure
Definition: pressure.cpp:9
rfield
t_rfield rfield
Definition: rfield.cpp:9
t_thermal::HeatLineMax
realnum HeatLineMax
Definition: thermal.h:181
atoms
t_atoms atoms
Definition: atoms.cpp:5
t_opac::lgCaseB
bool lgCaseB
Definition: opacity.h:174
FndLineHt
const TransitionProxy FndLineHt(long int *level)
Definition: lines_service.cpp:750
t_rfield::ConInterOut
realnum * ConInterOut
Definition: rfield.h:154
realnum
float realnum
Definition: cddefines.h:124
t_mole_local::species
valarray< class molezone > species
Definition: mole.h:468
t_rfield::uh
realnum uh
Definition: rfield.h:349
EXIT_FAILURE
#define EXIT_FAILURE
Definition: cddefines.h:168
INPUT_LINE_LENGTH
const int INPUT_LINE_LENGTH
Definition: cddefines.h:301
t_geometry::lgStatic
bool lgStatic
Definition: geometry.h:64
t_opac::thmin
realnum thmin
Definition: opacity.h:188
t_radius::lgdR2Small
bool lgdR2Small
Definition: radius.h:118
t_thermal::thist
realnum thist
Definition: thermal.h:68
t_opac::telec
realnum telec
Definition: opacity.h:188
hydro
t_hydro hydro
Definition: hydrogenic.cpp:5
notein
void notein(const char *chLine)
Definition: warnings.h:82
cdEXIT
#define cdEXIT(FAIL)
Definition: cddefines.h:482
t_pressure::RadBetaMax
realnum RadBetaMax
Definition: pressure.h:96
molecule::index
int index
Definition: mole.h:194
powi
double powi(double, long int)
Definition: service.cpp:594
t_thermal::totcol
double totcol
Definition: thermal.h:130
t_geometry::covrt
realnum covrt
Definition: geometry.h:61
h2
diatomics h2("h2", 4100.,&hmi.H2_total, Yan_H2_CS)
t_hydro::nbul
long int nbul
Definition: hydrogenic.h:141
wlAirVac
realnum wlAirVac(double wlAir)
Definition: lines_service.cpp:135
t_pressure::chLineRadPres
string chLineRadPres
Definition: pressure.h:109
iterations
t_iterations iterations
Definition: iterations.cpp:6
t_pressure::lgPradCap
bool lgPradCap
Definition: pressure.h:113
column
double column(const genericState &gs)
Definition: generic_state.cpp:15
nWindLine
long nWindLine
Definition: cdinit.cpp:19
t_atoms::xMg2Max
realnum xMg2Max
Definition: atoms.h:129
t_rfield::lgXRayOK
bool lgXRayOK
Definition: rfield.h:442
t_elementnames::chElementNameShort
char chElementNameShort[LIMELM][CHARS_ELEMENT_NAME_SHORT]
Definition: elementnames.h:21
t_struc::H2_abund
realnum * H2_abund
Definition: struc.h:71
t_hextra::cryden
realnum cryden
Definition: hextra.h:24
radius
t_radius radius
Definition: radius.cpp:5
timesc
t_timesc timesc
Definition: timesc.cpp:7
t_iso_sp::QuantumNumbers2Index
multi_arr< long, 3 > QuantumNumbers2Index
Definition: iso.h:490
t_StopCalc::TempLoStopZone
realnum TempLoStopZone
Definition: stopcalc.h:42
t_isoCTRL::lgCritDensLMix
bool lgCritDensLMix[NISO]
Definition: iso.h:422
t_dense::lgElmtOn
bool lgElmtOn[LIMELM]
Definition: dense.h:160
t_iso_sp::trans
TransitionProxy trans(const long ipHi, const long ipLo)
Definition: iso.h:473
t_elementnames::chElementSym
char chElementSym[LIMELM][CHARS_ELEMENT_SYM]
Definition: elementnames.h:25
t_dynamics::CoolMax
double CoolMax
Definition: dynamics.h:74
t_hydro::TeLyaMax
realnum TeLyaMax
Definition: hydrogenic.h:109
t_fudgec::nfudge
long int nfudge
Definition: fudgec.h:17
t_struc::nzlim
long int nzlim
Definition: struc.h:19
t_dense::gas_phase
realnum gas_phase[LIMELM]
Definition: dense.h:76
t_hmi::HeatH2DexcMax
realnum HeatH2DexcMax
Definition: hmi.h:57
t_iterations::itermx
long int itermx
Definition: iterations.h:37
ipH2p
const int ipH2p
Definition: iso.h:31
t_cpu_i::foundMD5Mismatch
bool foundMD5Mismatch() const
Definition: cpu.h:402
t_input::lgUnderscoreFound
bool lgUnderscoreFound
Definition: input.h:112
t_input::lgBracketFound
bool lgBracketFound
Definition: input.h:116
t_geometry::lgStaticNoIt
bool lgStaticNoIt
Definition: geometry.h:84
ASSERT
#define ASSERT(exp)
Definition: cddefines.h:613
t_geometry::covaper
realnum covaper
Definition: geometry.h:54
t_iterations::lgLastIt
bool lgLastIt
Definition: iterations.h:47
ioPrnErr
FILE * ioPrnErr
Definition: cddefines.cpp:9
t_dense::chDenseLaw
char chDenseLaw[5]
Definition: dense.h:176
ipH2s
const int ipH2s
Definition: iso.h:30
t_hydro::nZTLaMax
long int nZTLaMax
Definition: hydrogenic.h:114
ipCALCIUM
const int ipCALCIUM
Definition: cddefines.h:368
cdLine
long int cdLine(const char *chLabel, realnum wavelength, double *relint, double *absint)
Definition: cddrive.cpp:1067
t_thermal::ConstTemp
realnum ConstTemp
Definition: thermal.h:56
t_continuum::lgCon0
bool lgCon0
Definition: continuum.h:67
he
t_he he
Definition: he.cpp:5
ipH_LIKE
const int ipH_LIKE
Definition: iso.h:64
Wind::emdot
realnum emdot
Definition: wind.h:39
LIMELM
const int LIMELM
Definition: cddefines.h:308
diatomics::renorm_max
double renorm_max
Definition: h2_priv.h:343
t_atmdat::ipSpecIon
vector< vector< long > > ipSpecIon
Definition: atmdat.h:455
DEBUG_ENTRY
#define DEBUG_ENTRY(funcname)
Definition: cddefines.h:723
ipHELIUM
const int ipHELIUM
Definition: cddefines.h:350
t_abund::lgDepln
bool lgDepln
Definition: abund.h:251
t_conv::lgAutoIt
bool lgAutoIt
Definition: conv.h:249
t_rfield::lgHionRad
bool lgHionRad
Definition: rfield.h:452
t_struc::gas_phase
realnum ** gas_phase
Definition: struc.h:75
GrainVar::bin
vector< GrainBin * > bin
Definition: grainvar.h:583
t_conv::failmx
realnum failmx
Definition: conv.h:204
t_hyperfine::cooling_max
realnum cooling_max
Definition: hyperfine.h:65
t_rfield::lgMMok
bool lgMMok
Definition: rfield.h:442
t_he::frac_he0dest_23S
double frac_he0dest_23S
Definition: he.h:25
t_StopCalc::chReasonStop
char chReasonStop[nCHREASONSTOP]
Definition: stopcalc.h:130
t_input::lgSetNoBuffering
bool lgSetNoBuffering
Definition: input.h:122
t_opac::TauAbsGeo
realnum ** TauAbsGeo
Definition: opacity.h:91
lgTestCodeCalled
bool lgTestCodeCalled
Definition: cddefines.cpp:11
oxy
t_oxy oxy
Definition: oxy.cpp:5
cdSurprises
void cdSurprises(FILE *ioOUT)
Definition: cddrive.cpp:274
cdNwcns
void cdNwcns(bool *lgAbort_ret, long int *NumberWarnings, long int *NumberCautions, long int *NumberNotes, long int *NumberSurprises, long int *NumberTempFailures, long int *NumberPresFailures, long int *NumberIonFailures, long int *NumberNeFailures)
Definition: cddrive.cpp:1173
t_hydro::lgHiPop2
bool lgHiPop2
Definition: hydrogenic.h:88
MAX2
#define MAX2(a, b)
Definition: cddefines.h:824
fprintf
int fprintf(const Output &stream, const char *format,...)
Definition: service.cpp:1121
t_dense::lgEdenBad
bool lgEdenBad
Definition: dense.h:238
GrainVar::GrnElecHoldMax
realnum GrnElecHoldMax
Definition: grainvar.h:532
t_mole_global::list
MoleculeList list
Definition: mole.h:365
t_dynamics::lgRecom
bool lgRecom
Definition: dynamics.h:108
cdReasonGeo
void cdReasonGeo(FILE *ioOUT)
Definition: cddrive.cpp:170
t_co::codfrc
realnum codfrc
Definition: co.h:22
SDIV
sys_float SDIV(sys_float x)
Definition: cddefines.h:1002
lgPrnErr
bool lgPrnErr
Definition: cddefines.cpp:13
ipCARBON
const int ipCARBON
Definition: cddefines.h:354
t_rt::DoubleTau
realnum DoubleTau
Definition: rt.h:178
t_warnings::zero
void zero(void)
Definition: warnings.cpp:13
hcmap
t_hcmap hcmap
Definition: hcmap.cpp:23
Wind::lgStatic
bool lgStatic(void) const
Definition: wind.h:24
t_struc::testr
realnum * testr
Definition: struc.h:25
t_elementnames::chElementName
char chElementName[LIMELM][CHARS_ELEMENT_NAME]
Definition: elementnames.h:17
t_rfield::occmnu
realnum occmnu
Definition: rfield.h:455
t_thermal::lgEdnGTcm
bool lgEdnGTcm
Definition: thermal.h:77
gv
GrainVar gv
Definition: grainvar.cpp:5
hmi
t_hmi hmi
Definition: hmi.cpp:5
secondaries
t_secondaries secondaries
Definition: secondaries.cpp:5
t_radius::r1r0sq
double r1r0sq
Definition: radius.h:31
t_thermal::wlCoolHeatMax
realnum wlCoolHeatMax
Definition: thermal.h:126
cpu
static t_cpu cpu
Definition: cpu.h:427
prt_smooth_predictions
STATIC void prt_smooth_predictions(void)
Definition: prt_comment.cpp:2993
t_warnings::lgCautns
bool lgCautns
Definition: warnings.h:44
t_oxy::poiii3Max
realnum poiii3Max
Definition: oxy.h:19
ShowMe
void ShowMe(void)
Definition: service.cpp:205
dBaseTrans
vector< TransitionList > dBaseTrans
Definition: taulines.cpp:18
GrainVar::TotalDustHeat
realnum TotalDustHeat
Definition: grainvar.h:559
save
t_save save
Definition: save.cpp:5
t_phycon::te
double te
Definition: phycon.h:21
t_rfield::plsfrqmax
realnum plsfrqmax
Definition: rfield.h:430
GrainVar::dphmax
realnum dphmax
Definition: grainvar.h:559
t_timesc::time_therm_long
double time_therm_long
Definition: timesc.h:29
GrainVar::lgDustOn
bool lgDustOn() const
Definition: grainvar.h:475
t_thermal::CoolHeatMax
realnum CoolHeatMax
Definition: thermal.h:125
ipHYDROGEN
const int ipHYDROGEN
Definition: cddefines.h:349
t_rfield::lgComUndr
bool lgComUndr
Definition: rfield.h:282
t_conv::BigEdenError
realnum BigEdenError
Definition: conv.h:213
t_rfield::nflux
long int nflux
Definition: rfield.h:46
t_hydro::fbul
realnum fbul
Definition: hydrogenic.h:140
t_hmi::h2dfrc
realnum h2dfrc
Definition: hmi.h:57
t_colden::colden
realnum colden[NCOLD]
Definition: colden.h:32
DEPTH_OFFSET
const double DEPTH_OFFSET
Definition: cddefines.h:322
t_rfield::lgPlasNu
bool lgPlasNu
Definition: rfield.h:428
t_dynamics::HeatMax
double HeatMax
Definition: dynamics.h:74
t_struc::xIonDense
realnum *** xIonDense
Definition: struc.h:64
t_hydro::TLyaMax
realnum TLyaMax
Definition: hydrogenic.h:109
t_oxy::poimax
realnum poimax
Definition: oxy.h:19
t_rt::lgMaserCapHit
bool lgMaserCapHit
Definition: rt.h:205
ipH3p
const int ipH3p
Definition: iso.h:33
t_pressure::lgPres_radiation_ON
bool lgPres_radiation_ON
Definition: pressure.h:90
t_iso_sp::numLevels_local
long int numLevels_local
Definition: iso.h:529
t_iso_sp::nCollapsed_local
long int nCollapsed_local
Definition: iso.h:519
t_rfield::tbr4nu
realnum tbr4nu
Definition: rfield.h:455
t_hmi::h2pmax
realnum h2pmax
Definition: hmi.h:131
Wind::windv
realnum windv
Definition: wind.h:18
EmissionProxy::TauIn
realnum & TauIn() const
Definition: emission.h:470
t_warnings::lgWarngs
bool lgWarngs
Definition: warnings.h:44
called
t_called called
Definition: called.cpp:4
t_geometry::filpow
realnum filpow
Definition: geometry.h:29
lgAbort
bool lgAbort
Definition: cddefines.cpp:10
caunin
void caunin(const char *chLine)
Definition: warnings.h:98
t_rfield::occmax
realnum occmax
Definition: rfield.h:455
t_rfield::cinrat
double cinrat
Definition: rfield.h:275
rt
t_rt rt
Definition: rt.cpp:5
t_struc::pres_radiation_lines_curr
realnum * pres_radiation_lines_curr
Definition: struc.h:25
t_hmi::CoolH2DexcMax
realnum CoolH2DexcMax
Definition: hmi.h:57
t_rfield::tbrmnu
realnum tbrmnu
Definition: rfield.h:455
Generated by   doxygen 1.8.5