cloudy  trunk
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ion_photo.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 /*ion_photo fill array PhotoRate with photoionization rates for heavy elements */
4 #include "cddefines.h"
5 #include "yield.h"
6 #include "heavy.h"
7 #include "opacity.h"
8 #include "dense.h"
9 #include "thermal.h"
10 #include "conv.h"
11 #include "grainvar.h"
12 #include "elementnames.h"
13 #include "gammas.h"
14 #include "ionbal.h"
15 
16 void ion_photo(
17  /* nlem is atomic number on C scale, 0 for H */
18  long int nelem ,
19  /* debugging flag to turn on print */
20  bool lgPrintIt )
21 {
22  long int ion,
23  iphi,
24  iplow,
25  ipop,
26  limit_hi,
27  limit_lo,
28  ns;
29 
30  DEBUG_ENTRY( "ion_photo()" );
31 
32  /* IonLow(nelem) and IonHigh(nelem) are bounds for evaluation*/
33 
34  /*begin sanity checks */
35  ASSERT( nelem < LIMELM );
36  ASSERT( dense.IonLow[nelem] >= 0 );
37  ASSERT( dense.IonLow[nelem] <= nelem);
38  ASSERT( dense.IonHigh[nelem] <= nelem + 1);
39  /*end sanity checks */
40 
41  /* NB - in following, nelem is on c scale, so is 29 for Zn */
42 
43  /* min since iso-like atom rates produced in iso_photo.
44  * IonHigh and IonLow range from 0 to nelem+1, but for photo rates
45  * we want 0 to nelem since cannot destroy fully stripped ion.
46  * since iso-seq done elsewhere, want to actually do IonHigh-xxx.*/
47  /* >>chng 00 dec 07, logic on limit_hi now precisely identical to ion_solver */
48  /* >>chng 02 mar 31, change limit_hi to < in loop (had been <=) and
49  * also coded for arbitrary number of iso sequences */
50  /*limit_hi = MIN2(nelem,dense.IonHigh[nelem]);
51  limit_hi = MIN2(nelem-2,dense.IonHigh[nelem]-1);*/
52  limit_hi = MIN2( dense.IonHigh[nelem] , nelem+1-NISO );
53 
54  /* >>chng 03 sep 26, always do atom itself since may be needed for molecules */
55  limit_hi = MAX2( 1 , limit_hi );
56 
57  /* when grains are present want to use atoms as lower bound to number of stages of ionization,
58  * since atomic rates needed for species within grains */
59  if( !conv.nPres2Ioniz && gv.lgDustOn )
60  {
61  limit_lo = 0;
62  }
63  else
64  {
65  limit_lo = dense.IonLow[nelem];
66  }
67 
68  /* >>chng 01 dec 11, lower bound now limit_lo */
69  /* loop over all ions for this element */
70  /* >>chng 02 mar 31, now ion < limit_hi not <= */
71  for( ion=limit_lo; ion < limit_hi; ion++ )
72  {
73  /* loop over all shells for this ion */
74  for( ns=0; ns < Heavy.nsShells[nelem][ion]; ns++ )
75  {
76  /* always reevaluate the outer shell, and all shells if lgRedoStatic is set */
77  if( (ns==(Heavy.nsShells[nelem][ion]-1) || opac.lgRedoStatic) )
78  {
79  /* option to redo the rates only on occasion */
80  iplow = opac.ipElement[nelem][ion][ns][0];
81  iphi = opac.ipElement[nelem][ion][ns][1];
82  ipop = opac.ipElement[nelem][ion][ns][2];
83 
84  /* compute the photoionization rate, ionbal.lgPhotoIoniz_On is 1, set 0
85  * with "no photoionization" command */
86  ionbal.PhotoRate_Shell[nelem][ion][ns][0] =
87  GammaK(iplow,iphi,
88  ipop,t_yield::Inst().elec_eject_frac(nelem,ion,ns,0))*ionbal.lgPhotoIoniz_On;
89 
90  /* these three lines must be kept parallel with the lines
91  * in GammaK ion*/
92 
93  /* the heating rate */
96  }
97  }
98 
99  /* add on compton recoil ionization for atoms to outer shell */
100  /* >>chng 02 mar 24, moved here from ion_solver */
101  /* this is the outer shell */
102  ns = (Heavy.nsShells[nelem][ion]-1);
103  /* this must be moved to photoionize and have code parallel to iso_photo code */
104  ionbal.PhotoRate_Shell[nelem][ion][ns][0] += ionbal.CompRecoilIonRate[nelem][ion];
105  /* add the heat as secondary-ionization capable heating */
106  ionbal.PhotoRate_Shell[nelem][ion][ns][2] += ionbal.CompRecoilHeatRate[nelem][ion];
107  }
108 
109  /* option to print information about these rates for this element */
110  if( lgPrintIt )
111  {
112  /* option to print rates for particular shell */
113  ns = 5;
114  ion = 1;
115  GammaPrt(
116  opac.ipElement[nelem][ion][ns][0],
117  opac.ipElement[nelem][ion][ns][1],
118  opac.ipElement[nelem][ion][ns][2],
119  ioQQQ, /* io unit we will write to */
120  ionbal.PhotoRate_Shell[nelem][ion][ns][0],
121  0.05);
122 
123  /* outer loop is from K to most number of shells present in atom */
124  for( ns=0; ns < Heavy.nsShells[nelem][0]; ns++ )
125  {
126  fprintf( ioQQQ, "\n %s", elementnames.chElementNameShort[nelem] );
127  fprintf( ioQQQ, " %s" , Heavy.chShell[ns]);
128  /* MB hydrogenic photo rate may not be included in beow */
129  for( ion=0; ion < dense.IonHigh[nelem]; ion++ )
130  {
131  if( Heavy.nsShells[nelem][ion] > ns )
132  {
133  fprintf( ioQQQ, " %8.1e", ionbal.PhotoRate_Shell[nelem][ion][ns][0] );
134  }
135  else
136  {
137  break;
138  }
139  }
140  }
141  fprintf(ioQQQ,"\n");
142  }
143  return;
144 }

Generated for cloudy by doxygen 1.8.3.1