cloudy  trunk
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ion_iron.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 /*IonIron ionization balance for iron */
4 #include "cddefines.h"
5 #include "dense.h"
6 #include "fe.h"
7 #include "gammas.h"
8 #include "opacity.h"
9 #include "trace.h"
10 #include "grainvar.h"
11 #include "ionbal.h"
12 #include "mole.h"
13 
14 /*IonIron ionization balance for iron */
15 void IonIron(void)
16 {
17  const int NDIM = ipIRON+1;
18 
19  static const double dicoef[2][NDIM] = {
20  {1.58e-3,8.38e-3,1.54e-2,3.75e-2,0.117,0.254,0.291,0.150,0.140,0.100,0.200,0.240,
21  0.260,0.190,0.120,0.350,0.066,0.10,0.13,0.23,0.14,0.11,0.041,0.747,0.519,0.},
22  {.456,.323,.310,.411,.359,.0975,.229,4.20,3.30,5.30,1.50,0.700,.600,.5,1.,0.,7.8,
23  6.3,5.5,3.6,4.9,1.6,4.2,.284,.279,0.}
24  };
25  static const double dite[2][NDIM] = {
26  {6.00e4,1.94e5,3.31e5,4.32e5,6.28e5,7.50e5,7.73e5,2.62e5,2.50e5,2.57e5,2.84e5,
27  8.69e5,4.21e5,4.57e5,2.85e5,8.18e5,1.51e6,1.30e6,1.19e6,1.09e6,9.62e5,7.23e5,
28  4.23e5,5.84e7,6.01e7,0.},
29  {8.97e4,1.71e5,2.73e5,3.49e5,5.29e5,4.69e5,6.54e5,1.32e6,1.33e6,1.41e6,1.52e6,
30  1.51e6,1.82e6,1.84e6,2.31e6,0.,9.98e6,9.98e6,1.00e7,1.10e7,8.34e6,1.01e7,1.07e7,
31  1.17e7,9.97e6,0.}
32  };
33  static const double ditcrt[NDIM] = {1e11,1e11,1e11,1e11,1e11,1e11,1e11,
34  1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,
35  1e11,1e11,1e11,1e11,1e11,1e11,1e11};
36  static const double aa[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
37  0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
38  static const double bb[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
39  0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
40  static const double cc[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
41  0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
42  static const double dd[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
43  0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
44  /* >>chng 03 aug 15, zero out this array. having non-zero values in
45  * this array had the effect of using 0 for the guess to the low-T dr
46  * in makerecomb */
47  static const double ff[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
48  0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
49  static const double fyield[NDIM+1] = {.34,.34,.35,.35,.36,.37,.37,.38,.39,.40,
50  .41,.42,.43,.44,.45,.46,.47,.47,.48,.48,.49,.49,.11,.75,0.,0.,0.};
51 
52  long int i, limit, limit2;
53 
54  DEBUG_ENTRY( "IonIron()" );
55 
56  /* iron nelem=26
57  *
58  * Fe rates from Woods et al. Ap.J. 249, 399.
59  * rec from +23, 24 25 from Arnauld et al 85 */
60 
61  /* Pequignot and Aldrovandi Ast Ap 161, 169. */
62 
63  if( !dense.lgElmtOn[ipIRON] )
64  {
65  fe.fekcld = 0.;
66  fe.fekhot = 0.;
67  fe.fegrain = 0.;
68  return;
69  }
70 
72 
73  ion_photo(ipIRON,false);
74  /* debugging printout for shell photo rates */
75  /*GammaPrtRate( ioQQQ , 0 , ipIRON , true);*/
76  /*GammaPrtShells( ipIRON , 13 );
77  GammaPrtShells( ipIRON , 12 );
78  GammaPrtShells( ipIRON , 0 );*/
79  {
80  /*@-redef@*/
81  enum {DEBUG_LOC=false};
82  /*@+redef@*/
83  if( DEBUG_LOC )
84  {
85  long int iplow , iphi , ipop , ns , ion;
86  double rate;
87  ns = 6;
88  ion = 0;
89  /* show what some of the ionization sources are */
90  iplow = opac.ipElement[ipIRON][ion][ns][0];
91  iphi = opac.ipElement[ipIRON][ion][ns][1];
92  ipop = opac.ipElement[ipIRON][ion][ns][2];
93  rate = ionbal.PhotoRate_Shell[ipIRON][ion][ns][0];
94  GammaPrt( iplow , iphi , ipop , ioQQQ, rate , rate*0.1 );
95  }
96  }
97 
98  /* find collisional ionization rates */
100 
101  /* get recombination coefficients */
102  ion_recomb(false,(const double*)dicoef,(const double*)dite,ditcrt,aa,bb,cc,dd,ff,ipIRON);
103 
104  /* >>chng 06 Feb 28 -- NPA. Add in charge transfer ionization of Mg with S+, Si+, and C+ */
105  /* only include this if molecular network is enabled - otherwise no feedback onto
106  * Si+, S+, and C+ soln */
107  if( !co.lgNoCOMole )
108  {
109  ionbal.PhotoRate_Shell[ipIRON][0][6][0] +=
110  CO_findrk("S+,Fe=>S,Fe+")*dense.xIonDense[ipSULPHUR][1] +
111  CO_findrk("Si+,Fe=>Si,Fe+")*dense.xIonDense[ipSILICON][1] +
112  CO_findrk("C+,Fe=>C,Fe+")*dense.xIonDense[ipCARBON][1];
113  /* CO_sink_rate("Fe");*/
114 
115  }
116 
117  /* solve for ionization balance */
118  ion_solver(ipIRON,false);
119 
120  /* now find total Auger yield of K-alphas
121  * "cold" iron has M-shell electrons, up to Fe 18 */
122  fe.fekcld = 0.;
123  limit = MIN2(18,dense.IonHigh[ipIRON]);
124 
125  for( i=dense.IonLow[ipIRON]; i < limit; i++ )
126  {
127  ASSERT( i < NDIM + 1 );
128  fe.fekcld +=
130  fyield[i]);
131  }
132 
133  /* same sum for hot iron */
134  fe.fekhot = 0.;
135  limit = MAX2(18,dense.IonLow[ipIRON]);
136 
137  limit2 = MIN2(ipIRON+1,dense.IonHigh[ipIRON]);
138  ASSERT( limit2 <= LIMELM + 1 );
139 
140  for( i=limit; i < limit2; i++ )
141  {
142  ASSERT( i < NDIM + 1 );
143  fe.fekhot +=
145  fyield[i]);
146  }
147 
148  /* Fe Ka from grains - Fe in grains assumed to be atomic
149  * gv.elmSumAbund[ipIRON] is number density of iron added over all grain species */
150  i = 0;
151  /* fyield is 0.34 for atomic fe */
152  fe.fegrain = ( gv.lgWD01 ) ? 0.f : (realnum)(ionbal.PhotoRate_Shell[ipIRON][i][0][0]*fyield[i]*
154 
155  if( trace.lgTrace && trace.lgHeavyBug )
156  {
157  /* print densities of each stage of ionization */
158  fprintf( ioQQQ, " Fe" );
159  for( i=0; i < 15; i++ )
160  {
161  fprintf( ioQQQ, "\t%.1e", dense.xIonDense[ipIRON][i]/dense.gas_phase[ipIRON] );
162  }
163  fprintf( ioQQQ, "\n" );
164  }
165 
166  if( trace.lgTrace && trace.lgFeBug )
167  {
168  fprintf( ioQQQ, " IonIron-Abund:" );
169  for( i=0; i < 27; i++ )
170  {
171  fprintf( ioQQQ, "%10.2e", dense.xIonDense[ipIRON][i] );
172  }
173  fprintf( ioQQQ, "\n" );
174 
175  fprintf( ioQQQ, " IonIron - Ka(Auger)%10.2e\n", fe.fekcld +
176  fe.fekhot );
177  }
178  return;
179 }

Generated for cloudy by doxygen 1.8.1.2