cloudy  trunk
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
parse_commands.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 /*ParseCommands main command line parser, decode command, then call other routines to read */
4 #include "cddefines.h"
5 #include "parse.h"
6 #include "stopcalc.h"
7 #include "abund.h"
8 #include "geometry.h"
9 #include "dense.h"
10 #include "plot.h"
11 #include "grid.h"
12 #include "grainvar.h"
13 #include "dynamics.h"
14 #include "magnetic.h"
15 #include "trace.h"
16 #include "atmdat.h"
17 #include "h2.h"
18 #include "rt.h"
19 #include "thermal.h"
20 #include "opacity.h"
21 #include "called.h"
22 #include "wind.h"
23 #include "hextra.h"
24 #include "iterations.h"
25 #include "radius.h"
26 #include "input.h"
27 #include "monitor_results.h"
28 #include "phycon.h"
29 #include "fudgec.h"
30 #include "version.h"
31 #include "conv.h"
32 #include "cosmology.h"
33 #include "pressure.h"
34 #include "parser.h"
35 #include "dark_matter.h"
36 #include "iso.h"
37 #include "mole.h"
38 #include "parse_species.h"
39 #include "doppvel.h"
40 #include "rfield.h"
41 #include "prt.h"
42 
43 void ParseAperture(Parser &p);
45 void ParseCExtra(Parser &p);
46 void ParseCMBOuter(Parser &p);
47 void ParseCosm(Parser &p);
48 void ParseCovering(Parser &p);
49 void ParseCylinder(Parser &p);
50 void ParseDarkMatter(Parser &p);
51 void ParseDatabase(Parser &p);
53 void ParseDiffuse(Parser &p);
54 void ParseDistance(Parser &p);
55 void ParseDoubleTau(Parser &);
56 void ParseEden(Parser &p);
57 void ParseEnergy(Parser &p);
58 void ParseFail(Parser &p);
59 void ParseFill(Parser &p);
60 void ParseF_nuSpecific(Parser &p);
62 void ParseFudge(Parser &p);
63 void ParsePGrains(Parser &);
64 void ParseGravity(Parser &p);
65 void ParseHeLike(Parser &);
66 void ParseHelp(Parser &);
67 void ParseHExtra(Parser &p);
68 void ParseConvHighT(Parser &);
69 void ParseHydrogen(Parser &);
70 void ParseInitCount(Parser &p);
71 void ParseIntensity(Parser &p);
72 void ParseIterations(Parser &p);
73 void ParseL_nu(Parser &p);
74 void ParseLaser(Parser &p);
75 void ParseLuminosity(Parser &p);
76 void ParseNeutrons(Parser &p);
77 void ParseNuF_nu(Parser &p);
78 void ParseNuL_nu(Parser &p);
79 void ParsePhi(Parser &p);
80 void ParseQH(Parser &p);
81 void ParseRoberto(Parser &);
82 void ParseSpecial(Parser &);
83 void ParseTauMin(Parser &p);
84 void ParseTitle(Parser &);
85 void ParseTolerance(Parser &);
86 void ParseVLaw(Parser &p);
87 void ParseTurbulence(Parser &p);
88 
89 void ParseCommands(void)
90 {
91  bool lgStop ,
92  lgStop_not_enough_info;
93 
94  DEBUG_ENTRY( "ParseCommands()" );
95 
96  /* following says abundances are solar */
97  abund.lgAbnSolar = true;
98 
99  /* there are no plots desired yet */
100  plotCom.lgPlotON = false;
101  plotCom.nplot = 0;
102 
103  /* this flag remembers whether grains have ever been turned on. It is passed
104  * to routine ParseAbundances - there grains will be turned on with commands
105  * such as abundances ism, unless grains were previously set
106  * with a grains command. this way abundances will not clobber explicitly set
107  * grains. */
108 
109  radius.Radius = 0.;
110  radius.rinner = 0.;
111  rfield.nShape = 0;
112 
113  /* initialize some code variables in case assert command used in input stream */
115 
116  for( long int i=0; i < LIMSPC; i++ )
117  {
118  strcpy( rfield.chRSpec[i], "UNKN" );
119  }
120  optimize.nparm = 0;
121 
122  /* this is an option to turn on trace printout on the nth
123  * call from the optimizer */
124  if( optimize.lgTrOpt )
125  {
126  /* nTrOpt was set with "optimize trace" command,
127  * is iteration to turn on trace */
129  {
130  trace.lgTrace = true;
131  called.lgTalk = cpu.i().lgMPI_talk();
132  trace.lgTrOvrd = true;
133  fprintf( ioQQQ, " READR turns on trace from optimize option.\n" );
134  }
135  }
136 
137  if ( called.lgTalk && prt.lgPrintHTML )
138  {
139  fprintf( ioQQQ, "<!DOCTYPE html PUBLIC "
140  "\"-//W3C//DTD XHTML 1.0 Transitional//EN\" "
141  "\"http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd\">\n");
142  fprintf( ioQQQ,"<html xmlns=\"http://www.w3.org/1999/xhtml\">\n");
143  fprintf( ioQQQ,"<head>\n");
144  fprintf( ioQQQ,"<meta http-equiv=\"Content-Type\" "
145  "content=\"text/html; charset=UTF-8\" />\n");
146  fprintf( ioQQQ,"<title>Cloudy output</title>\n");
147  fprintf( ioQQQ,"</head>\n");
148  fprintf( ioQQQ,"<body>\n");
149  fprintf( ioQQQ,"<pre>\n");
150  }
151 
152  /* say this is a beta version if we are talking and it is the truth */
153  if( t_version::Inst().nBetaVer > 0 && called.lgTalk )
154  {
155  fprintf( ioQQQ,
156  "\n This is a beta release of Cloudy, and is intended for testing only.\n" );
157  fprintf( ioQQQ,
158  "Please help make Cloudy better by posing problems or suggestions on http://tech.groups.yahoo.com/group/cloudy_simulations/.\n\n" );
159  }
160 
161  if( called.lgTalk )
162  {
163  /* this code prints pretty lines at top of output box */
164  int indent = (int)((122 - strlen(t_version::Inst().chVersion))/2);
165  fprintf( ioQQQ, "%*cCloudy %s\n", indent, ' ', t_version::Inst().chVersion );
166 
167  fprintf( ioQQQ, "%57cwww.nublado.org\n\n", ' ' );
168 
169  /* now print box and date of version, before printing commands */
170  fprintf( ioQQQ, "%23c", ' ' );
171  fprintf( ioQQQ, "**************************************");
172  fprintf( ioQQQ, "%7.7s", t_version::Inst().chDate);
173  fprintf( ioQQQ, "**************************************\n");
174 
175  fprintf( ioQQQ, "%23c*%81c*\n", ' ', ' ' );
176  }
177 
178  /* read in commands and print them */
179 
180  /* following signals which read is in progress,
181  * 1 is forward read, of input commands
182  * -1 is reverse read from bottom of line array, of cloudy.ini file */
183  input.iReadWay = 1;
184 
185  /* initialize array reader, this sub does nothing but set
186  * initial value of a variable, depending on value of iReadWay */
187  input.init();
188 
189  static const CloudyCommand commands[] = {
190  {"ABSO",ParseAbsMag},
191  /* enter luminosity in absolute magnitudes, in reads2 */
192  {"AGE",ParseAge},
193  /* enter age of cloud so we can check for time-steady reads2 */
194  {"AGN ",ParseAgn},
195  /* enter generic style AGN continuum, in reads2 */
196  {"ABUN",ParseAbundances},
197  {"APER",ParseAperture},
198  {"BACK", ParseBackgrd},
199  /* cosmic background, in parse_backgrd */
200  {"BLAC", ParseBlackbody},
201  /* black body, in reads2 */
202  {"BREM", ParseBremsstrahlung},
203  {"CASE",ParseCaseB},
204  /* do Case A or Case B (usually Case B since most realistic) */
205  {"CEXT",ParseCExtra},
206  {"CHEM",ParseChemistry},
207  /* CMB command */
208  {"CLUM",ParseFill},
209  {"CMB", ParseCMBOuter},
210  /* cosmic thermal background radiation, argument is redshift */
211  /* if no number on line then (realnum)FFmtRead returns z=0; i.e., now */
212  {"COMP", ParseCompile},
213  /* compile ascii version of stellar atmosphere continua in volk */
214  {"CONS",ParseConstant},
215  /* constant temperature, pressure, density, or gas pressure
216  * in readsun */
217  {"CORO",ParseCoronal},
218  /* coronal equilibrium; set constant temperature to number on line
219  * in readsun */
220  {"COSMOLOGY",ParseCosmology},
221  {"COSMI", ParseCosmicRays},
222  {"COSM",ParseCosm}, // Backstop for ambiguity
223  {"COVE",ParseCovering},
224  {"CRAS",ParseCrashDo},
225  /* any of several tests to cause the code to crash */
226  {"CYLI",ParseCylinder},
227  {"DARK",ParseDarkMatter},
228  {"DATA",ParseDatabase},
229  {"ATOM",ParseDatabase},// accept the old ATOM command as a pseudonym for DATABASE
230  {"DIEL",ParseDielectronic},
231  {"DIFF",ParseDiffuse},
232  {"DIST",ParseDistance},
233  {"DLAW",ParseDLaw},
234  /* either use dense_fabden routine, or read in table of points */
235  {"DOUB",ParseDoubleTau},
236  /* double optical depth scale after each iteration */
237  {"DRIV",ParseDrive},
238  /* call one of several drivers, readsun */
239  {"EDEN",ParseEden},
240  /* option to add "extra" electrons, to test Compton limit
241  * for very low T(star) - option is log of eden */
242  {"ELEM",ParseElement},
243  /* element command;
244  * scale or abundance options, to change abundance of specific element
245  * read option to change order of elements
246  * in reads2.f */
247  {"ENER",ParseEnergy},
248  {"EXTI",ParseExtinguish},
249  /* extinguish ionizing continuum by absorbing column AFTER
250  * setting luminosity or Q(H). First number is the column
251  * density (log), second number is leakage (def=0%)
252  * last number is lowest energy (ryd), last two may be omitted
253  * from right to left
254  *
255  * extinction is actually done in extin, which is called by ContSetIntensity */
256  {"FAIL",ParseFail},
257  /* reset number of temp failures allowed, default=20 */
258  {"FILL",ParseFill},
259  {"FLUC",ParseFluc},
260  /* rapid density fluctuations, in readsun */
261  {"F(NU",ParseF_nuSpecific},
262  /* this is the specific flux at nu
263  * following says F(nu) not nuF(nu) */
264  {"FORC",ParseForceTemperature},
265  /* force temperature of first zone, but don't keep constant
266  * allow to then go to nearest equilibrium
267  * log if < 10 */
268  {"FUDG",ParseFudge},
269  {"GLOB",ParseGlobule},
270  /* globule with density increasing inward
271  * parameters are outer density, radius of globule, and density power */
272  {"GRAI",ParseGrain},
273  /* read parameters dealing with grains, in reads2 */
274  {"PGRA",ParsePGrains},
275  {"GRAV",ParseGravity},
276  /* (self-)Gravity forces: Yago Ascasibar (UAM, Spring 2009) */
277  {"GRID",ParseGrid},
278  /* option to run grid of models by varying certain parameters
279  * in reads2 */
280  {"HDEN",ParseHDEN},
281  /* parse the hden command to set the hydrogen density, in reads2 */
282  {"HELI",ParseHeLike},
283  {"HELP",ParseHelp},
284  {"HEXT",ParseHExtra},
285  /* "extra" heating rate, so that first= log(erg(cm-3, s-1),
286  * second optional number is scale radius, so that HXTOT = TurbHeat*SEXP(DEPTH/SCALE)
287  * if missing then constant heating.
288  * third option is depth from shielded face, to mimic irradiation from both sides*/
289  {"HIGH",ParseConvHighT},
290  /* approach equilibrium from high te */
291  {"HYDROGEN",ParseHydrogen},
292  {"ILLU",ParseIlluminate},
293  // illuminate command
294  {"INIT",ParseInitCount},
295  {"INTEN",ParseIntensity},
296  {"INTER",ParseInterp},
297  /* interpolate on input tables for continuum, set of power laws used
298  * input ordered pairs nu( ryd or log(Hz),, log(fnu)
299  * additional lines begin CONTINUE
300  * first check that this is the one and only INTERP command
301  * in readsun */
302  {"IONI",ParseIonParI},
303  /* inter ionization parameter U=Q/12 R*R N C;
304  * defined per hydrogen, not per electron (as before)
305  * radius must also be entered if spherical, not needed if plane */
306  {"ITER",ParseIterations},
307  {"L(NU",ParseL_nu},
308  {"LASE",ParseLaser},
309  {"LUMI",ParseLuminosity},
310  {"MAGN", ParseMagnet},
311  /* parse the magnetic field command, routine in magnetic.c */
312  {"MAP ",ParseMap},
313  /* do cooling space map for specified zones, in reads2 */
314  {"META",ParseMetal},
315  /* read depletion for metals, all elements heavier than He
316  * in reads2 */
317  {"MONI",ParseMonitorResults},
318  /* monitor that code predicts certain results, in monitor_results.h */
319  {"NEUT",ParseNeutrons},
320  {"NO ", ParseDont},
321  /* don't do something, in readsun */
322  {"NORM",ParseNorm},
323  /* normalize lines to this rather than h-b, sec number is scale factor */
324  {"NUF(",ParseNuF_nu},
325  {"NUL(",ParseNuL_nu},
326  {"OPTI",ParseOptimize},
327  /* option to optimize the model by varying certain parameters
328  * in reads2 */
329  {"PHI(",ParsePhi},
330  {"PLOT", ParsePlot},
331  /* make plot of log(nu.f(nu), vs log(nu) or opacity
332  * in file plot */
333  {"POWE", ParsePowerlawContinuum},
334  /* power law with cutoff, in reads2 */
335  {"PRIN",ParsePrint},
336  /* adjust print schedule, in readsun */
337  {"PUNC",ParseSave},
338  /* save something, in save */
339  {"Q(H)",ParseQH},
340  {"RATI",ParseRatio},
341  /* enter a continuum luminosity as a ratio of
342  * nuFnu for this continuum relative to a previous continuum
343  * format; first number is ratio of second to first continuum
344  * second number is energy for this ratio
345  * if third number on line, then 2nd number is energy of
346  * first continuum, while 3rd number is energy of second continuum
347  * in reads2 */
348  {"RADI", ParseRadius},
349  /* log of inner and outer radii, default second=infinity,
350  * if R2<R1 then R2=R1+R2
351  * there is an optional keyword, "PARSEC" on the line
352  * to use PC as units, reads2 */
353  {"ROBE",ParseRoberto},
354  {"SAVE",ParseSave},
355  {"SET ",ParseSet},
356  {"SPECIAL", ParseSpecial},
357  {"SPECIES",ParseSpecies},
358  {"SPHE", ParseSphere},
359  /* compute a spherical model, diffuse field from other side in
360  * in reads2 */
361  {"STAT",ParseState},
362  /* either get or put the code's state as a file on disk */
363  {"STOP",ParseStop},
364  /* stop model at desired zone, temperature, column density or tau-912
365  * in readsun */
366  {"TABL",ParseTable},
367  /* interpolate on input tables for continuum, set of power laws used
368  * input stored in big BLOCK data
369  * first check that this is the one and only INTERP command
370  * in readsun */
371  {"TAUM",ParseTauMin},
372  {"TEST",ParseTest},
373  /* parse the test command and its options */
374  {"TIME",ParseDynaTime},
375  /* parse the time dependent command, in dynamics.c */
376  {"TITL",ParseTitle},
377  {"TLAW", ParseTLaw},
378  /* some temperature vs depth law */
379  {"TOLE", ParseTolerance},
380  {"TRAC", ParseTrace},
381  /* turn on trace output, in reads2 */
382  {"VLAW",ParseVLaw},
383  {"TURB",ParseTurbulence},
384  {"WIND",ParseDynaWind},
385  /* NB - advection and wind commands are now a single command */
386  /* parse the wind command, in dynamics.c */
387  {"XI",ParseIonParX},
388  {NULL,NULL}}; // {NULL,NULL} sentinel must come last
389 
390  Parser p(commands);
391 
392  p.m_nqh = 0;
393  p.m_nInitFile=0;/* used to count number of init files read in */
394  p.m_lgDSet = false;
395  p.m_lgEOF = false;
396 
397  // set default solar abundances
398  char chDUMMY[INPUT_LINE_LENGTH];
399  sprintf(chDUMMY, "abundances \"default.abn\"" );
400  p.setline(chDUMMY);
401  ParseAbundances( p );
402 
403  // set default isotopic abundances
404  sprintf( chDUMMY, "abundances isotopes \"default-iso.abn\"" );
405  p.setline( chDUMMY );
406  ParseAbundances( p );
407 
408  /* read until eof or blank line, then return control back to main program */
409  while (p.getline())
410  {
411  /* line starting with blank is taken as end of commands */
412  if ( p.last( ) )
413  break;
414 
415  /* echo the line but only if it does not contain the keyword HIDE */
416  p.echo( );
417 
418  /* check whether VARY is on line */
419  if( p.nMatch("VARY") )
420  {
421  optimize.lgVarOn = true;
422  if( optimize.nparm + 1 > LIMPAR )
423  {
424  fprintf( ioQQQ, " Too many VARY lines entered; the limit is%4ld\n",
425  LIMPAR );
427  }
428  }
429  else
430  {
431  optimize.lgVarOn = false;
432  }
433 
434  if( p.isCommandComment() )
435  {
436  ((void)0); // do nothing for comments
437  }
438  else if( p.isVar() )
439  {
440  p.doSetVar();
441  }
442  else
443  {
444  long int i;
445  for (i=0; commands[i].name != NULL; ++i)
446  {
447  if (p.Command(commands[i].name,commands[i].action))
448  break;
449  }
450  if (commands[i].name == NULL)
451  {
452  p.CommandError(); // No command was recognized
453  }
454  }
455  }
456 
457  /*------------------------------------------------------------------- */
458  /* fall through - hit lgEOF or blank line */
459 
460  if( called.lgTalk )
461  {
462  fprintf( ioQQQ, "%23c*%81c*\n", ' ', ' ' );
463  fprintf( ioQQQ, "%23c***********************************************************************************\n\n\n\n", ' ' );
464  }
465 
466  /* set R to large value if U specified but R is not */
467  /* set radius to very large value if not already set */
468  /* >>chng 01 jul 24, from Radius == 0 to this, as per PvH comments */
469  if( !radius.lgRadiusKnown )
470  {
472  }
473 
474  /* this is an option to turn on trace printout on the nth
475  * call from the optimizer - only allow trace if
476  * this is the case and nOptimiz 1 below nTrOpt */
477  if( optimize.lgTrOpt )
478  {
479  /* nTrOpt was set with "optimize trace" command,
480  * is iteration to turn on trace */
482  {
483  trace.lgTrace = false;
484  /* following overrides turning on trace elsewhere */
485  trace.lgTrOvrd = false;
486  }
487  else
488  {
489  trace.lgTrace = true;
490  called.lgTalk = cpu.i().lgMPI_talk();
491  trace.lgTrOvrd = true;
492  fprintf( ioQQQ, " READR turns on trace from optimize option.\n" );
493  }
494  }
495 
496  /* set density from DLAW command, must be done here since it may depend on later input commands */
497  if( strcmp(dense.chDenseLaw,"DLW1") == 0 )
498  {
500 
501  if( dense.gas_phase[ipHYDROGEN] <= 0. )
502  {
503  fprintf( ioQQQ, " PROBLEM DISASTER Hydrogen density set by DLAW must be > 0.\n" );
505  }
506  }
507  else if( strcmp(dense.chDenseLaw,"DLW2") == 0 )
508  {
510 
511  if( dense.gas_phase[ipHYDROGEN] <= 0. )
512  {
513  fprintf( ioQQQ, " PROBLEM DISASTER Hydrogen density set by DLAW must be > 0.\n" );
515  }
516  }
517  else if( strcmp(dense.chDenseLaw,"DLW3") == 0 )
518  {
520 
521  if( dense.gas_phase[ipHYDROGEN] <= 0. )
522  {
523  fprintf( ioQQQ, " PROBLEM DISASTER Hydrogen density set by DLAW must be > 0.\n" );
525  }
526  }
527 
528  /* start checks on parameters set properly - this begins with same line saying start .. */
529 
530  /* lgStop_not_enough_info says that not enough info for model, so stop
531  * set true in following tests if anything missing */
532  lgStop_not_enough_info = false;
533  lgStop = false;
534 
535  /* check whether hydrogen density has been set - this value was set to 0 in zero */
536  if( dense.gas_phase[ipHYDROGEN] <= 0. )
537  {
538  fprintf( ioQQQ, " PROBLEM DISASTER Hydrogen density MUST be specified.\n" );
539  lgStop_not_enough_info = true;
540  lgStop = true;
541  /* need to set something since used below - will abort
542  * since lgStop is set */
544  }
545 
546  /* the SAVE XSPEC command cannot be combined with negative increments on the GRID command */
548  {
549  if( called.lgTalk )
550  {
551  fprintf( ioQQQ, " PROBLEM DISASTER The SAVE XSPEC command cannot be combined with negative grid increments.\n" );
552  fprintf( ioQQQ, " PROBLEM DISASTER Please check your GRID commands.\n\n\n" );
554  }
555  }
556 
557  if( geometry.covaper < 0.f || geometry.iEmissPower == 2 )
559 
561 
562  /* mass flux for wind model - used for mass conservation */
564 
565  /* set converge criteria - limit number of iterations and zones */
567  {
569  for( long int j=0; j < iterations.iter_malloc; j++ )
570  {
572  }
573  }
574 
575  /* NSAVE is number of lines saved, =0 if no commands entered */
576  if( input.nSave < 0 )
577  {
578  fprintf( ioQQQ, " PROBLEM DISASTER No commands were entered - whats up?\n" );
580  }
581 
582  /* iterate to convergence and wind models are mutually exclusive */
583  if( wind.lgBallistic() && conv.lgAutoIt )
584  {
585  if( called.lgTalk )
586  {
587  fprintf( ioQQQ, " NOTE PROBLEM Due to the nature of the Sobolev approximation, it makes no sense to converge a windy model.\n" );
588  fprintf( ioQQQ, " NOTE Iterate to convergence is turned off\n\n\n" );
589  }
590  conv.lgAutoIt = false;
591  iterations.itermx = 0;
592  }
593 
594  /* iterate to convergence and case b do not make sense together */
595  /* WJH 22 May 2004: unless we are doing i-front dynamics (negative wind.windv0) */
596  if( opac.lgCaseB && conv.lgAutoIt && (wind.lgBallistic() || wind.lgStatic()) )
597  {
598  if( called.lgTalk )
599  {
600  fprintf( ioQQQ, " NOTE Case B is an artificial test, it makes no sense to converge this model.\n" );
601  fprintf( ioQQQ, " NOTE Iterate to convergence is turned off.\n\n\n" );
602  }
603  conv.lgAutoIt = false;
604  iterations.itermx = 0;
605  }
606 
607  /* specifying a density power and constant pressure makes no sense */
608  if( dense.DensityPower!=0. && strcmp( dense.chDenseLaw, "CPRE" )==0 )
609  {
610  if( called.lgTalk )
611  {
612  fprintf( ioQQQ, " NOTE Specifying both a density power law and constant pressure is impossible.\n" );
613  }
614  lgStop = true;
615  }
616 
617  if( !rfield.lgIonizReevaluate && strcmp( dense.chDenseLaw, "CDEN" )!=0 )
618  {
619  if( called.lgTalk )
620  {
621  fprintf( ioQQQ, " NOTE NO REEVALUATE IONIZATION can only be used with constant density.\n" );
622  fprintf( ioQQQ, " NOTE Resetting to reevaluate ionization.\n\n" );
623  }
624  rfield.lgIonizReevaluate = true;
625  }
626 
627  if( !rfield.lgOpacityReevaluate && strcmp( dense.chDenseLaw, "CDEN" )!=0 )
628  {
629  if( called.lgTalk )
630  {
631  fprintf( ioQQQ, " NOTE NO REEVALUATE OPACITY can only be used with constant density.\n" );
632  fprintf( ioQQQ, " NOTE Resetting to reevaluate opacity.\n\n" );
633  }
635  }
636 
637  /* check that a symmetry is specified if gravity from an external mass has been added */
638  if( pressure.external_mass[0].size()>0 && pressure.gravity_symmetry==-1 )
639  {
640  if( called.lgTalk )
641  {
642  fprintf( ioQQQ, " NOTE Gravity from an external mass has been added, but no symmetry (spherical/mid-plane) was specified.\n" );
643  fprintf( ioQQQ, " NOTE It will be ignored.\n\n\n" );
644  }
645  }
646 
647  /* check if the combination of stopping column density and H density are physically plausible */
648  double thickness = min( MIN3( StopCalc.tauend, StopCalc.colpls, StopCalc.colnut ),
650  if( thickness < COLUMN_INIT )
651  {
652  /* a stop column density was specified - check on physical thickness this corresponds to */
653  thickness /= (dense.gas_phase[ipHYDROGEN]*geometry.FillFac);
654  /* don't complain if outer radius set small with `stop thickness' command. */
655  if( thickness > 1e25 && iterations.StopThickness[0] > 1e25 )
656  {
657  fprintf( ioQQQ,
658  "NOTE The specified column density and hydrogen density correspond to a thickness of %.2e cm.\n",
659  thickness);
660  fprintf( ioQQQ,
661  "NOTE This seems large to me.\n");
662  fprintf(ioQQQ,"NOTE a very large radius may cause overflow.\n\n");
663  }
664  }
665 
667  {
668  /* warn if constant grain temperature but gas-grain thermal effects
669  * are still included */
670  fprintf( ioQQQ,
671  "NOTE The grain temperatures are set to a constant value with the "
672  "CONSTANT GRAIN TEMPERATURE command, but "
673  "energy exchange \n");
674  fprintf( ioQQQ,
675  "NOTE is still included. The grain-gas heating-cooling will be incorrect. "
676  "Consider turning off gas-grain collisional energy\n");
677  fprintf( ioQQQ,
678  "NOTE exchange with the NO GRAIN GAS COLLISIONAL ENERGY EXCHANGE command.\n\n\n");
679  }
680 
682  {
683  if( called.lgTalk )
684  {
685  fprintf( ioQQQ, " NOTE NO LINE TRANSER set but fine opacities still computed.\n" );
686  fprintf( ioQQQ, " NOTE Turning off fine opacities.\n\n" );
687  }
688  rfield.lgOpacityFine = false;
689  }
690 
692  {
693  if( called.lgTalk )
694  {
695  fprintf( ioQQQ, " NOTE Large H2 molecule turned on but line transfer and fine opacities are not.\n" );
696  fprintf( ioQQQ, " NOTE Turning on line transfer and fine opacities.\n\n" );
697  }
698  rfield.lgOpacityFine = true;
699  rfield.lgDoLineTrans = true;
700  }
701 
703  {
704  /* one of the input continua needs to have H-ionizing radiation
705  * blocked with extinguish command, but it was not done */
706  fprintf( ioQQQ, "\n NOTE\n"
707  " NOTE One of the incident continuum is a form used when no H-ionizing radiation is produced.\n" );
708  fprintf( ioQQQ, " NOTE You must also include the EXTINGUISH command to make sure this is done.\n" );
709  fprintf( ioQQQ, " NOTE The EXTINGUISH command was not included.\n" );
710  fprintf( ioQQQ, " NOTE YOU MAY BE MAKING A BIG MISTAKE!!\n NOTE\n\n\n\n" );
711  }
712 
713  /* if stop temp set below default then we are going into cold and possibly molecular
714  * gas - check some parameters in this case */
716  /* thermal.ConstTemp def is zero, set pos when constant temperature is set */
718  {
719  /* print warning if temperature set below default but cosmic rays not turned on
720  * do not print if molecules are off */
721  if( (hextra.cryden == 0.) && !mole_global.lgNoMole )
722  {
723  fprintf( ioQQQ, "\n NOTE\n"
724  " NOTE The simulation is going into neutral gas but cosmic rays are not included.\n" );
725  fprintf( ioQQQ, " NOTE Ion-molecule chemistry will not occur without a source of ionization.\n" );
726  fprintf( ioQQQ, " NOTE The chemistry network may collapse deep in molecular regions.\n" );
727  fprintf( ioQQQ, " NOTE Consider adding galactic background cosmic rays with the COSMIC RAYS BACKGROUND command.\n" );
728  fprintf( ioQQQ, " NOTE You may be making a BIG mistake.\n NOTE\n\n\n\n" );
729  }
730  }
731 
732  /* dense.gas_phase[ipHYDROGEN] is linear hydrogen density (cm-3) */
733  /* test for hydrogen density properly set has already been done above */
734  if( called.lgTalk && dense.gas_phase[ipHYDROGEN] < 1e-4 )
735  {
736  fprintf( ioQQQ, " NOTE Is the entered value of the hydrogen density (%.2e) reasonable?\n",
738  fprintf( ioQQQ, " NOTE It seems pretty low to me.\n\n\n" );
739  }
740  else if( called.lgTalk && dense.gas_phase[ipHYDROGEN] > 1e15 )
741  {
742  fprintf( ioQQQ, " NOTE Is this value of the hydrogen density reasonable?\n" );
743  fprintf( ioQQQ, " NOTE It seems pretty high to me.\n\n\n" );
744  }
745 
746  /* is the model going to crash because of extreme density? */
747  if( called.lgTalk && !lgStop && !lgStop_not_enough_info )
748  {
749  if( dense.gas_phase[ipHYDROGEN] < 1e-6 || dense.gas_phase[ipHYDROGEN] > 1e19 )
750  {
751  fprintf( ioQQQ, " NOTE Simulation may crash because of extreme "
752  "density. The value was %.2e\n\n" ,
754  }
755  }
756 
757  if( rfield.nShape == 0 && (p.m_nqh) == 0 )
758  {
759  // no incident radiation field, at all
760  // This is ok if temperature or heating is specified
761  if( thermal.ConstTemp <=0 && hextra.TurbHeat<=0. )
762  {
763  fprintf( ioQQQ, " PROBLEM DISASTER No incident radiation field was specified - "
764  "at least put in the CMB.\n" );
765  lgStop = true;
766  lgStop_not_enough_info = true;
767  }
768 
769  }
770  else if( rfield.nShape == 0 )
771  {
772  fprintf( ioQQQ, " PROBLEM DISASTER No incident radiation field was specified - "
773  "at least put in the CMB.\n" );
774  lgStop = true;
775  lgStop_not_enough_info = true;
776  }
777  else if( (p.m_nqh) == 0 )
778  {
779  fprintf( ioQQQ, " PROBLEM DISASTER Luminosity of continuum MUST be specified.\n" );
780  lgStop = true;
781  lgStop_not_enough_info = true;
782  }
783 
784  /* Testing on 0 is safe since this it is initialized that way
785  * only print comment if continuum has been specified,
786  * if continuum not given then we are aborting anyway */
787  if( radius.Radius == 0. && rfield.nShape> 0)
788  {
789  fprintf( ioQQQ, " PROBLEM DISASTER Starting radius MUST be specified.\n" );
790  lgStop = true;
791  lgStop_not_enough_info = true;
792  }
793 
794  if( rfield.nShape != (p.m_nqh) )
795  {
796  fprintf( ioQQQ, " PROBLEM DISASTER There were not the same number of continuum shapes and luminosities entered.\n" );
797  lgStop = true;
798  }
799 
800  /* we only want to do this test on the first call to the command
801  * parser - it will be called many more times but with no grid command
802  * during the actual grid calculation */
803  static bool lgFirstPass = true;
804 
805  /* the NO VARY option sets this flag, and can be used to turn off
806  * the grid command, as well as the optimizer */
807  if( optimize.lgNoVary && grid.lgGrid )
808  {
809  /* ignore grids */
810  grid.lgGrid = false;
811  optimize.nparm = 0;
812  grid.nGridCommands = 0;
813  }
814 
815  if( lgFirstPass && grid.lgGrid && (optimize.nparm!=grid.nGridCommands) )
816  {
817  /* number of grid vary options do match */
818  fprintf( ioQQQ, " PROBLEM DISASTER The GRID command was entered "
819  "but there were %li GRID commands and %li commands with a VARY option.\n" ,
821  fprintf( ioQQQ, " There must be the same number of GRIDs and VARY.\n" );
822  lgStop = true;
823  }
824  lgFirstPass = false;
825 
826  if( lgStop_not_enough_info )
827  {
828  fprintf( ioQQQ, " PROBLEM DISASTER I do not have enough information to do the simulation, I cannot go on.\n" );
829  fprintf( ioQQQ, "\n\n Sorry.\n\n\n" );
831  }
832 
833  {
834  bool lgParserTest = false;
835  if ( lgParserTest )
836  {
837  // Quit after parse phase, to allow quick verification that
838  // there are no immediate syntax handling failures.
839  fprintf(ioQQQ,"Parser phase PASSED\n");
840  exit(0);
841  }
842  }
843 
844  if( lgStop )
845  {
847  }
848 
849  /* end checks on parameters set properly - this begins with same line saying start .. */
850  return;
851 }
852 
854 {
855  DEBUG_ENTRY( "ParseAperture()" );
856  /* aperture command to simulate pencil beam or long slit */
857 
858  /* if the "BEAM" or "SLIT" option is specified then only part
859  * of the geometry is observed, and intensities
860  * should not be weighted by r^2. There are two limiting cases, SLIT in which
861  * the slit is longer than the diameter of the nebula and the contribution to the
862  * detected luminosity goes as r^1, and BEAM when the contribution is r^0,
863  * the same as plane parallel
864  *
865  * default value of geometry.iEmissPower is 2 (set in zero.c) for full geometry
866  */
867  if( p.nMatch("SLIT") )
868  {
869  /* long slit is case where slit is longer than diameter, so emissivity contributes
870  * r^1 to the observed luminosity */
871  geometry.iEmissPower = 1;
872  }
873  else if( p.nMatch("BEAM") )
874  {
875  /* pencil beam is case where we view the nebula through a narrow square
876  * centered on the central source, this gives r^0 dependence */
877  geometry.iEmissPower = 0;
878  }
879  else if( p.nMatch("SIZE") )
880  {
881  /* set the aperture size: slit width or suface area of the pencil beam */
882  /* units are arcsec for slit width, or arcsec^2 for pencil beam */
883  geometry.size = realnum(p.FFmtRead());
884  if( p.lgEOL() )
885  {
886  p.NoNumb("aperture size");
887  }
888  if( geometry.size <= 0.f )
889  {
890  fprintf( ioQQQ, " The aperture size must be positive. Sorry.\n" );
892  }
893  geometry.lgSizeSet = true;
894  }
895  else if( p.nMatch("COVE") )
896  {
897  /* set the aperture covering factor, see Hazy 1 for a discussion
898  * this is a dimensionless fraction between 0 and 1 */
900  if( p.lgEOL() )
901  {
902  p.NoNumb("aperture covering factor");
903  }
904  if( geometry.covaper <= 0.f || geometry.covaper > 1.f )
905  {
906  fprintf( ioQQQ, " The aperture covering factor must be > 0 and <= 1. Sorry.\n" );
908  }
909  }
910  else
911  {
912  fprintf( ioQQQ, " One of the keywords SLIT, BEAM, SIZE or COVEring factor must appear.\n" );
913  fprintf( ioQQQ, " Sorry.\n" );
915  }
916 }
917 
919 {
920  DEBUG_ENTRY( "ParseDatabase()" );
921 
922  string chString_quotes_original;
923  bool lgQuotesFound = true;
924  if (p.GetQuote(chString_quotes_original))
925  lgQuotesFound = false;
926 
927  /* accept both forms of feii */
928  if( p.nMatch("FEII") || p.nMatch("FE II") )
929  {
930  fprintf( ioQQQ, " Warning: The 'atom feii' command is obsolete. "
931  " Instead, please use 'species \"Fe+\" levels=all'.\n Sorry.\n\n" );
932  cdEXIT( EXIT_FAILURE );
933  }
934 
935  else if( p.nMatch("H-LI") )
936  {
937  /* parse the atom h-like command */
939  }
940 
941  else if( p.nMatch("HE-L") )
942  {
943  /* parse the atom he-like command */
945  }
946 
947  else if( p.nMatch("ROTO") || p.nMatch(" CO ") )
948  {
949  fprintf(ioQQQ," The old CO models no longer exist, and this command is no longer supported.\n" );
950  fprintf( ioQQQ, " Sorry.\n" );
952  }
953 
954  else if( p.nMatch(" H2 ") )
955  {
956  ParseDatabaseH2(p);
957  }
958 
959  /* enable models from CHIANTI database */
960  else if (p.nMatch("CHIANTI"))
961  {
962  // option to specify different CloudyChianti.ini file, was initialized
963  // with string CloudyChianti.ini
964  if (lgQuotesFound == true)
965  strcpy(atmdat.chCloudyChiantiFile, chString_quotes_original.c_str());
966 
967  atmdat.lgChiantiOn = true;
968  if (p.nMatch(" OFF"))
969  {
970  atmdat.lgChiantiOn = false;
971  atmdat.lgChiantiHybrid = false;
972  }
973 
974  // hybrid, chianti with OP for higher energy lines
975  // Turn off hybrid, use Chianti only
976  if (p.nMatch("NO HYBR"))
977  atmdat.lgChiantiHybrid = false;
978 
979  // Print which species are being used in output and # of levels
980  if (p.nMatch("PRINT"))
981  atmdat.lgChiantiPrint = true;
982 
983  // Use Experimental energies exclusively. Default use experimental.
984  if (p.nMatch("THEOR"))
985  atmdat.lgChiantiExp = false;
986 
987  // Input the maximum number of Chianti levels to use
988  if (p.nMatch("LEVEL"))
989  {
990  if (p.nMatch(" MAX"))
991  {
992  atmdat.nChiantiMaxLevels = LONG_MAX;
993  atmdat.nChiantiMaxLevelsFe = LONG_MAX;
994  }
995  else
996  {
997  atmdat.nChiantiMaxLevelsFe = (long)p.FFmtRead();
998  atmdat.nChiantiMaxLevels = (long)p.FFmtRead();
999  if (p.lgEOL())
1000  {
1001  p.NoNumb("two numbers, the maximum number of levels in Fe, and in other elements, or the keyword MAX,");
1002  }
1004  {
1005  fprintf(ioQQQ,
1006  " \nPROBLEM The maximum number of chianti levels should be two or greater.\n");
1007  fprintf(ioQQQ, " To turn off the Chianti data use \"atom Chianti off\" instead.\n");
1008  fprintf(ioQQQ, " See Hazy 1 for details.\n");
1009  cdEXIT( EXIT_FAILURE );
1010  }
1011  }
1012  atmdat.lgChiantiLevelsSet = true;
1013  }
1014  }
1015 
1016  /* enable models from STOUT database */
1017  else if (p.nMatch("STOUT"))
1018  {
1019 
1020  // option to specify different Stout.ini file, was initialized
1021  // with string Stout.ini
1022  if (lgQuotesFound == true)
1023  strcpy(atmdat.chStoutFile, chString_quotes_original.c_str());
1024 
1025  atmdat.lgStoutOn = true;
1026  // Print which species are being used in output and # of levels
1027  if (p.nMatch("PRINT"))
1028  atmdat.lgStoutPrint = true;
1029 
1030  if (p.nMatch(" OFF"))
1031  {
1032  atmdat.lgStoutOn = false;
1033  atmdat.lgStoutHybrid = false;
1034  }
1035 
1036  // hybrid, Stout with OP for higher energy lines
1037  // Turn off hybrid, use Stout only
1038  if (p.nMatch("NO HYBR"))
1039  atmdat.lgStoutHybrid = false;
1040 
1041  // Input the maximum number of Stout levels to use
1042  if (p.nMatch("LEVEL"))
1043  {
1044  if (p.nMatch(" MAX"))
1045  {
1046  atmdat.nStoutMaxLevels = LONG_MAX;
1047  atmdat.nStoutMaxLevelsFe = LONG_MAX;
1048  }
1049  else
1050  {
1051  atmdat.nStoutMaxLevelsFe = (long)p.FFmtRead();
1052  atmdat.nStoutMaxLevels = (long)p.FFmtRead();
1053  if (p.lgEOL())
1054  {
1055  p.NoNumb("two numbers, the maximum number of levels in Fe, and in other elements, or the keyword MAX,");
1056  }
1058  {
1059  fprintf(ioQQQ,
1060  " \nPROBLEM The maximum number of stout levels should be two or greater.\n");
1061  fprintf(ioQQQ, " To turn off the Stout data use \"atom Stout off\" instead.\n");
1062  fprintf(ioQQQ, " See Hazy 1 for details.\n");
1063  cdEXIT( EXIT_FAILURE );
1064  }
1065  }
1066  atmdat.lgStoutLevelsSet = true;
1067  }
1068  }
1069 
1070  /* enable models from LAMDA (Leiden Atomic and Molecular Database) */
1071  else if (p.nMatch("LAMDA"))
1072  {
1073  // option to specify different Lamda.ini file, was initialized
1074  // with string Lamda.ini
1075  if (lgQuotesFound == true)
1076  strcpy(atmdat.chLamdaFile, chString_quotes_original.c_str());
1077 
1078  atmdat.lgLamdaOn = true;
1079  // Print which species are being used in output and # of levels
1080  if (p.nMatch("PRINT"))
1081  atmdat.lgLamdaPrint = true;
1082 
1083  if(p.nMatch(" OFF"))
1084  {
1085  atmdat.lgLamdaOn = false;
1086  }
1087 
1088  // Input the maximum number of Lamda levels to use
1089  if (p.nMatch("LEVEL"))
1090  {
1091  if(p.nMatch(" MAX"))
1092  {
1093  atmdat.nLamdaMaxLevels = LONG_MAX;
1094  }
1095  else
1096  {
1097  atmdat.nLamdaMaxLevels = (long)p.FFmtRead();
1098  if( p.lgEOL())
1099  {
1100  p.NoNumb("the maximum number of levels,");
1101  }
1102  if( atmdat.nLamdaMaxLevels < 2 )
1103  {
1104  fprintf(ioQQQ,
1105  " \nPROBLEM The maximum number of Lamda levels should be two or greater.\n");
1106  fprintf(ioQQQ, " To turn off the Lamda data use \"atom lamda off\" instead.\n");
1107  fprintf(ioQQQ, " See Hazy 1 for details.\n");
1108  cdEXIT( EXIT_FAILURE );
1109  }
1110  }
1111  atmdat.lgLamdaLevelsSet = true;
1112  }
1113  }
1114 
1115  /* no single species, but print summaries of all */
1116  else if (p.nMatch("PRINT"))
1117  {
1118  atmdat.lgChiantiPrint = true;
1119  atmdat.lgLamdaPrint = true;
1120  atmdat.lgStoutPrint = true;
1122  }
1123 
1124  else
1125  {
1126  fprintf( ioQQQ, " I could not recognize a keyword on this species command.\n");
1127  fprintf( ioQQQ, " The available keys are FeII, H-Like, He-like, H2, Chianti, Lamda, and Stout.\n");
1128  fprintf( ioQQQ, " Sorry.\n" );
1130  }
1131 }
1133 {
1134  DEBUG_ENTRY( "ParseBremsstrahlung()" );
1135  /* bremsstrahlung continuum from central object */
1136  strcpy( rfield.chSpType[rfield.nShape], "BREMS" );
1138  (realnum)p.FFmtRead();
1139  if( p.lgEOL() )
1140  {
1141  p.NoNumb("temperature");
1142  }
1143 
1144  /* temperature is interpreted as log if <= 10, or if keyword is present */
1145  if( rfield.slope[rfield.nShape] <= 10. || p.nMatch(" LOG") )
1146  {
1149  }
1150  rfield.cutoff[rfield.nShape][0] = 0.;
1151 
1152  /* option for vary keyword */
1153  if( optimize.lgVarOn )
1154  {
1155  /* only one parameter */
1157  strcpy( optimize.chVarFmt[optimize.nparm], "BREMS, T=%f LOG" );
1158  /* pointer to where to write */
1160  /* log of temp will be pointer */
1162  optimize.vincr[optimize.nparm] = 0.5;
1163  ++optimize.nparm;
1164  }
1165  ++rfield.nShape;
1166  if( rfield.nShape >= LIMSPC )
1167  {
1168  /* too many continua were entered */
1169  fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
1171  }
1172 }
1174 {
1175  /* add "extra" cooling, power law temp dependence */
1176  thermal.lgCExtraOn = true;
1178  if( p.lgEOL() )
1179  {
1180  p.NoNumb("extra cooling");
1181  }
1182  thermal.cextpw = (realnum)p.FFmtRead();
1183 }
1185 {
1188 
1189  /* >>chng 06 mar 22, add time option to vary only some continua with time */
1190  if( p.nMatch( "TIME" ) )
1191  rfield.lgTimeVary[(p.m_nqh)] = true;
1192 
1194 
1195  /* option to also set the hydrogen density at given redshift. */
1196  if( p.nMatch("DENS") && !p.lgEOL() )
1197  {
1198  char chStuff[INPUT_LINE_LENGTH];
1199 
1200  /* hydrogen density */
1201  sprintf( chStuff , "HDEN %.2e LINEAR", GetDensity( cosmology.redshift_start ) );
1202  p.setline(chStuff);
1203  p.set_point(4);
1204  ParseHDEN(p);
1205  }
1206 }
1208 {
1209  DEBUG_ENTRY( "ParseCosm()" );
1210  fprintf(ioQQQ,
1211  "This command is now ambiguous -- please specify either COSMIC RAYS or COSMOLOGY.\nSorry.\n");
1213 }
1215 {
1216  DEBUG_ENTRY( "ParseCovering()" );
1217  /* covering factor for gas */
1218  /* The geometric covering factor accounts for how much of 4\pi is
1219  * covered by gas, and so linearly multiplies the predicted intensities */
1221 
1222  if( p.lgEOL() )
1223  {
1224  p.NoNumb("covering factor");
1225  }
1226 
1227  /* if negative then log, convert to linear */
1228  if( geometry.covgeo <= 0. )
1229  {
1231  }
1232 
1233  if( geometry.covgeo > 1. )
1234  {
1235  fprintf( ioQQQ, " A covering factor greater than 1 makes no physical sense. Sorry.\n" );
1237  }
1238 
1239  /* radiative transfer covering factor will be equal to the geometric one */
1241 }
1243 {
1244  /* height of cylinder in cm */
1245  radius.lgCylnOn = true;
1247  if( p.lgEOL() )
1248  {
1249  p.NoNumb("height");
1250  }
1251 }
1253 {
1254  DEBUG_ENTRY( "ParseDarkMatter()" );
1255 
1256  if( p.nMatch(" NFW") )
1257  {
1258  /* specify an NFW profile */
1259  /* >>refer dark matter Navarro, Frenk, & White, 1996, ApJ, 462, 563 */
1260 
1261  dark.r_200 = (realnum)p.getNumberCheckAlwaysLog("NFW r_200");
1262  /* Let r_s have a default corresponding to c_200 = 10. */
1263  dark.r_s = (realnum)p.getNumberDefaultAlwaysLog("NFW r_s", log10(dark.r_200)-1. );
1264  dark.lgNFW_Set = true;
1265 
1266  /* option for vary keyword */
1267  if( optimize.lgVarOn )
1268  {
1269  /* only one parameter */
1271  strcpy( optimize.chVarFmt[optimize.nparm], "DARK NFW %f" );
1272  /* pointer to where to write */
1274  /* log of temp will be pointer */
1275  optimize.vparm[0][optimize.nparm] = (realnum)log10(dark.r_200);
1276  optimize.vincr[optimize.nparm] = 0.5;
1277  ++optimize.nparm;
1278  }
1279  }
1280  else
1281  {
1282  fprintf( ioQQQ, " Did not recognize a valid option for this DARK command.\nSorry.\n\n" );
1284  }
1285 
1286  return;
1287 }
1289 {
1290  DEBUG_ENTRY( "ParseDielectronic()" );
1291  /* >>chng 05 dec 21, change dielectronic command to set dielectronic recombination command */
1292  fprintf( ioQQQ, " The DIELectronic command has been replaced with the SET DIELectronic recombination command.\n" );
1293  fprintf( ioQQQ, " Please have a look at Hazy.\n Sorry.\n\n" );
1295 }
1297 {
1298  DEBUG_ENTRY( "ParseDiffuse()" );
1299  /* set method of transferring diffuse fields,
1300  * default is outward only, cdDffTrns label "OU2", set in zero.c */
1301  if( p.nMatch(" OTS") )
1302  {
1303  if( p.nMatch("SIMP") )
1304  {
1305  /* this is simple ots, which never allows photons to escape */
1306  strcpy( rfield.chDffTrns, "OSS" );
1307  }
1308  else
1309  {
1310  /* this is regular ots, which allows photons to escape */
1311  strcpy( rfield.chDffTrns, "OTS" );
1312  }
1313  rfield.lgOutOnly = false;
1314  }
1315  else if( p.nMatch(" OUT") )
1316  {
1317  rfield.lgOutOnly = true;
1318  long int j = (long int)p.FFmtRead();
1319  if( p.lgEOL() )
1320  {
1321  /* this is the default set in zero */
1322  strcpy( rfield.chDffTrns, "OU2" );
1323  }
1324  else
1325  {
1326  if( j > 0 && j < 10 )
1327  {
1328  sprintf( rfield.chDffTrns, "OU%1ld", j );
1329  }
1330  else
1331  {
1332  fprintf( ioQQQ, " must be between 1 and 9 \n" );
1334  }
1335  }
1336  }
1337 
1338  else
1339  {
1340  fprintf( ioQQQ, " There should have been OUTward or OTS on this line. Sorry.\n" );
1342  }
1343 }
1345 {
1346  /* distance to the object */
1347  radius.distance = p.FFmtRead();
1348  if( p.lgEOL() )
1349  {
1350  p.NoNumb("distance");
1351  }
1352  /* default is for quantity to be log of distance, linear keyword
1353  * overrides this - is LINEAR is not on line then exp */
1354  if( !p.nMatch("LINE" ) )
1355  {
1357  }
1358  /* default is radius in cm - if parsecs appears then must
1359  * convert to cm */
1360  if( p.nMatch("PARS") )
1361  {
1362  radius.distance *= PARSEC;
1363  }
1364 
1365  /* vary option */
1366  if( optimize.lgVarOn )
1367  {
1368  strcpy( optimize.chVarFmt[optimize.nparm], "DISTANCE %f LOG" );
1371  optimize.vincr[optimize.nparm] = 0.3f;
1373  ++optimize.nparm;
1374  }
1375 }
1377 {
1378  rt.DoubleTau = 2.;
1379 }
1381 {
1383  if( p.lgEOL() )
1384  {
1385  p.NoNumb("electron density");
1386  }
1387  /* warn that this model is meaningless */
1388  phycon.lgPhysOK = false;
1389 }
1391 {
1392  DEBUG_ENTRY( "ParseEnergy()" );
1393  /* energy density (degrees K) of this continuum source */
1394  if( (p.m_nqh) >= LIMSPC )
1395  {
1396  /* too many continua were entered */
1397  fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
1399  }
1400  /* silly, but calms down the lint */
1401  ASSERT( (p.m_nqh) < LIMSPC );
1402  strcpy( rfield.chRSpec[(p.m_nqh)], "SQCM" );
1403  realnum teset = (realnum)p.FFmtRead();
1404  if( p.lgEOL() )
1405  {
1406  p.NoNumb("energy density");
1407  }
1408  /* set radius to very large value if not already set */
1409  /* >>chng 01 jul 24, from Radius == 0 to this, as per PvH comments */
1410  if( !radius.lgRadiusKnown )
1411  {
1412  /* set radius to large value */
1414  }
1415 
1416  /* convert temp to log, recognize linear option */
1417  if( !p.nMatch(" LOG") && (p.nMatch("LINE") || teset > 10.) )
1418  {
1419  /* option for linear temperature, must store log */
1420  teset = (realnum)log10(teset);
1421  }
1422 
1423  if( teset > 5. )
1424  {
1425  fprintf( ioQQQ, " This intensity may be too large. The code may crash due to overflow. Was log intended?\n" );
1426  }
1427 
1428  /* teset is now log of temp, now get log of total luminosity */
1429  strcpy( rfield.chSpNorm[(p.m_nqh)], "LUMI" );
1430 
1431  /* full range of continuum will be used */
1432  rfield.range[(p.m_nqh)][0] = rfield.emm();
1433  rfield.range[(p.m_nqh)][1] = rfield.egamry();
1434  rfield.totpow[(p.m_nqh)] = (4.*(teset) - 4.2464476 + 0.60206);
1435 
1436  /* >>chng 06 mar 22, add time option to vary only some continua with time */
1437  if( p.nMatch( "TIME" ) )
1438  rfield.lgTimeVary[(p.m_nqh)] = true;
1439 
1440  /* vary option */
1441  if( optimize.lgVarOn )
1442  {
1443  strcpy( optimize.chVarFmt[optimize.nparm], "ENERGY DENSITY %f LOG" );
1444  if( rfield.lgTimeVary[(p.m_nqh)] )
1445  strcat( optimize.chVarFmt[optimize.nparm], " TIME" );
1446  /* pointer to where to write */
1448  optimize.vparm[0][optimize.nparm] = teset;
1449  optimize.vincr[optimize.nparm] = 0.1f;
1451  ++optimize.nparm;
1452  }
1453 
1454  /* finally increment number of continua */
1455  ++(p.m_nqh);
1456 }
1458 {
1459  /* save previous value */
1460  long int j = conv.LimFail;
1461  conv.LimFail = (long int)p.FFmtRead();
1462  if( p.lgEOL() )
1463  {
1464  p.NoNumb("limit");
1465  }
1466 
1467  /* >>chng 01 mar 14, switch logic on maps */
1468  /* ' map' flag, make map when failure, default is no map,
1469  * second check is so that no map does not trigger a map */
1470  if( p.nMatch(" MAP") && !p.nMatch(" NO ") )
1471  {
1472  conv.lgMap = true;
1473  }
1474 
1475  /* complain if failures was increased above default */
1476  if( conv.LimFail > j )
1477  {
1478  fprintf( ioQQQ, " This command should not be necessary.\n" );
1479  fprintf( ioQQQ, " Please show this input stream to Gary Ferland if this command is really needed for this simulation.\n" );
1480  }
1481 }
1483 {
1484  /* filling factor, power law exponent (default=1., 0.) */
1485  realnum a = (realnum)p.FFmtRead();
1486  if( p.lgEOL() )
1487  {
1488  p.NoNumb("filling factor");
1489  }
1490 
1491  if( a <= 0. || p.nMatch(" LOG") )
1492  {
1493  /* number less than or equal to 0, must have been entered as log */
1494  geometry.FillFac = exp10(a);
1495  }
1496  else
1497  {
1498  /* number greater than zero, must have been the real thing */
1499  geometry.FillFac = a;
1500  }
1501  if( geometry.FillFac > 1.0 )
1502  {
1503  if( called.lgTalk )
1504  fprintf( ioQQQ, " Filling factor > 1, reset to 1\n" );
1505  geometry.FillFac = 1.;
1506  }
1508 
1509  /* option to have power law dependence, filpow set to 0 in zerologic */
1511 
1512  /* vary option */
1513  if( optimize.lgVarOn )
1514  {
1515  strcpy( optimize.chVarFmt[optimize.nparm], "FILLING FACTOR= %f LOG power= %f" );
1516 
1517  /* pointer to where to write */
1520  optimize.vincr[optimize.nparm] = 0.5f;
1521 
1522  /* power law dependence here, but cannot be varied */
1525 
1526  /* do not allow filling factor to go positive */
1527  optimize.varang[optimize.nparm][0] = -FLT_MAX;
1528  optimize.varang[optimize.nparm][1] = 0.f;
1529  ++optimize.nparm;
1530  }
1531 }
1533 {
1534  bool lgNu2 = false;
1535  /* in reads2 */
1536  ParseF_nu(p,"SQCM",lgNu2);
1537 }
1539 {
1541  if( p.lgEOL() )
1542  {
1543  p.NoNumb("temperature");
1544  }
1545 
1546  if( p.nMatch(" LOG") || (thermal.ConstTemp <= 10. &&
1547  !p.nMatch("LINE")) )
1548  {
1550  }
1551 
1552  /* low energy bounds of continuum array do not permit too-low a temp */
1553  if( thermal.ConstTemp < 3. )
1554  {
1555  fprintf( ioQQQ, " TE reset to 3K: entered number too small.\n" );
1556  thermal.ConstTemp = 3.;
1557  }
1558 }
1560 {
1561  /* enter a fudge factor */
1562  fudgec.nfudge = 0;
1563  for( long int j=0; j < NFUDGC; j++ )
1564  {
1565  fudgec.fudgea[j] = (realnum)p.FFmtRead();
1566  /* fudgec.nfudge is number of factors on the line */
1567  if( !p.lgEOL() )
1568  fudgec.nfudge = j+1;
1569  }
1570  if( fudgec.nfudge == 0 )
1571  p.NoNumb("fudge factor");
1572 
1573  /* vary option */
1574  if( optimize.lgVarOn )
1575  {
1576  /* no luminosity options on vary */
1578  strcpy( optimize.chVarFmt[optimize.nparm], "FUDGE= %f" );
1579  /* enter remaining parameters as constants */
1580  char chHold[1000];
1581  for( long int j=1; j<fudgec.nfudge; ++j )
1582  {
1583  sprintf( chHold , " %f" , fudgec.fudgea[j] );
1584  strcat( optimize.chVarFmt[optimize.nparm] , chHold );
1585  }
1587 
1588  /* pointer to where to write */
1590  /* fudge factor stored here */
1592  /* the increment in the first steps away from the original value */
1593  optimize.vincr[optimize.nparm] = abs(0.2f*fudgec.fudgea[0]);
1594  /* we have no clue what to use when initial estimate is 0... */
1595  if( optimize.vincr[optimize.nparm] == 0.f )
1596  optimize.vincr[optimize.nparm] = 1.f;
1597  ++optimize.nparm;
1598  }
1599 }
1601 {
1602  DEBUG_ENTRY( "ParsePGrains()" );
1603  fprintf(ioQQQ," Sorry, this command is obsolete, you can now use the normal GRAINS command.\n");
1605 }
1607 {
1608  DEBUG_ENTRY( "ParseGravity()" );
1609 
1610  if( p.nMatch("EXTE") )
1611  {
1612  /* external mass: M/Msun or Sigma/(Msun/pc^2) depending on symmetry */
1613  /* if no number is read, FFmtRead returns 0 */
1614  /* default is linear, unless "log" keyword is present */
1615  double M_i = p.FFmtRead();
1616 
1617  if( !p.lgEOL() && p.nMatch("LOG") )
1618  M_i = exp10( M_i );
1619  pressure.external_mass[0].push_back( M_i );
1620 
1621  /* extent of the mass distribution, in pc */
1622  /* default is linear, unless "log" keyword is present */
1623  double x_i = p.FFmtRead();
1624 
1625  if( !p.lgEOL() && p.nMatch("LOG") )
1626  x_i = exp10( x_i );
1627  pressure.external_mass[1].push_back( x_i * PARSEC );
1628 
1629  /* exponential index of the mass distribution */
1630  pressure.external_mass[2].push_back( p.FFmtRead() );
1631  }
1632  else
1633  {
1634  /* choose spherical or mid-plane symmetry for the gas mass distribution */
1635  if( p.nMatch("SPHE") )
1636  {
1638  }
1639  else if( p.nMatch("PLAN") )
1640  {
1642  }
1643  else
1644  {
1645  fprintf( ioQQQ, " The symmetry of the gravitational mass must be specified explicitly. Sorry.\n" );
1647  }
1648 
1649  /* external mass, proportional to the gas density */
1650  /* default is linear, unless "log" keyword is present */
1652 
1653  if( p.lgEOL() )
1655  else if( p.nMatch("LOG") )
1656  {
1658  }
1659  }
1660 }
1662 {
1663  DEBUG_ENTRY( "ParseHeLike()" );
1664  fprintf(ioQQQ,"Sorry, this command is replaced with SPECIES HE-LIKE\n");
1666 }
1668 {
1669  DEBUG_ENTRY( "ParseHelp()" );
1670  p.help(ioQQQ);
1671 }
1673 {
1674  DEBUG_ENTRY( "ParseHExtra()" );
1676  if( p.lgEOL() )
1677  p.NoNumb("extra heating first parameter" );
1678 
1679  /* save initial value in case reset with time dependent command */
1681 
1682  /* remember which type of scale dependence we will use */
1683  const char *chHextraScale;
1684  /* these are optional numbers on depth or density option */
1685  realnum par1=0. , par2=0.;
1686  long int nPar;
1687  if( p.nMatch("DEPT") )
1688  {
1689  /* depend on depth */
1690  hextra.lgHextraDepth = true;
1691  chHextraScale = "DEPTH";
1692  /* optional scale radius */
1693  realnum a = (realnum)p.FFmtRead();
1694  if( p.lgEOL() )
1695  p.NoNumb("depth");
1696  else
1697  hextra.turrad = exp10(a);
1698 
1699  /* depth from shielded face, to mimic illumination from both sides
1700  * may not be specified */
1701  realnum b = (realnum)p.FFmtRead();
1702  if( p.lgEOL() )
1703  {
1704  hextra.turback = 0.;
1705  nPar = 2;
1706  }
1707  else
1708  {
1709  hextra.turback = exp10(b);
1710  nPar = 3;
1711  }
1712  par1 = a;
1713  par2 = b;
1714  }
1715  else if( p.nMatch("DENS") )
1716  {
1717  /* depends on density */
1718  chHextraScale = "DENSITY";
1719  hextra.lgHextraDensity = true;
1720 
1721  /* the log of the density */
1722  realnum a = (realnum)p.FFmtRead();
1723  /* if number not entered then unity is used, heating
1724  * proportional to density */
1726  par1 = a;
1727  par2 = 0;
1728  nPar = 2;
1729  }
1730  else if( p.nMatch("SS") )
1731  {
1732  /* alpha disk model, specify alpha, mass of black hole, and distance */
1733  chHextraScale = "SS";
1734  hextra.lgHextraSS = true;
1735 
1736  /* the parameter alpha of alpha-disk model */
1738 
1739  /* mass in solar masses of center black hole */
1740  realnum a = (realnum)p.FFmtRead();
1741  if( p.lgEOL() )
1742  p.NoNumb("hextraSS Mass");
1743  hextra.HextraSS_M = a*SOLAR_MASS;
1744 
1745  /* radius (cm) from center */
1746  realnum b = (realnum)p.FFmtRead();
1747  if( p.lgEOL() )
1748  p.NoNumb("hextraSS radius");
1749  hextra.HextraSSradius = b;
1750  par1 = a;
1751  par2 = 0;
1752  nPar = 2;
1753  }
1754  else
1755  {
1756  /* constant heating */
1757  chHextraScale = "";
1758  nPar = 1;
1759  }
1760 
1761  /* option to have heating vary with time */
1762  if( p.nMatch( "TIME" ) )
1763  hextra.lgTurbHeatVaryTime = true;
1764 
1765  /* vary option */
1766  if( optimize.lgVarOn )
1767  {
1768  if( hextra.lgHextraSS )
1769  {
1770  fprintf(ioQQQ,"Sorry, HEXTRA SS command does not now support vary option.\n");
1772  }
1773  /* 1 to 3 options on hextra vary */
1774  optimize.nvarxt[optimize.nparm] = nPar;
1775  optimize.vparm[0][optimize.nparm] = log10(hextra.TurbHeat);
1776  optimize.vparm[1][optimize.nparm] = par1;
1777  optimize.vparm[2][optimize.nparm] = par2;
1778 
1779  /* the keyword LOG is not used above, but is checked elsewhere */
1780  strcpy( optimize.chVarFmt[optimize.nparm], "HEXTra %f LOG " );
1781  strcat( optimize.chVarFmt[optimize.nparm], chHextraScale );
1782  while( nPar > 1 )
1783  {
1784  /* add extra spots to write parameters */
1785  --nPar;
1786  strcat( optimize.chVarFmt[optimize.nparm], " %f" );
1787  }
1789  strcat( optimize.chVarFmt[optimize.nparm], " TIME" );
1790  /* pointer to where to write */
1792  /* the increment in the first steps away from the original value */
1793  optimize.vincr[optimize.nparm] = 0.1f;
1794  ++optimize.nparm;
1795  }
1796 }
1797 
1799 {
1800  thermal.lgTeHigh = true;
1801 }
1803 {
1804  DEBUG_ENTRY( "ParseHydrogen()" );
1805  fprintf(ioQQQ," Sorry, this command has been replaced with the SPECIES H-LIKE command.\n");
1807 }
1809 {
1810  DEBUG_ENTRY( "ParseInitCount()" );
1811  /* read cloudy.ini initialization file
1812  * need to read in file every time, since some vars
1813  * get reset in zero - would be unsafe not to read in again */
1814  ParseInit(p);
1815 
1816  /* check that only one init file was in the input stream -
1817  * we cannot now read more than one */
1818  ++p.m_nInitFile;
1819  if( p.m_nInitFile > 1 )
1820  {
1821  fprintf( ioQQQ,
1822  " This is the second init file, I can only handle one.\nSorry.\n" );
1824  }
1825 
1826  /* we will put the data for the ini file at the end of the array of
1827  * line images, this is the increment for stuffing the lines in - negative */
1828  input.iReadWay = -1;
1829 
1830  /* initialize array reader, this sub does nothing but set
1831  * initial value of a variable, depending on value of iReadWay */
1832  input.init();
1833 }
1835 {
1836  DEBUG_ENTRY( "ParseIntensity()" );
1837  /* intensity of this continuum source */
1838  if( (p.m_nqh) >= LIMSPC )
1839  {
1840  /* too many continua were entered */
1841  fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
1843  }
1844 
1845  /* silly, but calms down the lint */
1846  ASSERT( (p.m_nqh) < LIMSPC );
1847  strcpy( rfield.chRSpec[(p.m_nqh)], "SQCM" );
1848  rfield.totpow[(p.m_nqh)] = p.FFmtRead();
1849  if( p.lgEOL() )
1850  p.NoNumb("intensity");
1851 
1852  /* set radius to very large value if not already set */
1853  /* >>chng 01 jul 24, from Radius == 0 to this, as per PvH comments */
1854  if( !radius.lgRadiusKnown )
1855  {
1856  /* set radius to large value */
1858  }
1859  if( p.nMatch("LINE") )
1860  {
1861  /* silly, but calms down the lint */
1862  ASSERT( (p.m_nqh) < LIMSPC );
1863  /* option for linear input parameter */
1864  rfield.totpow[(p.m_nqh)] = log10(rfield.totpow[(p.m_nqh)]);
1865  }
1866  strcpy( rfield.chSpNorm[(p.m_nqh)], "LUMI" );
1867  /* ParseRangeOption in readsun */
1868  ParseRangeOption(p);
1869 
1870  /* >>chng 06 mar 22, add time option to vary only some continua with time */
1871  if( p.nMatch( "TIME" ) )
1872  rfield.lgTimeVary[(p.m_nqh)] = true;
1873 
1874  /* vary option */
1875  if( optimize.lgVarOn )
1876  {
1877  /* create the format we will use to write the command */
1878  strcpy( optimize.chVarFmt[optimize.nparm], "INTENSITY %f LOG range %f %f" );
1879  if( rfield.lgTimeVary[(p.m_nqh)] )
1880  strcat( optimize.chVarFmt[optimize.nparm], " TIME" );
1881  /* array index for this command */
1884  optimize.vincr[optimize.nparm] = 0.5;
1885  /* range option, but cannot be varied */
1886  optimize.vparm[1][optimize.nparm] = (realnum)log10(rfield.range[(p.m_nqh)][0]);
1887  optimize.vparm[2][optimize.nparm] = (realnum)log10(rfield.range[(p.m_nqh)][1]);
1889  ++optimize.nparm;
1890  }
1891  ++(p.m_nqh);
1892 }
1894 {
1895  /* iterate to get optical depths and diffuse field properly */
1896  iterations.itermx = (long int)p.FFmtRead() - 1;
1898  /* >>chng 06 mar 19, malloc itrdim arrays
1899  * iterations.iter_malloc is -1 before space allocated and set to
1900  * itermx if not previously set */
1902  {
1903  long int iter_malloc_save = iterations.iter_malloc;
1904  /* add one for mismatch between array dim and iterations defn,
1905  * another few for safety */
1907  iterations.IterPrnt.resize((size_t)iterations.iter_malloc);
1908  iterations.nend.resize((size_t)iterations.iter_malloc);
1910  iterations.StopRadius.resize((size_t)iterations.iter_malloc);
1911  /* >>chng 06 jun 07, need to set IterPrnt, thickness, and nend */
1912  for(long int j=iter_malloc_save; j<iterations.iter_malloc; ++j )
1913  {
1914  iterations.IterPrnt[j] = iterations.IterPrnt[iter_malloc_save-1];
1915  iterations.nend[j] = iterations.nend[iter_malloc_save-1];
1916  iterations.StopThickness[j] = iterations.StopThickness[iter_malloc_save-1];
1917  iterations.StopRadius[j] = iterations.StopRadius[iter_malloc_save-1];
1918  }
1919  }
1920 
1921  if( p.nMatch("CONV") )
1922  {
1923  /* option to keep iterating until it converges, checks on convergence
1924  * are done in update, and checked again in prtcomment */
1925  conv.lgAutoIt = true;
1926  /* above would have been legitimate setting of ITERMX, set to default 10 */
1927  if( p.lgEOL() )
1928  {
1929  iterations.itermx = 10 - 1;
1930  }
1931  realnum a = (realnum)p.FFmtRead();
1932  /* change convergence criteria, preset in zero */
1933  if( !p.lgEOL() )
1934  {
1935  conv.autocv = a;
1936  }
1937  if( p.nMatch("ALL") )
1938  {
1939  conv.lgAllTransitions = true;
1940  }
1941  }
1942 }
1944 {
1945  /* this is the specific luminosity at nu
1946  * following says n(nu) not nuF(nu) */
1947  bool lgNu2 = false;
1948  /* in reads2 */
1949  ParseF_nu(p,"4 PI",lgNu2);
1950 }
1952 {
1953  DEBUG_ENTRY( "ParseLaser()" );
1954  /* mostly one frequency (a laser) to check gamma's,*/
1955  /* say the continuum type */
1956  strcpy( rfield.chSpType[rfield.nShape], "LASER" );
1957 
1958  /* scan off the laser's energy ion Rydbergs */
1960 
1961  /* negative energies are logs */
1962  if( rfield.slope[rfield.nShape] <= 0. )
1963  {
1966  }
1967  if( p.lgEOL() )
1968  {
1969  p.NoNumb("energy");
1970  }
1971 
1972  /* next number is optional relative width of laser */
1973  rfield.cutoff[rfield.nShape][0] = p.FFmtRead();
1974  if( p.lgEOL() )
1975  {
1976  /* default width is 0.05 */
1977  rfield.cutoff[rfield.nShape][0] = 0.05;
1978  }
1979 
1980  /* increment counter making sure we still have space enough */
1981  ++rfield.nShape;
1982  if( rfield.nShape >= LIMSPC )
1983  {
1984  /* too many continua were entered */
1985  fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
1987  }
1988 }
1990 {
1991  DEBUG_ENTRY( "ParseLuminosity()" );
1992  /* luminosity of this continuum source */
1993  if( (p.m_nqh) >= LIMSPC )
1994  {
1995  /* too many continua were entered */
1996  fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
1998  }
1999 
2000  strcpy( rfield.chRSpec[(p.m_nqh)], "4 PI" );
2001  rfield.totpow[(p.m_nqh)] = p.FFmtRead();
2002  if( p.lgEOL() )
2003  {
2004  p.NoNumb("luminosity");
2005  }
2006  if( p.nMatch("LINE") )
2007  {
2008  /* option for linear input parameter */
2009  rfield.totpow[(p.m_nqh)] = log10(rfield.totpow[(p.m_nqh)]);
2010  }
2011  if( rfield.totpow[(p.m_nqh)] > 200. )
2012  {
2013  fprintf( ioQQQ, " The log of the luminosity is very large: %g\n", rfield.totpow[(p.m_nqh)] );
2014  fprintf( ioQQQ, " Did you omit the keyword LINEAR?\n" );
2016  }
2017 
2018  strcpy( rfield.chSpNorm[(p.m_nqh)], "LUMI" );
2019 
2020  /* the solar option - number is log total solar luminosity */
2021  if( p.nMatch("SOLA") )
2022  {
2023  /* option to use log of total luminosity in solar units */
2024  rfield.range[(p.m_nqh)][0] = rfield.emm();
2025  rfield.range[(p.m_nqh)][1] = rfield.egamry();
2026  /* log of solar luminosity in erg s-1 */
2027  rfield.totpow[(p.m_nqh)] += 33.5827f;
2028  }
2029  else
2030  {
2031  /* check if range is present and parse it if it is - ParseRangeOption in readsun
2032  * this includes TOTAL keyword for total luminosity */
2033  ParseRangeOption(p);
2034  }
2035 
2036  /* >>chng 06 mar 22, add time option to vary only some continua with time */
2037  if( p.nMatch( "TIME" ) )
2038  rfield.lgTimeVary[(p.m_nqh)] = true;
2039 
2040  /* vary option */
2041  if( optimize.lgVarOn )
2042  {
2043  strcpy( optimize.chVarFmt[optimize.nparm], "LUMINOSITY %f LOG range %f %f" );
2044  if( rfield.lgTimeVary[(p.m_nqh)] )
2045  strcat( optimize.chVarFmt[optimize.nparm], " TIME" );
2046  /* pointer to where to write */
2049  optimize.vincr[optimize.nparm] = 0.5;
2050  /* range option in, but cannot be varied */
2051  optimize.vparm[1][optimize.nparm] = (realnum)log10(rfield.range[(p.m_nqh)][0]);
2052  optimize.vparm[2][optimize.nparm] = (realnum)log10(rfield.range[(p.m_nqh)][1]);
2054  ++optimize.nparm;
2055  }
2056  ++(p.m_nqh);
2057 }
2059 {
2060  /* heating and ionization due to fast neutrons */
2061  hextra.lgNeutrnHeatOn = true;
2062 
2063  /* first number is neutron luminosity
2064  * relative to bolometric luminosity */
2065  hextra.frcneu = (realnum)p.FFmtRead();
2066  if( p.lgEOL() )
2067  p.NoNumb("neutron luminosity");
2068 
2069  /* save as log of fraction */
2070  if( hextra.frcneu > 0. )
2071  hextra.frcneu = (realnum)log10(hextra.frcneu);
2072 
2073  /* second number is efficiency */
2074  hextra.effneu = (realnum)p.FFmtRead();
2075  if( p.lgEOL() )
2076  {
2077  hextra.effneu = 1.0;
2078  }
2079  else
2080  {
2081  if( hextra.effneu <= 0. )
2083  }
2084 }
2086 {
2087  /* flux density of this continuum source, at optional frequency
2088  * expressed as product nu*f_nu */
2089  bool lgNu2 = true;
2090  /* in reads2 */
2091  ParseF_nu(p,"SQCM",lgNu2);
2092 }
2094 {
2095  /* specific luminosity density of this continuum source, at opt freq
2096  * expressed as product nu*f_nu */
2097  bool lgNu2 = true;
2098  /* in reads2 */
2099  ParseF_nu(p,"4 PI",lgNu2);
2100 }
2102 {
2103  DEBUG_ENTRY( "ParsePhi()" );
2104  /* enter phi(h), the number of h-ionizing photons/cm2 */
2105  if( (p.m_nqh) >= LIMSPC )
2106  {
2107  /* too many continua were entered */
2108  fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
2110  }
2111  /* silly, but calms down the lint */
2112  ASSERT( (p.m_nqh) < LIMSPC );
2113  strcpy( rfield.chRSpec[(p.m_nqh)], "SQCM" );
2114  strcpy( rfield.chSpNorm[(p.m_nqh)], "PHI " );
2115  rfield.totpow[(p.m_nqh)] = p.FFmtRead();
2116  if( p.lgEOL() )
2117  {
2118  p.NoNumb("number of h-ionizing photons");
2119  }
2120  /* set radius to very large value if not already set */
2121  /* >>chng 01 jul 24, from Radius == 0 to this, as per PvH comments */
2122  if( !radius.lgRadiusKnown )
2123  {
2124  /* set radius to large value */
2126  }
2127  /* check if continuum so intense that crash is likely */
2128  if( rfield.totpow[(p.m_nqh)] > 29. )
2129  {
2130  fprintf( ioQQQ, " Is the flux for this continuum correct?\n" );
2131  fprintf( ioQQQ, " It appears too bright to me.\n" );
2132  }
2133  /* ParseRangeOption in readsun xx parse_rangeoption*/
2134  ParseRangeOption(p);
2135 
2136  /* >>chng 06 mar 22, add time option to vary only some continua with time */
2137  if( p.nMatch( "TIME" ) )
2138  rfield.lgTimeVary[(p.m_nqh)] = true;
2139 
2140  /* vary option */
2141  if( optimize.lgVarOn )
2142  {
2143  strcpy( optimize.chVarFmt[optimize.nparm], "phi(h) %f LOG range %f %f" );
2144  if( rfield.lgTimeVary[(p.m_nqh)] )
2145  strcat( optimize.chVarFmt[optimize.nparm], " TIME" );
2146  /* pointer to where to write */
2149  optimize.vincr[optimize.nparm] = 0.5;
2150  /* range option in, but cannot be varied */
2151  optimize.vparm[1][optimize.nparm] = (realnum)log10(rfield.range[(p.m_nqh)][0]);
2152  optimize.vparm[2][optimize.nparm] = (realnum)log10(rfield.range[(p.m_nqh)][1]);
2154  ++optimize.nparm;
2155  }
2156  ++(p.m_nqh);
2157 }
2158 void ParseQH(Parser &p)
2159 {
2160  DEBUG_ENTRY( "ParseQH()" );
2161  /* log of number of ionizing photons */
2162  if( (p.m_nqh) >= LIMSPC )
2163  {
2164  /* too many continua were entered */
2165  fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
2167  }
2168 
2169  /* silly, but calms down the lint */
2170  ASSERT( (p.m_nqh) < LIMSPC );
2171  strcpy( rfield.chRSpec[(p.m_nqh)], "4 PI" );
2172  strcpy( rfield.chSpNorm[(p.m_nqh)], "Q(H)" );
2173  rfield.totpow[(p.m_nqh)] = p.FFmtRead();
2174  if( rfield.totpow[(p.m_nqh)] > 100. && called.lgTalk )
2175  {
2176  fprintf( ioQQQ, " Is this reasonable?\n" );
2177  }
2178  if( p.lgEOL() )
2179  {
2180  p.NoNumb("number of ionizing photons");
2181  }
2182  /* ParseRangeOption in readsun */
2183  ParseRangeOption(p);
2184 
2185  /* >>chng 06 mar 22, add time option to vary only some continua with time */
2186  if( p.nMatch( "TIME" ) )
2187  rfield.lgTimeVary[(p.m_nqh)] = true;
2188 
2189  /* vary option */
2190  if( optimize.lgVarOn )
2191  {
2192  strcpy( optimize.chVarFmt[optimize.nparm], "Q(H) %f LOG range %f %f" );
2193  if( rfield.lgTimeVary[(p.m_nqh)] )
2194  strcat( optimize.chVarFmt[optimize.nparm], " TIME" );
2195  /* pointer to where to write */
2198  optimize.vincr[optimize.nparm] = 0.5;
2199  /* range option in, but cannot be varied */
2200  optimize.vparm[1][optimize.nparm] = (realnum)log10(rfield.range[(p.m_nqh)][0]);
2201  optimize.vparm[2][optimize.nparm] = (realnum)log10(rfield.range[(p.m_nqh)][1]);
2203  ++optimize.nparm;
2204  }
2205  /* increment number of luminosity sources specified */
2206  ++(p.m_nqh);
2207 }
2209 {
2210  /* this is the Roberto Terlivich command, no telling if it still works */
2211  radius.dRadSign = -1.;
2212 }
2214 {
2215  DEBUG_ENTRY( "ParseSpecial()" );
2216  /* special key, can do anything */
2218 }
2220 {
2221  /* taumin command minimum optical depths for lines dafault 1e-20 */
2222  opac.taumin = (realnum)exp10(p.FFmtRead());
2223  if( p.lgEOL() )
2224  p.NoNumb("minimum optical depth");
2225 }
2227 {
2228  /* read in title of model starting in col 5 -- prefer to get string
2229  in quotes, but for the moment if it's not present take what you
2230  can get */
2231  string title;
2232  if (p.GetQuote(title) != 0)
2233  strcpy( input.chTitle , p.getRawTail().c_str() );
2234  else
2235  strcpy( input.chTitle , title.c_str() );
2236 }
2238 {
2239  DEBUG_ENTRY( "ParseTolerance()" );
2240  fprintf(ioQQQ,
2241  "Sorry, this command has been replaced with the SET TEMPERATURE TOLERANCE command.\n");
2243 }
2245 {
2246  /* velocity power law as a function of radius. */
2248 
2249  DoppVel.lgTurbLawOn = true;
2251  /* for now, insist upon non-positive power,
2252  * so that velocity decreases with increasing radius. */
2253  ASSERT( DoppVel.TurbVelLaw <= 0.f );
2254 }
2256 {
2257  DEBUG_ENTRY( "ParseTurbulence()" );
2258  string ExtraPars;
2259 
2260  if( p.nMatch("EQUIPART") )
2261  {
2262  /* turbulence equipartition option */
2263  DoppVel.lgTurbEquiMag = true;
2264 
2265  /* this is optional F parameter from equation 34 of
2266  *>>refer pressure turb Heiles, C. & Troland, T.H. 2005, ApJ, 624, 773
2267  * turbulent energy density will be rho F v^2 / 2 */
2269  if( p.lgEOL() )
2270  {
2271  /* this is the default value of 3 for isotropic turbulence */
2272  DoppVel.Heiles_Troland_F = 3.f;
2273  }
2274  }
2275  else
2276  {
2277  /* line has turbulent velocity in km/s */
2278  DoppVel.lgTurbEquiMag = false;
2280  if( p.lgEOL() )
2281  p.NoNumb("microturbulent velocity");
2282 
2283  if( p.nMatch(" LOG") )
2284  {
2285  if( DoppVel.TurbVel > 32. )
2286  {
2287  fprintf( ioQQQ, "PROBLEM the log of the turbulence is "
2288  "%.2e - I cannot handle a number this big.\n",
2289  DoppVel.TurbVel );
2290  fprintf( ioQQQ, " The line image was\n" );
2291  p.PrintLine(ioQQQ);
2292  fprintf( ioQQQ, " Sorry.\n" );
2294  }
2296  }
2297  /* now convert from km/s to cm/s */
2298  DoppVel.TurbVel *= 1e5;
2299 
2300  if( DoppVel.TurbVel < 0. )
2301  {
2302  fprintf( ioQQQ, " PROBLEM: the turbulent velocity needs to be > 0, but this was entered: %e\n",
2303  DoppVel.TurbVel );
2304  fprintf( ioQQQ, " Bailing out. Sorry.\n" );
2306  }
2307  else if( DoppVel.TurbVel >= SPEEDLIGHT )
2308  {
2309  fprintf( ioQQQ, " PROBLEM: A turbulent velocity greater than speed of light is not allowed, this was entered: %e\n",
2310  DoppVel.TurbVel );
2311  fprintf( ioQQQ, " Bailing out. Sorry.\n" );
2313  }
2314 
2315  /* this is optional F parameter from equation 34 of
2316  *>>refer pressure turb Heiles, C. & Troland, T.H. 2005, ApJ, 624, 773
2317  * turbulent energy density will be rho F v^2 / 2 */
2319  if( p.lgEOL() )
2320  {
2321  /* this is the default value of 3 for isotropic turbulence */
2322  DoppVel.Heiles_Troland_F = 3.f;
2323  }
2324 
2325  /* option to have turbulence be dissipative, keyword is dissipate,
2326  * argument is log of scale length in cm */
2327  if( p.nMatch("DISS") )
2328  {
2330  if( p.lgEOL() )
2331  p.NoNumb("turbulence dissipation scale");
2332  ExtraPars += " DISSIPATE %f";
2333  }
2334  }
2335 
2336  /* include turbulent pressure in equation of state?
2337  * >>chng 06 mar 24, include turbulent pressure in gas equation of state by default,
2338  * but not if NO PRESSURE occurs */
2339  if( p.nMatch(" NO ") && p.nMatch("PRES") )
2340  {
2341  DoppVel.lgTurb_pressure = false;
2342  ExtraPars += " NO PRESSURE";
2343  }
2344  else
2345  {
2346  /* default is to include turbulent pressure in equation of state */
2347  DoppVel.lgTurb_pressure = true;
2348  }
2349 
2350  /* vary option */
2351  if( optimize.lgVarOn && !p.nMatch("EQUIPART") )
2352  {
2353  /* only one parameter to vary */
2355  // the line image
2356  strcpy( optimize.chVarFmt[optimize.nparm], "TURBULENCE= %f LOG %f" );
2357  strcat( optimize.chVarFmt[optimize.nparm], ExtraPars.c_str() );
2358 
2359  /* pointer to where to write */
2361  /* turbulent velocity */
2362  optimize.vparm[0][optimize.nparm] = log10(DoppVel.TurbVel/1e5f);
2364  if( p.nMatch("DISS") )
2365  {
2368  }
2369  optimize.vincr[optimize.nparm] = 0.1f;
2370  ++optimize.nparm;
2371  }
2372 
2373  /* set initial turbulence */
2375 }
void ParseState(Parser &p)
Definition: parse_state.cpp:9
long int nSave
Definition: input.h:62
realnum col_h2
Definition: stopcalc.h:74
void ParseF_nuSpecific(Parser &p)
void ParseL_nu(Parser &p)
long int lim_iter
Definition: iterations.h:62
bool lgStoutLevelsSet
Definition: atmdat.h:415
realnum cextpw
Definition: thermal.h:152
double emm() const
Definition: mesh.h:84
bool nMatch(const char *chKey) const
Definition: parser.h:140
void ParseHDEN(Parser &p)
Definition: parse_hden.cpp:10
t_fudgec fudgec
Definition: fudgec.cpp:5
void ParseAperture(Parser &p)
t_mole_global mole_global
Definition: mole.cpp:7
long int iter_malloc
Definition: iterations.h:40
char chLamdaFile[FILENAME_PATH_LENGTH]
Definition: atmdat.h:395
void ParseAgn(Parser &p)
Definition: parse_agn.cpp:10
double Radius
Definition: radius.h:31
realnum colnut
Definition: stopcalc.h:69
double depth
Definition: radius.h:31
void ParseHydrogen(Parser &)
realnum GetDensity(realnum z)
Definition: cosmology.cpp:28
t_atmdat atmdat
Definition: atmdat.cpp:6
void InitMonitorResults(void)
t_thermal thermal
Definition: thermal.cpp:6
void SetGasPhaseDensity(const long nelem, const realnum density)
Definition: dense.cpp:106
void ParseDont(Parser &p)
Definition: parse_dont.cpp:28
bool lgTurbLawOn
Definition: doppvel.h:33
bool lgStoutHybrid
Definition: atmdat.h:404
void ParseInitCount(Parser &p)
realnum DensityPower
Definition: dense.h:250
bool lgChiantiLevelsSet
Definition: atmdat.h:388
realnum size
Definition: geometry.h:77
void echo(void) const
Definition: parser.cpp:167
double FFmtRead(void)
Definition: parser.cpp:438
double exp10(double x)
Definition: cddefines.h:1383
void ParseDynaWind(Parser &p)
Definition: dynamics.cpp:1816
const int ipHE_LIKE
Definition: iso.h:65
void ParseAge(Parser &p)
Definition: parse_age.cpp:38
char chStoutFile[FILENAME_PATH_LENGTH]
Definition: atmdat.h:408
void ParsePlot(Parser &p)
Definition: parse_plot.cpp:16
void ParseDatabaseISO(long ipISO, Parser &p)
t_input input
Definition: input.cpp:12
void ParseCylinder(Parser &p)
bool lgGrid
Definition: grid.h:42
void ParseBlackbody(Parser &p)
t_opac opac
Definition: opacity.cpp:5
realnum Heiles_Troland_F
Definition: doppvel.h:37
void setline(const char *const card)
Definition: parser.h:79
realnum turrad
Definition: hextra.h:51
long int nvfpnt[LIMPAR]
Definition: optimize.h:198
void ParseBackgrd(Parser &p)
realnum turback
Definition: hextra.h:51
bool Command(const char *name, OptionParser doOpts)
Definition: parser.h:183
bool lgStoutOn
Definition: atmdat.h:402
long int m_nqh
Definition: parser.h:52
t_isoCTRL iso_ctrl
Definition: iso.cpp:9
t_cpu_i & i()
Definition: cpu.h:415
void ParseFail(Parser &p)
double totpow[LIMSPC]
Definition: rfield.h:286
void ParseEnergy(Parser &p)
bool lgCylnOn
Definition: radius.h:126
char chDffTrns[4]
Definition: rfield.h:219
vector< long int > nend
Definition: iterations.h:71
void ParseDatabaseH2(Parser &p)
realnum windv0
Definition: wind.h:11
char chRSpec[LIMSPC][5]
Definition: rfield.h:333
void ParseMonitorResults(Parser &p)
void ParseNeutrons(Parser &p)
double CylindHigh
Definition: radius.h:125
void ParseExtinguish(Parser &p)
realnum redshift_start
Definition: cosmology.h:26
long int iEmissPower
Definition: geometry.h:71
void ParseLuminosity(Parser &p)
void ParseTable(Parser &p)
Definition: parse_table.cpp:94
void ParseSet(Parser &p)
Definition: parse_set.cpp:38
realnum HextraSSradius
Definition: hextra.h:73
bool lgTrOpt
Definition: optimize.h:256
long nStoutMaxLevelsFe
Definition: atmdat.h:410
t_StopCalc StopCalc
Definition: stopcalc.cpp:7
t_conv conv
Definition: conv.cpp:5
long int nOptimiz
Definition: optimize.h:250
t_hextra hextra
Definition: hextra.cpp:5
bool lgNeutrnHeatOn
Definition: hextra.h:81
long nLamdaMaxLevels
Definition: atmdat.h:397
bool lgOpacityFine
Definition: rfield.h:402
double distance
Definition: radius.h:70
int GetQuote(string &chLabel)
Definition: parser.cpp:184
realnum varang[LIMPAR][2]
Definition: optimize.h:201
realnum fiscal
Definition: geometry.h:29
long int nRead
Definition: input.h:62
t_phycon phycon
Definition: phycon.cpp:6
realnum effneu
Definition: hextra.h:85
realnum TurbVel
Definition: doppvel.h:21
void ParseFill(Parser &p)
bool isCommandComment(void) const
Definition: parser.cpp:117
void ParseRatio(Parser &p)
Definition: parse_ratio.cpp:10
void ParseHExtra(Parser &p)
void ParseConstant(Parser &p)
realnum TurbHeatSave
Definition: hextra.h:42
void ParseTitle(Parser &)
double getNumberDefaultAlwaysLog(const char *chDesc, double fdef)
Definition: parser.cpp:412
bool lgIonizReevaluate
Definition: rfield.h:112
bool lgTeHigh
Definition: thermal.h:72
void ParseHelp(Parser &)
bool lgChiantiPrint
Definition: atmdat.h:378
#define NFUDGC
Definition: fudgec.h:9
void ParseFudge(Parser &p)
void ParseCrashDo(Parser &p)
char chVarFmt[LIMPAR][FILENAME_PATH_LENGTH_2]
Definition: optimize.h:267
vector< double > StopThickness
Definition: iterations.h:77
realnum covgeo
Definition: geometry.h:45
FILE * ioQQQ
Definition: cddefines.cpp:7
void ParseStop(Parser &p)
Definition: parse_stop.cpp:17
realnum vparm[LIMEXT][LIMPAR]
Definition: optimize.h:192
realnum FillFac
Definition: geometry.h:29
bool lgAllTransitions
Definition: conv.h:254
char chTitle[INPUT_LINE_LENGTH]
Definition: input.h:48
bool lgCExtraOn
Definition: thermal.h:151
bool lgTalk
Definition: called.h:12
void ParseCosm(Parser &p)
void ParseSave(Parser &p)
Definition: parse_save.cpp:73
void ParseGrid(Parser &p)
Definition: parse_grid.cpp:10
t_DoppVel DoppVel
Definition: doppvel.cpp:5
void ParseDoubleTau(Parser &)
bool lgConverge_set
Definition: iterations.h:60
const char * name
Definition: parser.h:27
void ParseDLaw(Parser &p)
Definition: parse_dlaw.cpp:10
#define MIN2(a, b)
Definition: cddefines.h:807
Definition: parser.h:42
bool lgTurbEquiMag
Definition: doppvel.h:46
realnum colpls
Definition: stopcalc.h:69
realnum TurbVelZero
Definition: doppvel.h:25
void ParseRangeOption(Parser &p)
double cutoff[LIMSPC][3]
Definition: rfield.h:286
realnum frcneu
Definition: hextra.h:83
void ParseVLaw(Parser &p)
void ParsePrint(Parser &p)
bool lgVarOn
Definition: optimize.h:207
bool lgPhysOK
Definition: phycon.h:111
t_dense dense
Definition: global.cpp:15
void ParseSphere(Parser &p)
Definition: parse_sphere.cpp:9
void ParseTurbulence(Parser &p)
bool lgTimeVary[LIMSPC]
Definition: rfield.h:292
realnum autocv
Definition: conv.h:260
void ParseMap(Parser &p)
Definition: parse_map.cpp:9
static t_version & Inst()
Definition: cddefines.h:209
double range[LIMSPC][2]
Definition: rfield.h:329
void ParseIlluminate(Parser &p)
long int m_nInitFile
Definition: parser.h:52
void ParseGravity(Parser &p)
const int LIMSPC
Definition: rfield.h:20
void ParseQH(Parser &p)
void ParseCommands(void)
Wind wind
Definition: wind.cpp:5
t_trace trace
Definition: trace.cpp:5
bool lgAbnSolar
Definition: abund.h:202
vector< long int > IterPrnt
Definition: iterations.h:43
void ParseCosmology(Parser &p)
bool lgBallistic(void) const
Definition: wind.h:31
OptionParser action
Definition: parser.h:28
void ParseCompile(Parser &p)
double rinner
Definition: radius.h:31
t_abund abund
Definition: abund.cpp:5
t_geometry geometry
Definition: geometry.cpp:5
void ParseTest(Parser &p)
Definition: parse_test.cpp:9
bool lgOutOnly
Definition: rfield.h:224
const double TEMP_STOP_DEFAULT
Definition: phycon.h:119
bool lgChiantiHybrid
Definition: atmdat.h:376
t_dark_matter dark
Definition: dark_matter.cpp:5
bool lgMap
Definition: conv.h:233
realnum CoolExtra
Definition: thermal.h:152
void ParseSpecies(Parser &p)
realnum HColStop
Definition: stopcalc.h:69
realnum redshift_current
Definition: cosmology.h:26
double slope[LIMSPC]
Definition: rfield.h:286
bool lgLamdaPrint
Definition: atmdat.h:393
bool lgSaveXspec
Definition: grid.h:39
void ParseNuL_nu(Parser &p)
void ParseHeLike(Parser &)
void ParseGlobule(Parser &p)
void ParseCaseB(Parser &p)
Definition: parse_caseb.cpp:9
void ParseTolerance(Parser &)
string getRawTail()
Definition: parser.h:211
void ParseMagnet(Parser &p)
Definition: magnetic.cpp:156
bool lgChiantiOn
Definition: atmdat.h:374
bool lgTrace
Definition: trace.h:12
void ParsePGrains(Parser &)
bool m_lgDSet
Definition: parser.h:53
void set_point(long int ipnt)
Definition: parser.h:87
void ParseDiffuse(Parser &p)
bool lgEnabled
Definition: h2_priv.h:352
void ParseF_nu(Parser &p, const char *chType, bool lgNU2)
Definition: parse_f_nu.cpp:9
void ParseIonParX(Parser &p)
char chCloudyChiantiFile[FILENAME_PATH_LENGTH]
Definition: atmdat.h:382
realnum TurbHeat
Definition: hextra.h:42
long int nparm
Definition: optimize.h:204
t_plotCom plotCom
Definition: plot.cpp:20
t_pressure pressure
Definition: pressure.cpp:9
t_rfield rfield
Definition: rfield.cpp:9
long nStoutMaxLevels
Definition: atmdat.h:413
realnum TurbVelLaw
Definition: doppvel.h:29
bool lgCaseB
Definition: opacity.h:173
void ParseAbsMag(Parser &p)
Definition: parse_absmag.cpp:9
realnum DispScale
Definition: doppvel.h:51
long int nTrOpt
Definition: optimize.h:255
float realnum
Definition: cddefines.h:124
#define EXIT_FAILURE
Definition: cddefines.h:168
bool lgPrintNumberOfLevels
Definition: iso.h:346
bool last(void) const
Definition: parser.h:196
const realnum COLUMN_INIT
Definition: stopcalc.h:14
const int INPUT_LINE_LENGTH
Definition: cddefines.h:301
realnum HextraScaleDensity
Definition: hextra.h:60
void ParseFluc(Parser &p)
Definition: parse_fluc.cpp:9
void ParseNorm(Parser &p)
Definition: parse_norm.cpp:10
bool lgTurbHeatVaryTime
Definition: hextra.h:76
#define cdEXIT(FAIL)
Definition: cddefines.h:484
DepthTable DLW
Definition: dense.h:198
bool lgMPI_talk() const
Definition: cpu.h:394
void ParseCMB(double z, long int *nqh)
Definition: parse_CMB.cpp:10
double tabval(double r0, double depth) const
Definition: depth_table.cpp:8
NORETURN void NoNumb(const char *chDesc) const
Definition: parser.cpp:318
bool lgRadiusKnown
Definition: radius.h:121
realnum covrt
Definition: geometry.h:61
long min(int a, long b)
Definition: cddefines.h:766
void ParsePhi(Parser &p)
bool lgMustBlockHIon
Definition: rfield.h:94
diatomics h2("h2", 4100.,&hmi.H2_total, Yan_H2_CS)
const long LIMPAR
Definition: optimize.h:61
long nChiantiMaxLevels
Definition: atmdat.h:386
int gravity_symmetry
Definition: pressure.h:84
t_iterations iterations
Definition: iterations.cpp:6
void ParseLaser(Parser &p)
void ParseDynaTime(Parser &p)
Definition: dynamics.cpp:1669
void ParseConvHighT(Parser &)
void ParseIntensity(Parser &p)
t_optimize optimize
Definition: optimize.cpp:6
void ParseDrive(Parser &p)
Definition: parse_drive.cpp:31
char chSpNorm[LIMSPC][5]
Definition: rfield.h:333
t_grid grid
Definition: grid.cpp:5
realnum cryden
Definition: hextra.h:24
bool isVar(void) const
Definition: parser.cpp:122
t_radius radius
Definition: radius.cpp:5
double dRadSign
Definition: radius.h:73
realnum vincr[LIMPAR]
Definition: optimize.h:195
realnum TempLoStopZone
Definition: stopcalc.h:42
long int lim_zone
Definition: iterations.h:61
void ParseTrace(Parser &p)
Definition: parse_trace.cpp:11
t_prt prt
Definition: prt.cpp:14
void ParseTauMin(Parser &p)
bool lgTrOvrd
Definition: trace.h:121
bool lgDoLineTrans
Definition: rfield.h:101
void ParseInit(Parser &p)
Definition: parse_init.cpp:9
void ParseGrain(Parser &p)
Definition: parse_grain.cpp:12
long int nfudge
Definition: fudgec.h:17
bool lgHextraSS
Definition: hextra.h:64
realnum gas_phase[LIMELM]
Definition: dense.h:76
long int itermx
Definition: iterations.h:37
void ParsePowerlawContinuum(Parser &p)
double HextraSS_M
Definition: hextra.h:70
bool lgOptimizeAsLinear[LIMPAR]
Definition: optimize.h:184
void help(FILE *fp) const
Definition: parser.cpp:267
long int LimFail
Definition: conv.h:230
bool lgNegativeIncrements
Definition: grid.h:38
void ParseIonParI(Parser &p)
#define ASSERT(exp)
Definition: cddefines.h:617
realnum covaper
Definition: geometry.h:54
long int iReadWay
Definition: input.h:62
char chDenseLaw[5]
Definition: dense.h:176
bool lgOpacityReevaluate
Definition: rfield.h:105
bool lgPrintHTML
Definition: prt.h:277
realnum ConstTemp
Definition: thermal.h:56
bool getline()
Definition: parser.cpp:249
const int ipH_LIKE
Definition: iso.h:64
realnum emdot
Definition: wind.h:39
bool lgHextraDensity
Definition: hextra.h:57
#define DEBUG_ENTRY(funcname)
Definition: cddefines.h:729
void ParseCovering(Parser &p)
t_cosmology cosmology
Definition: cosmology.cpp:8
void ParseSpecial(Parser &)
void ParseNuF_nu(Parser &p)
bool lgAutoIt
Definition: conv.h:251
double getNumberCheckAlwaysLog(const char *chDesc)
Definition: parser.cpp:393
#define MIN3(a, b, c)
Definition: cddefines.h:812
void ParseRadius(Parser &p)
realnum ConstGrainTemp
Definition: thermal.h:59
double egamry() const
Definition: mesh.h:88
void ParseTLaw(Parser &p)
Definition: parse_tlaw.cpp:13
void ParseRoberto(Parser &)
realnum tauend
Definition: stopcalc.h:23
bool lgEOL(void) const
Definition: parser.h:103
#define MAX2(a, b)
Definition: cddefines.h:828
void ParseIterations(Parser &p)
void ParseCosmicRays(Parser &p)
bool lgStoutPrint
Definition: atmdat.h:406
int fprintf(const Output &stream, const char *format,...)
Definition: service.cpp:1217
NORETURN void CommandError(void) const
Definition: parser.cpp:239
double dense_fabden(double radius, double depth)
bool lgTurb_pressure
Definition: doppvel.h:42
void ParseBremsstrahlung(Parser &p)
vector< double > external_mass[3]
Definition: pressure.h:86
realnum DoubleTau
Definition: rt.h:178
long int nShape
Definition: rfield.h:308
void ParseCMBOuter(Parser &p)
double dense_parametric_wind(double rad)
realnum HextraSSalpha
Definition: hextra.h:67
void ParseOptimize(Parser &p)
int PrintLine(FILE *fp) const
Definition: parser.h:200
bool lgStatic(void) const
Definition: wind.h:24
vector< double > StopRadius
Definition: iterations.h:80
double self_mass_factor
Definition: pressure.h:85
bool lgSizeSet
Definition: geometry.h:80
bool lgPlotON
Definition: plot.h:21
double rdfalt
Definition: radius.h:129
GrainVar gv
Definition: grainvar.cpp:5
bool lgLamdaOn
Definition: atmdat.h:391
static t_cpu cpu
Definition: cpu.h:423
bool m_lgEOF
Definition: parser.h:53
bool lgNoMole
Definition: mole.h:325
bool lgNoVary
Definition: optimize.h:175
void ParseMetal(Parser &p)
Definition: parse_metal.cpp:12
long nChiantiMaxLevelsFe
Definition: atmdat.h:384
const int ipHYDROGEN
Definition: cddefines.h:348
long int nGridCommands
Definition: grid.h:55
bool lgChiantiExp
Definition: atmdat.h:380
void doSetVar(void)
Definition: parser.cpp:139
void ParseCExtra(Parser &p)
bool lgBlockHIon
Definition: rfield.h:98
void ParseDielectronic(Parser &)
void ParseDatabase(Parser &p)
realnum EdenExtra
Definition: dense.h:217
realnum col_h2_nut
Definition: stopcalc.h:80
double r_200
Definition: dark_matter.h:19
void ParseCoronal(Parser &p)
long int nvarxt[LIMPAR]
Definition: optimize.h:198
char chSpType[LIMSPC][6]
Definition: rfield.h:333
void ParseElement(Parser &p)
realnum taumin
Definition: opacity.h:166
t_called called
Definition: called.cpp:4
void ParseForceTemperature(Parser &p)
bool lgLamdaLevelsSet
Definition: atmdat.h:399
long int nplot
Definition: plot.h:27
realnum filpow
Definition: geometry.h:29
void ParseInterp(Parser &p)
void ParseEden(Parser &p)
void ParseDistance(Parser &p)
void ParseDarkMatter(Parser &p)
realnum fudgea[NFUDGC]
Definition: fudgec.h:15
void init(void)
Definition: input.cpp:76
void ParseAbundances(Parser &p)
void ParseChemistry(Parser &p)
Definition: mole.cpp:144
bool lgDColOn
Definition: grainvar.h:496
bool lgHextraDepth
Definition: hextra.h:49
t_rt rt
Definition: rt.cpp:5