cloudy  trunk
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
opacity_addtotal.cpp
Go to the documentation of this file.
1 /* This file is part of Cloudy and is copyright (C)1978-2022 by Gary J. Ferland and
2  * others. For conditions of distribution and use see copyright notice in license.txt */
3 /*OpacityAddTotal derive total opacity for this position,
4  * called by ConvBase */
5 #include "cddefines.h"
6 #include "iso.h"
7 #include "ipoint.h"
8 #include "grainvar.h"
9 #include "ca.h"
10 #include "rfield.h"
11 #include "oxy.h"
12 #include "h2.h"
13 #include "hmi.h"
14 #include "atoms.h"
15 #include "conv.h"
16 #include "ionbal.h"
17 #include "trace.h"
18 #include "phycon.h"
19 #include "opacity.h"
20 #include "mole.h"
21 #include "freebound.h"
22 #include "dense.h"
23 #include "atmdat_gaunt.h"
24 #include "cool_eval.h"
25 
26 void OpacityAddTotal(void)
27 {
28  long int i,
29  ion,
30  ipHi,
31  ipISO,
32  ipop,
33  limit,
34  low,
35  n;
36  double DepartCoefInv ,
37  fac,
38  sum;
39  realnum SaveOxygen1 ,
40  SaveCarbon1;
41 
42  DEBUG_ENTRY( "OpacityAddTotal()" );
43 
44  /* OpacityZero will zero out scattering and absorption opacities,
45  * and set OldOpacSave to opac to save it */
46  OpacityZero();
47 
48  /* free electron scattering opacity, Compton recoil energy loss */
49  for( i=0; i < rfield.nflux; i++ )
50  {
51  /* scattering part of total opacity */
53  dense.eden;
54  }
55 
56  /* opacity due to Compton bound recoil ionization */
57  for( long nelem=ipHYDROGEN; nelem<LIMELM; ++nelem )
58  {
59  if( dense.lgElmtOn[nelem] )
60  {
61  for( ion=0; ion<nelem+1; ++ion )
62  {
63  realnum factor = dense.xIonDense[nelem][ion];
64  /*>>chng 05 nov 26, add molecular hydrogen assuming same
65  * as two free hydrogen atoms - hence mult density by two
66  *>>KEYWORD H2 bound Compton ionization */
67  if( nelem==ipHYDROGEN )
68  factor += hmi.H2_total*2.f;
69  if( factor > 0. )
70  {
71  // loop_min and loop_max are needed to work around a bug in icc 10.0
72  long loop_min = ionbal.ipCompRecoil[nelem][ion]-1;
73  long loop_max = rfield.nflux;
74  /* ionbal.nCompRecoilElec number of electrons in valence shell
75  * that can compton recoil ionize */
76  factor *= ionbal.nCompRecoilElec[nelem-ion];
77  for( i=loop_min; i < loop_max; i++ )
78  {
79  /* add in bound hydrogen electron scattering, treated as absorption */
80  opac.opacity_abs[i] += opac.OpacStack[i-1+opac.iopcom]*factor;
81  }
82  }
83  }
84  }
85  }
86 
87  /* opacity due to pair production - does not matter what form these
88  * elements are in */
92 
93  // only update ee brems spectrum when temperature changes
94  // the electron density dependence is explicitly factored in further below
96  {
99  }
100 
101  /* free free brems emission for all ions */
102 
103  /* Precomputed opacity value in OpacStack is inflated by 30 dex to avoid underflows;
104  * deflation by 1e-30 corrects for that. */
105  double bfac = 1e-30 * dense.eden / phycon.sqrte;
106  /* >>chng 02 jul 21, use full set of ions and gaunt factor */
107  /* gaunt factors depend only on photon energy and ion charge, so do
108  * sum of ions here before entering into loop over photon energy */
109  t_brems_den bsum;
111  // there is no need to keep these separate here...
112  bsum.den_ion[1] += bsum.den_Hp + bsum.den_Hep;
113  bsum.den_ion[2] += bsum.den_Hepp;
114 
115  for( i=0; i < rfield.nflux; i++ )
116  opac.FreeFreeOpacity[i] = 0.;
117 
119  for( ion=1; ion < LIMELM+1; ++ion )
120  if( bsum.den_ion[ion] > 0. )
121  t_gaunt::Inst().brems_opac( ion, phycon.te, bfac*bsum.den_ion[ion], opac.FreeFreeOpacity );
122 
123  for( i=0; i < rfield.nflux; i++ )
124  {
125  if( rfield.ContBoltz[i] < 0.995 )
127  else
130  }
131 
132  /* H minus absorption, with correction for stimulated emission */
133  if( hmi.hmidep > SMALLFLOAT )
134  {
135  DepartCoefInv = 1./hmi.hmidep;
136  }
137  else
138  {
139  /* the hmidep departure coef can become vastly small in totally
140  * neutral gas (no electrons) */
141  DepartCoefInv = 1.;
142  }
143 
144  limit = iso_sp[ipHE_LIKE][ipHELIUM].fb[0].ipIsoLevNIonCon;
145  if(limit > rfield.nflux)
146  limit = rfield.nflux;
147 
148  const molezone *MHm = findspecieslocal("H-");
149  for( i=hmi.iphmin-1; i < limit; i++ )
150  {
151  double factor;
152  factor = 1. - rfield.ContBoltz[i]*DepartCoefInv;
153  if(factor > 0)
155  MHm->den*factor;
156  }
157 
158  /* H2P h2plus bound free opacity */
159  {
160  const molezone *MH2p = findspecieslocal("H2+");
161  limit = opac.ih2pnt[1];
162  if(limit > rfield.nflux)
163  limit = rfield.nflux;
164  double frac = MAX2( 0., 1.-hmi.h2plus_exc_frac );
165  for( i=opac.ih2pnt[0]-1; i < limit; i++ )
166  {
167  opac.opacity_abs[i] += MH2p->den*frac*opac.OpacStack[i-opac.ih2pnt[0]+
168  opac.ih2pof];
169  ASSERT( !isnan( opac.opacity_abs[i] ) );
170  }
171 
172  frac = hmi.h2plus_exc_frac;
173  limit = opac.ih2pnt_ex[1];
174  if(limit > rfield.nflux)
175  limit = rfield.nflux;
176  for( i=opac.ih2pnt_ex[0]-1; i < limit; i++ )
177  {
178  opac.opacity_abs[i] += MH2p->den*frac*opac.OpacStack[i-opac.ih2pnt_ex[0]+
179  opac.ih2pof_ex];
180  ASSERT( !isnan( opac.opacity_abs[i] ) );
181  }
182  }
183 
184  /* H2 continuum dissociation opacity */
185  for( diatom_iter diatom = diatoms.begin(); diatom != diatoms.end(); ++diatom )
186  {
187  if( (*diatom)->lgEnabled && mole_global.lgStancil )
188  {
189  for( vector< diss_tran >::iterator tran = (*diatom)->Diss_Trans.begin(); tran != (*diatom)->Diss_Trans.end(); ++tran )
190  {
191  long lower_limit = ipoint(tran->energies[0]);
192  long upper_limit = ipoint(tran->energies.back());
193  upper_limit = MIN2( upper_limit, rfield.nflux-1 );
194  for(i = lower_limit; i <= upper_limit; ++i)
195  {
196  opac.opacity_abs[i] += (*diatom)->MolDissocOpacity( *tran, rfield.anu(i));
197  }
198  }
199  }
200  }
201 
202  /* get total population of hydrogen ground to do Rayleigh scattering */
203  if( dense.xIonDense[ipHYDROGEN][1] <= 0. )
204  {
205  fac = dense.xIonDense[ipHYDROGEN][0];
206  }
207  else
208  {
209  fac = iso_sp[ipH_LIKE][ipHYDROGEN].st[ipH1s].Pop();
210  }
211 
212  /* Ly a damp wing opac (Rayleigh scattering) */
213  limit = iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH1s].ipIsoLevNIonCon;
214  if(limit > rfield.nflux)
215  limit = rfield.nflux;
216  for( i=0; i < limit; i++ )
217  {
218  opac.opacity_sct[i] += (fac*opac.OpacStack[i-1+opac.ipRayScat]);
219  }
220 
221  /* remember largest correction for stimulated emission */
222  if( iso_sp[ipH_LIKE][ipHYDROGEN].st[ipH1s].DepartCoef() > 1e-30 && !conv.lgSearch )
223  {
224  realnum factor;
225  factor = (realnum)(rfield.ContBoltz[iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH1s].ipIsoLevNIonCon-1]/iso_sp[ipH_LIKE][ipHYDROGEN].st[ipH1s].DepartCoef());
226  if(opac.stimax[0] < factor)
227  opac.stimax[0] = factor;
228  }
229 
230  if( iso_sp[ipH_LIKE][ipHYDROGEN].st[ipH2p].DepartCoef() > 1e-30 && !conv.lgSearch )
231  {
232  realnum factor;
233  factor = (realnum)(rfield.ContBoltz[iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH2p].ipIsoLevNIonCon-1]/iso_sp[ipH_LIKE][ipHYDROGEN].st[ipH2p].DepartCoef());
234  if(opac.stimax[1] < factor)
235  opac.stimax[1] = factor;
236  }
237 
238 # if 0
239  /* check whether hydrogen or Helium singlets mased, if not in search mode */
240  if( !conv.lgSearch )
241  {
242  if( iso_sp[ipH_LIKE][ipHYDROGEN].trans(ipH2p,ipH1s).PopOpc() < 0. )
243  {
244  hydro.lgHLyaMased = true;
245  }
246  }
247 # endif
248 
249  /* >>chng 05 nov 25, use Yan et al. H2 photo cs
250  * following reference gives cross section for all energies
251  * >>refer H2 photo cs Yan, M., Sadeghpour, H.R., & Dalgarno, A., 1998, ApJ, 496, 1044
252  * Wilms, J., Allen, A., & McCray, R. 2000, ApJ, 542, 914 */
253  /* >>chng 02 jan 16, approximate inclusion of H_2 photoelectric opacity
254  * include H_2 in total photoelectric opacity as twice H0 cs */
255  /* set lower and upper limits to this range */
256  /*>>KEYWORD H2 photoionization opacity */
257  low = h2.ip_photo_opac_thresh;
258  ipHi = rfield.nflux_with_check;
259  ipop = h2.ip_photo_opac_offset;
260  /* OpacityAdd1Subshell just returns for static opacities if opac.lgRedoStatic not set*/
261  /* >>chng 05 nov 27, change on nov 25 had left 2*density from H0, so
262  * twice the H2 density was used - * also changed to static opacity
263  * this assumes that all v,J levels contribute the same opacity */
264  OpacityAdd1Subshell( ipop , low , ipHi , hmi.H2_total , 's' );
265 
266  /*>>KEYWORD CO photoionization opacity */
267  /* include photoionization of CO - assume C and O in CO each have
268  * same photo cs as atom - this should only be significant in highly
269  * shielded regions where only very hard photons penetrate
270  * also H2O condensed onto grain surfaces - very important deep in cloud */
271  SaveOxygen1 = dense.xIonDense[ipOXYGEN][0];
272  SaveCarbon1 = dense.xIonDense[ipCARBON][0];
273  /* atomic C and O will include CO during the heating sum loop */
274  fixit("CO opacity hack breaks the invariant for total mass");
275  /*
276  * In any case, is CO the only species which contributes?
277  * -- might expect all other molecular species to do so,
278  * i.e. the item added should perhaps be
279  * dense.xMolecules[nelem]) -- code duplicated in heatsum
280  */
281  dense.xIonDense[ipOXYGEN][0] += (realnum) (findspecieslocal("CO")->den + findspecieslocal("H2Ogrn")->den);
283 
284  /* following loop adds standard opacities for first 30 elements
285  * most heavy element opacity added here */
286  for( long nelem=ipHYDROGEN; nelem < LIMELM; nelem++ )
287  {
288  /* this element may be turned off */
289  if( dense.lgElmtOn[nelem] )
290  {
291  OpacityAdd1Element(nelem);
292  }
293  }
294 
295  /* now reset the abundances */
296  dense.xIonDense[ipOXYGEN][0] = SaveOxygen1;
297  dense.xIonDense[ipCARBON][0] = SaveCarbon1;
298 
299  /* following are opacities due to specific excited levels */
300 
301  /* nitrogen opacity
302  * excited level of N+ */
304  dense.xIonDense[ipNITROGEN][0]*atoms.p2nit , 'v' );
305 
306  /* oxygen opacity
307  * excited level of Oo */
309  dense.xIonDense[ipOXYGEN][0]*oxy.poiexc,'v');
310 
311  /* O2+ excited states */
313  dense.xIonDense[ipOXYGEN][2]*oxy.poiii2,'v');
314 
316  dense.xIonDense[ipOXYGEN][2]*oxy.poiii3,'v');
317 
318  /* magnesium opacity
319  * excited level of Mg+ */
322 
323  /* calcium opacity
324  * photoionization of excited levels of Ca+ */
326  ca.popca2ex,'v');
327 
328  /*******************************************************************
329  *
330  * complete evaluation of total opacity by adding in the static part and grains
331  *
332  *******************************************************************/
333 
334  /* this loop defines the variable iso_sp[ipH_LIKE][nelem].fb[n].ConOpacRatio,
335  * the ratio of not H to Hydrogen opacity. for grain free environments
336  * at low densities this is nearly zero. The correction includes
337  * stimulated emission correction */
338  for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
339  {
340  for( long nelem=ipISO; nelem < LIMELM; nelem++ )
341  {
342  /* this element may be turned off */
343  if( dense.lgElmtOn[nelem] )
344  {
345  /* this branch is for startup only */
346  if( nzone < 1 )
347  {
348  /* >>chng 06 aug 17, should go to numLevels_local instead of _max */
349  for( n=0; n < iso_sp[ipISO][nelem].numLevels_local; n++ )
350  {
351  if(iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon < rfield.nflux )
352  {
353  /*>>chng 04 dec 12, had tested against 1e-30, set to zero if not
354  * greater than this - caused oscillations as opacity fell below
355  * and around this value - change to greater than 0 */
356  /*if( opac.opacity_abs[iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon-1] > 1e-30 )*/
357  if( opac.opacity_abs[iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon-1] > 0. )
358  {
359  /* >>chng 02 may 06, use general form of threshold cs */
360  /*double t1 = atmdat_sth(n)/(POW2(nelem+1.-ipISO));*/
361  long int ip = iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon;
362  double t2 = csphot(
363  ip ,
364  ip ,
365  iso_sp[ipISO][nelem].fb[n].ipOpac );
366 
367  iso_sp[ipISO][nelem].fb[n].ConOpacRatio = 1.-
368  (iso_sp[ipISO][nelem].st[n].Pop()*t2)/
369  opac.opacity_abs[ip-1];
370  }
371  else
372  {
373  iso_sp[ipISO][nelem].fb[n].ConOpacRatio = 0.;
374  }
375  }
376  }
377  }
378  /* end branch for startup only, start branch for all zones including startup */
379  /* >>chng 06 aug 17, should go to numLevels_local instead of _max */
380  for( n=0; n < iso_sp[ipISO][nelem].numLevels_local; n++ )
381  {
382  /* ratios of other to total opacity for continua of all atoms done with iso model */
383  if(iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon < rfield.nflux )
384  {
385  /*>>chng 04 dec 12, had tested against 1e-30, set to zero if not
386  * greater than this - caused oscillations as opacity fell below
387  * and around this value - change to greater than 0 */
388  /*if( opac.opacity_abs[iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon-1] > 1e-30 )*/
389  if( opac.opacity_abs[iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon-1] > 0. )
390  {
391  /* first get departure coef */
392  if( iso_sp[ipISO][nelem].st[n].DepartCoef() > 1e-30 && (!conv.lgSearch ) )
393  {
394  /* this is the usual case, use inverse of departure coef */
395  fac = 1./iso_sp[ipISO][nelem].st[n].DepartCoef();
396  }
397  else if( conv.lgSearch )
398  {
399  /* do not make correction for stim emission during search
400  * for initial temperature solution, since trys are very non-equil */
401  fac = 0.;
402  }
403  else
404  {
405  fac = 1.;
406  }
407 
410  /* now get opaicty ratio with correction for stimulated emission */
411  /*>>chng 04 dec 12, had tested against 1e-30, set to zero if not
412  * greater than this - caused oscillations as opacity fell below
413  * and around this value - change to greater than 0 */
414  /*if( opac.opacity_abs[iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon-1] > 1e-30 )*/
415  if( opac.opacity_abs[iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon-1] > 0. )
416  {
417  /* >>chng 02 may 06, use general form of threshold cs */
418  long int ip = iso_sp[ipISO][nelem].fb[n].ipIsoLevNIonCon;
419 
420  double t2 = csphot(
421  ip ,
422  ip ,
423  iso_sp[ipISO][nelem].fb[n].ipOpac );
424 
425  double opacity_this_species =
426  iso_sp[ipISO][nelem].st[n].Pop()*t2*
427  (1. - fac*rfield.ContBoltz[ip-1]);
428 
429  double opacity_fraction = 1. - opacity_this_species / opac.opacity_abs[ip-1];
430  if(opacity_fraction < 0)
431  opacity_fraction = 0.;
432 
433  /* use mean of old and new ratios */
434  iso_sp[ipISO][nelem].fb[n].ConOpacRatio =
435  iso_sp[ipISO][nelem].fb[n].ConOpacRatio* 0.75 + 0.25*opacity_fraction;
436 
437  if(iso_sp[ipISO][nelem].fb[n].ConOpacRatio < 0.)
438  iso_sp[ipISO][nelem].fb[n].ConOpacRatio = 0.;
439  }
440  else
441  {
442  iso_sp[ipISO][nelem].fb[n].ConOpacRatio = 0.;
443  }
444  }
445  else
446  {
447  iso_sp[ipISO][nelem].fb[n].ConOpacRatio = 0.;
448  }
449  }
450  else
451  {
452  iso_sp[ipISO][nelem].fb[n].ConOpacRatio = 0.;
453  }
454  }
455  }
456  }
457  }
458 
459  /* add dust grain opacity if dust present */
460  if( gv.lgDustOn() )
461  {
462  /* generate current grain opacities since may be function of depth */
463  /* >>chng 01 may 11, removed code to update grain opacities, already done by GrainChargeTemp */
464  for( i=0; i < rfield.nflux; i++ )
465  {
466  /* units cm-1 */
469  }
470  }
471 
472  /* check that opacity is sane */
473  for( i=0; i < rfield.nflux; i++ )
474  {
475  /* OpacStatic was zeroed in OpacityZero, incremented in opacityadd1subshell */
476  opac.opacity_abs[i] += opac.OpacStatic[i];
477  /* make sure that opacity is positive */
478  /*ASSERT( opac.opacity_abs[i] > 0. );*/
479  }
480 
481  /* compute gas albedo here */
482  for( i=0; i < rfield.nflux; i++ )
483  {
484  opac.albedo[i] = opac.opacity_sct[i]/
485  (opac.opacity_sct[i] + opac.opacity_abs[i]);
486  }
487 
488  /* during search phase set old opacity array to current value */
489  if( conv.lgSearch )
490  OpacityZeroOld();
491 
492  if( trace.lgTrace )
493  {
494  fprintf( ioQQQ, " OpacityAddTotal returns; grd rec eff (opac) for Hn=1,4%10.2e%10.2e%10.2e%10.2e%10.2e%10.2e%10.2e%10.2e%10.2e%10.2e HeI,II:%10.2e%10.2e\n",
495  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH1s].ConOpacRatio,
496  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH2s].ConOpacRatio,
497  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH2p].ConOpacRatio,
498  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH3s].ConOpacRatio,
499  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH3p].ConOpacRatio,
500  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH3d].ConOpacRatio,
501  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH4s].ConOpacRatio,
502  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH4p].ConOpacRatio,
503  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH4d].ConOpacRatio,
504  iso_sp[ipH_LIKE][ipHYDROGEN].fb[ipH4f].ConOpacRatio,
505  iso_sp[ipHE_LIKE][ipHELIUM].fb[ipHe1s1S].ConOpacRatio,
506  iso_sp[ipH_LIKE][ipHELIUM].fb[ipH1s].ConOpacRatio );
507  }
508 
509  {
510  /* following should be set true to print strongest ots contributors */
511  enum {DEBUG_LOC=false};
512  if( DEBUG_LOC && (nzone>=378) )
513  {
514  if( nzone > 380 )
515  cdEXIT( EXIT_FAILURE );
516  for( i=0; i<rfield.nflux; ++i )
517  {
518  fprintf(ioQQQ,"rtotsbugggg\t%li\t%.3e\t%.3e\t%.3e\t%.3e\n",
520  rfield.anu(i),
521  opac.opacity_abs[i],
522  rfield.otscon[i],
523  rfield.otslin[i]);
524  }
525  }
526  }
527  return;
528 }
void OpacityAdd1Element(long int ipZ)
realnum popca2ex
Definition: ca.h:43
long int iphmin
Definition: hmi.h:128
t_mole_global mole_global
Definition: mole.cpp:7
bool lgStancil
Definition: mole.h:337
vector< double > dstab
Definition: grainvar.h:525
vector< double, allocator_avx< double > > ContBoltz
Definition: rfield.h:124
long int ipRayScat
Definition: opacity.h:223
long int nCompRecoilElec[LIMELM]
Definition: ionbal.h:184
const int ipH4d
Definition: iso.h:37
const int ipMAGNESIUM
Definition: cddefines.h:360
double * OpacStack
Definition: opacity.h:164
double * opacity_abs
Definition: opacity.h:104
qList st
Definition: iso.h:482
double * albedo
Definition: opacity.h:113
const int ipHE_LIKE
Definition: iso.h:65
t_opac opac
Definition: opacity.cpp:5
long ip_photo_opac_thresh
Definition: h2_priv.h:320
double * OpacStatic
Definition: opacity.h:123
const realnum SMALLFLOAT
Definition: cpu.h:246
const int NISO
Definition: cddefines.h:311
long int ipo1exc[3]
Definition: opacity.h:223
long ipmgex
Definition: opacity.h:300
double den_Hepp
Definition: atmdat_gaunt.h:40
const int ipOXYGEN
Definition: cddefines.h:356
double den_Hep
Definition: atmdat_gaunt.h:39
t_conv conv
Definition: conv.cpp:5
t_phycon phycon
Definition: phycon.cpp:6
long int ih2pof
Definition: opacity.h:223
realnum stimax[2]
Definition: opacity.h:317
long int ipo3exc[3]
Definition: opacity.h:223
FILE * ioQQQ
Definition: cddefines.cpp:7
double * opacity_sct
Definition: opacity.h:107
molezone * findspecieslocal(const char buf[])
long int nzone
Definition: cddefines.cpp:14
const int ipHe1s1S
Definition: iso.h:43
Definition: mole.h:378
void OpacityAdd1Subshell(long int ipOpac, long int ipLowLim, long int ipUpLim, realnum abundance, char chStat)
long int ipBrems
Definition: opacity.h:223
vector< freeBound > fb
Definition: iso.h:481
#define MIN2(a, b)
Definition: cddefines.h:803
double * eeFreeFreeOpacity
Definition: opacity.h:127
double anu(size_t i) const
Definition: mesh.h:120
long int ih2pof_ex
Definition: opacity.h:223
t_ca ca
Definition: ca.cpp:5
t_dense dense
Definition: global.cpp:15
static t_gaunt & Inst()
Definition: cddefines.h:209
realnum popMg2
Definition: atoms.h:129
t_iso_sp iso_sp[NISO][LIMELM]
Definition: iso.cpp:11
long int nflux_with_check
Definition: rfield.h:49
vector< double > dstsc
Definition: grainvar.h:526
double xIonDense[LIMELM][LIMELM+1]
Definition: dense.h:135
realnum * otslin
Definition: rfield.h:183
bool lgSearch
Definition: conv.h:168
t_trace trace
Definition: trace.cpp:5
bool fp_equal(sys_float x, sys_float y, int n=3)
Definition: cddefines.h:854
t_ionbal ionbal
Definition: ionbal.cpp:8
long ipoint(double energy_ryd)
Definition: cont_ipoint.cpp:15
realnum poiii3
Definition: oxy.h:19
const int ipH1s
Definition: iso.h:29
long int nPres2Ioniz
Definition: conv.h:145
bool lgTrace
Definition: trace.h:12
void OpacityZeroOld(void)
t_rfield rfield
Definition: rfield.cpp:9
void OpacityAddTotal(void)
t_atoms atoms
Definition: atoms.cpp:5
float realnum
Definition: cddefines.h:124
#define EXIT_FAILURE
Definition: cddefines.h:168
const int ipH4p
Definition: iso.h:36
long int iopcom
Definition: opacity.h:223
long ipOpMgEx
Definition: opacity.h:300
vector< diatomics * > diatoms
Definition: h2.cpp:8
t_hydro hydro
Definition: hydrogenic.cpp:5
realnum * otscon
Definition: rfield.h:183
#define cdEXIT(FAIL)
Definition: cddefines.h:482
long ip_photo_opac_offset
Definition: h2_priv.h:321
void eeBremsSpectrum(double Te, realnum jnu[], double knu[])
Definition: cool_eval.cpp:1798
double eeFreeFreeTemp
Definition: opacity.h:130
diatomics h2("h2", 4100.,&hmi.H2_total, Yan_H2_CS)
const int ipH3s
Definition: iso.h:32
realnum * eeBremsDif
Definition: rfield.h:273
double den_ion[LIMELM+1]
Definition: atmdat_gaunt.h:42
long int iphmop
Definition: opacity.h:223
const int ipH3d
Definition: iso.h:34
long int ioppr
Definition: opacity.h:223
bool lgElmtOn[LIMELM]
Definition: dense.h:160
realnum gas_phase[LIMELM]
Definition: dense.h:76
const int ipH2p
Definition: iso.h:31
long int ipo3exc3[3]
Definition: opacity.h:223
#define ASSERT(exp)
Definition: cddefines.h:613
long int ippr
Definition: opacity.h:223
const int ipH2s
Definition: iso.h:30
realnum p2nit
Definition: atoms.h:126
double den_Hm
Definition: atmdat_gaunt.h:37
const int ipNITROGEN
Definition: cddefines.h:355
const int ipH_LIKE
Definition: iso.h:64
double csphot(long int inu, long int ithr, long int iofset)
Definition: service.cpp:1655
const int LIMELM
Definition: cddefines.h:308
realnum poiii2
Definition: oxy.h:19
T pow2(T a)
Definition: cddefines.h:981
double hmidep
Definition: hmi.h:42
double den
Definition: mole.h:421
realnum poiexc
Definition: oxy.h:19
#define DEBUG_ENTRY(funcname)
Definition: cddefines.h:723
#define isnan
Definition: cddefines.h:659
const int ipHELIUM
Definition: cddefines.h:350
double * FreeFreeOpacity
Definition: opacity.h:126
long int in1[3]
Definition: opacity.h:223
double H2_total
Definition: hmi.h:25
double eden
Definition: dense.h:201
t_oxy oxy
Definition: oxy.cpp:5
#define MAX2(a, b)
Definition: cddefines.h:824
int fprintf(const Output &stream, const char *format,...)
Definition: service.cpp:1121
const int ipH4s
Definition: iso.h:35
const int ipCARBON
Definition: cddefines.h:354
long ica2ex[2]
Definition: opacity.h:300
double sqrte
Definition: phycon.h:58
#define fixit(a)
Definition: cddefines.h:417
long int ih2pnt[2]
Definition: opacity.h:223
GrainVar gv
Definition: grainvar.cpp:5
t_hmi hmi
Definition: hmi.cpp:5
void brems_sum_ions(t_brems_den &sum) const
double te
Definition: phycon.h:21
void brems_opac(long ion, double Te, double abun, double arr[])
double frac(double d)
bool lgDustOn() const
Definition: grainvar.h:475
const int ipHYDROGEN
Definition: cddefines.h:349
long int nflux
Definition: rfield.h:46
long int ih2pnt_ex[2]
Definition: opacity.h:223
const int ipH3p
Definition: iso.h:33
long int numLevels_local
Definition: iso.h:529
const int ipH4f
Definition: iso.h:38
vector< diatomics * >::iterator diatom_iter
Definition: h2.h:13
double den_Hp
Definition: atmdat_gaunt.h:38
long ica2op
Definition: opacity.h:300
double h2plus_exc_frac
Definition: hmi.h:50
void OpacityZero(void)
Definition: opacity_zero.cpp:8
long int ** ipCompRecoil
Definition: ionbal.h:159