cloudy  trunk
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
parse_abundances.cpp
Go to the documentation of this file.
1 /* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and
2  * others. For conditions of distribution and use see copyright notice in license.txt */
3 /*ParseAbundances parse and read in composition as set by abundances command */
4 #include "cddefines.h"
5 #include "grains.h"
6 #include "abund.h"
7 #include "phycon.h"
8 #include "called.h"
9 #include "elementnames.h"
10 #include "input.h"
11 #include "parse.h"
12 
13 void ParseAbundances(char *chCard,
14  /* following set true by grains command,
15  * so this will not set more grains if grains already on. */
16  bool lgDSet)
17 {
18  bool lgEOF,
19  lgEOL,
20  lgLog;
21  long int i,
22  j,
23  k;
24  double absav[LIMELM],
25  chk;
26  GrainPar gp;
27 
28  DEBUG_ENTRY( "ParseAbundances()" );
29 
30  j = 5;
31  absav[0] = FFmtRead(chCard,&j,INPUT_LINE_LENGTH,&lgEOL);
32  /* check whether it scanned off end of line on first try
33  * at least one number on the line if it passes this test
34  * no numbers means a keyword */
35  if( !lgEOL )
36  {
37  absav[1] = FFmtRead(chCard,&j,INPUT_LINE_LENGTH,&lgEOL);
38  if( lgEOL )
39  {
40  /* this branch, we found one, but not two, numbers -
41  * must be one of a few special cases */
42  if( nMatch(" ALL",chCard) )
43  {
44  /* special option, all abundances will be this number */
45  if( absav[0] <= 0. )
46  {
47  absav[0] = pow(10.,absav[0]);
48  }
49  for( i=1; i < LIMELM; i++ )
50  {
51  abund.solar[i] = (realnum)absav[0];
52  }
53 
54  }
55  else if( nMatch("OLD ",chCard) && nMatch("SOLA",chCard) )
56  {
57  i = (int)absav[0];
58  /* old solar - better be the number "84" on the line */
59  if( i!=84 )
60  {
61  fprintf( ioQQQ,
62  " The only old abundance set I have is for version 84 - %3ld was requested. Sorry.\n",
63  i );
64  cdEXIT(EXIT_FAILURE);
65  }
66  for( i=1; i < LIMELM; i++ )
67  {
68  /* these are the old abundances */
70  abund.solar[i] = abund.OldSolar84[i];
71  }
72  }
73 
74  else if( nMatch("STAR",chCard) )
75  {
76  /* Fred Hamonn's star burst galaxy mixture */
77  abund_starburst(chCard);
78  }
79  else
80  {
81  fprintf( ioQQQ,
82  " I did not recognize a sub-keyword - options are ALL and OLD SOLAR 84. Sorry.\n");
83  cdEXIT(EXIT_FAILURE);
84  }
85 
86  /* normal return */
87  return;
88  }
89 
90  /* we get here if there is a second number - read in all abundances */
91  for( i=2; i < abund.npSolar; i++ )
92  {
93  absav[i] = FFmtRead(chCard,&j,INPUT_LINE_LENGTH,&lgEOL);
94  if( lgEOL )
95  {
96  /* read CONTINUE line if J scanned off end before reading all abund */
97  input_readarray(chCard,&lgEOF);
98  if( lgEOF )
99  {
100  fprintf( ioQQQ, " There MUST be%3ld abundances entered, there were only%3ld. Sorry.\n",
101  abund.npSolar, i );
102  cdEXIT(EXIT_FAILURE);
103  }
104 
105  /* option to ignore all lines starting with #, *, or % */
106  while( lgInputComment(chCard) )
107  {
108  if( lgEOF )
109  {
110  fprintf( ioQQQ, " There MUST be%3ld abundances entered, there were only%3ld. Sorry.\n",
111  abund.npSolar, i );
112  cdEXIT(EXIT_FAILURE);
113  }
114  /* read in another line */
115  input_readarray(chCard,&lgEOF);
116  }
117 
118  /* we have the line image, print it */
119  if( called.lgTalk )
120  {
121  fprintf( ioQQQ, " * ");
122  k=0;
123  while( chCard[k]!='\0' )
124  {
125  fprintf(ioQQQ,"%c",chCard[k]);
126  ++k;
127  }
128  while( k<80 )
129  {
130  fprintf(ioQQQ,"%c",' ');
131  ++k;
132  }
133  fprintf( ioQQQ,"*\n");
134  }
135 
136  /* now convert to caps */
137  caps(chCard);
138  if( strncmp(chCard,"CONT",4) != 0 )
139  {
140  fprintf( ioQQQ, " There MUST be%3ld abundances entered, there were only%3ld. Sorry.\n",
141  abund.npSolar, i );
142  cdEXIT(EXIT_FAILURE);
143  }
144  else
145  {
146  j = 1;
147  absav[i] = FFmtRead(chCard,&j,INPUT_LINE_LENGTH,&lgEOL);
148  if( lgEOF )
149  {
150  fprintf( ioQQQ, " There MUST be%3ld abundances entered, there were only%3ld. Sorry.\n",
151  abund.npSolar, i);
152  cdEXIT(EXIT_FAILURE);
153  }
154  }
155  }
156  }
157 
158  /* fell through to here after reading in N abundances for N elements
159  * check that there are no more abundances on the line - that would
160  * be an error - a typo */
161  chk = FFmtRead(chCard,&j,INPUT_LINE_LENGTH,&lgEOL);
162  if( !lgEOL || (chk!=0.) )
163  {
164  /* got another number, not lgEOL, so too many numbers */
165  fprintf( ioQQQ, " There were more than %ld abundances entered\n",
166  abund.npSolar );
167  fprintf( ioQQQ, " Could there have been a typo somewhere?\n" );
168  }
169 
170  /* are numbers scale factors, or log of abund rel to H?? */
171  lgLog = false;
172  for( i=0; i < abund.npSolar; i++ )
173  {
174  if( absav[i] < 0. )
175  lgLog = true;
176  }
177 
178  if( lgLog )
179  {
180  /* entered as log of number rel to hydrogen */
181  for( i=0; i < abund.npSolar; i++ )
182  {
183  abund.solar[abund.ipSolar[i]-1] = (realnum)pow(10.,absav[i]);
184  }
185  }
186  else
187  {
188  /* scale factors relative to solar */
189  for( i=0; i < abund.npSolar; i++ )
190  {
191  abund.solar[abund.ipSolar[i]-1] *= (realnum)absav[i];
192  }
193  }
194 
195  /* check whether the abundances are reasonable */
196  if( abund.solar[1] > 0.2 )
197  {
198  fprintf( ioQQQ, " Is an abundance of%10.3e for %2.2s reasonable?\n",
200  }
201 
202  for( i=2; i < LIMELM; i++ )
203  {
204  if( abund.solar[i] > 0.1 )
205  {
206  fprintf( ioQQQ, " Is an abundance of%10.3e for %2.2s reasonable?\n",
208  }
209  }
210  return;
211  }
212 
213  /* this branch if no numbers on the line - option to use one of several
214  * special stored abundance sets */
215  if( nMatch(" AGB",chCard) || nMatch("AGB ",chCard) || nMatch("PLAN",chCard) )
216  {
217  /* AGB/planetary nebula abundances */
218  /* only turn on grains if "no grains" is not present */
219  if( !nMatch("O GR",chCard) )
220  {
221  /* turn on grains if not already done */
222  if( !lgDSet )
223  {
224  /* return bins allocated by previous abundances ... commands */
225  ReturnGrainBins();
226  /* now allocate new grain bins */
227  gp.dep = 1.;
228  gp.lgAbunVsDepth = false;
229  gp.lgRequestQHeating = true;
230  gp.lgForbidQHeating = false;
231  gp.lgGreyGrain = false;
232 
233  /* NO QHEAT option to turn off quantum heating for grains */
234  if( nMatch("O QH",chCard) )
235  {
236  gp.lgForbidQHeating = true;
237  gp.lgRequestQHeating = false;
238  phycon.lgPhysOK = false;
239  }
240 
241  /* actually set up the grains */
242  mie_read_opc("graphite_ism_10.opc",gp);
243  mie_read_opc("silicate_ism_10.opc",gp);
244  }
245  }
246 
247  for( i=0; i < LIMELM; i++ )
248  {
249  abund.solar[i] = abund.apn[i];
250  }
251  }
252 
253  else if( nMatch("HII ",chCard) || nMatch("H II",chCard) || nMatch("ORIO",chCard) )
254  {
255  /* H II region abundances */
256  /* only turn on grains if "no grains" is not present */
257  if( !nMatch("O GR",chCard) )
258  {
259  /* option to turn on grains */
260  if( !lgDSet )
261  {
262  /* return bins allocated by previous abundances ... commands */
263  ReturnGrainBins();
264  /* now allocate new grain bins */
265  gp.dep = 1.;
266  gp.lgAbunVsDepth = false;
267  gp.lgRequestQHeating = true;
268  gp.lgForbidQHeating = false;
269  gp.lgGreyGrain = false;
270 
271  /* NO QHEAT option to turn off quantum heating for grains */
272  if( nMatch("O QH",chCard) )
273  {
274  gp.lgForbidQHeating = true;
275  gp.lgRequestQHeating = false;
276  phycon.lgPhysOK = false;
277  }
278  /* This scales the Orion grain abundance so that the observed
279  * dust to gas ratio that Cloudy predicts is in agreement with
280  * that observed in the Veil,
281  *>>refer grain Abel, N., Brogan, C., Ferland, G., O'Dell, C.R.,
282  *>>refercon Shaw, G., Troland, T., 2004, ApJ, submitted */
283  gp.dep *= 0.85;
284 
285  mie_read_opc("graphite_orion_10.opc",gp);
286  mie_read_opc("silicate_orion_10.opc",gp);
287  }
288  }
289 
290  for( i=0; i < LIMELM; i++ )
291  {
292  abund.solar[i] = abund.ahii[i];
293  }
294  }
295 
296  else if( nMatch("ISM ",chCard) || nMatch(" ISM",chCard) )
297  {
298  /* ISM abundances from Cowie and Songaila Ann Rev '86 */
299  /* only turn on grains if "no grains" is not present */
300  if( !nMatch("O GR",chCard) )
301  {
302  if( !lgDSet )
303  {
304  /* return bins allocated by previous abundances ... commands */
305  ReturnGrainBins();
306  /* now allocate new grain bins */
307  gp.dep = 1.;
308  gp.lgAbunVsDepth = false;
309  gp.lgRequestQHeating = true;
310  gp.lgForbidQHeating = false;
311  gp.lgGreyGrain = false;
312 
313  /* NO QHEAT option to turn off quantum heating */
314  if( nMatch("O QH",chCard) )
315  {
316  gp.lgForbidQHeating = true;
317  gp.lgRequestQHeating = false;
318  phycon.lgPhysOK = false;
319  }
320 
321  /* actually set up the grains */
322  mie_read_opc("graphite_ism_10.opc",gp);
323  mie_read_opc("silicate_ism_10.opc",gp);
324  }
325  }
326 
327  for( i=0; i < LIMELM; i++ )
328  {
329  abund.solar[i] = abund.aism[i];
330  }
331  }
332 
333  else if( nMatch("NOVA",chCard) )
334  {
335  /* Nova Cyg abundances */
336  for( i=0; i < LIMELM; i++ )
337  {
338  abund.solar[i] = abund.anova[i];
339  }
340  }
341 
342  else if( nMatch("PRIM",chCard) )
343  {
344  /* roughly primordial abundances: He/H=.072 */
345  for( i=0; i < LIMELM; i++ )
346  {
347  abund.solar[i] = abund.aprim[i];
348  }
349 
350  /* now turn off the heavy elements */
351  for( i=4; i < LIMELM; i++ )
352  {
353  /* turn off heavy elements - do this way to make sure,
354  * that H-like and He-like level limits are handled properly */
355  char chDUMMY[LIMELM];
356  sprintf(chDUMMY,"element %s off ", elementnames.chElementName[i] );
357  /* now convert to caps */
358  caps(chDUMMY);
359  ParseElement( chDUMMY );
360  }
361  }
362 
363  else if( nMatch("CAME",chCard) )
364  {
365  /* mix from Cameron 1982, "Essays on Nuclear Astrophysics" */
366  for( i=0; i < LIMELM; i++ )
367  {
368  abund.solar[i] = abund.camern[i];
369  }
370  }
371 
372  else
373  {
374  fprintf( ioQQQ,
375  " ABUND must have PLAN, H II, CAMERON, NOVA, ALL, STARBURST, OLD SOLAR 84 or PRIMORDIAL. Sorry.\n" );
376  cdEXIT(EXIT_FAILURE);
377  }
378  return;
379 }

Generated for cloudy by doxygen 1.8.3.1