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