00001
00002
00003
00004 #include "cddefines.h"
00005 #include "dense.h"
00006 #include "fe.h"
00007 #include "thermal.h"
00008 #include "gammas.h"
00009 #include "opacity.h"
00010 #include "trace.h"
00011 #include "grainvar.h"
00012 #include "ionbal.h"
00013 #include "mole.h"
00014
00015
00016 void IonIron(void)
00017 {
00018 const int NDIM = ipIRON+1;
00019
00020 static const double dicoef[2][NDIM] = {
00021 {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,
00022 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.},
00023 {.456,.323,.310,.411,.359,.0975,.229,4.20,3.30,5.30,1.50,0.700,.600,.5,1.,0.,7.8,
00024 6.3,5.5,3.6,4.9,1.6,4.2,.284,.279,0.}
00025 };
00026 static const double dite[2][NDIM] = {
00027 {6.00e4,1.94e5,3.31e5,4.32e5,6.28e5,7.50e5,7.73e5,2.62e5,2.50e5,2.57e5,2.84e5,
00028 8.69e5,4.21e5,4.57e5,2.85e5,8.18e5,1.51e6,1.30e6,1.19e6,1.09e6,9.62e5,7.23e5,
00029 4.23e5,5.84e7,6.01e7,0.},
00030 {8.97e4,1.71e5,2.73e5,3.49e5,5.29e5,4.69e5,6.54e5,1.32e6,1.33e6,1.41e6,1.52e6,
00031 1.51e6,1.82e6,1.84e6,2.31e6,0.,9.98e6,9.98e6,1.00e7,1.10e7,8.34e6,1.01e7,1.07e7,
00032 1.17e7,9.97e6,0.}
00033 };
00034 static const double ditcrt[NDIM] = {1e11,1e11,1e11,1e11,1e11,1e11,1e11,
00035 1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,1e11,
00036 1e11,1e11,1e11,1e11,1e11,1e11,1e11};
00037 static const double aa[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
00038 0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
00039 static const double bb[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
00040 0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
00041 static const double cc[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
00042 0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
00043 static const double dd[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
00044 0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
00045
00046
00047
00048 static const double ff[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,
00049 0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.};
00050 static const double fyield[NDIM+1] = {.34,.34,.35,.35,.36,.37,.37,.38,.39,.40,
00051 .41,.42,.43,.44,.45,.46,.47,.47,.48,.48,.49,.49,.11,.75,0.,0.,0.};
00052
00053 long int i, limit, limit2;
00054
00055 DEBUG_ENTRY( "IonIron()" );
00056
00057
00058
00059
00060
00061
00062
00063
00064 if( !dense.lgElmtOn[ipIRON] )
00065 {
00066 fe.fekcld = 0.;
00067 fe.fekhot = 0.;
00068 fe.fegrain = 0.;
00069 return;
00070 }
00071
00072 ion_zero(ipIRON);
00073
00074 ion_photo(ipIRON,false);
00075
00076
00077
00078
00079
00080 {
00081
00082 enum {DEBUG_LOC=false};
00083
00084 if( DEBUG_LOC )
00085 {
00086 long int iplow , iphi , ipop , ns , ion;
00087 double rate;
00088 ns = 6;
00089 ion = 0;
00090
00091 iplow = opac.ipElement[ipIRON][ion][ns][0];
00092 iphi = opac.ipElement[ipIRON][ion][ns][1];
00093 ipop = opac.ipElement[ipIRON][ion][ns][2];
00094 rate = ionbal.PhotoRate_Shell[ipIRON][ion][ns][0];
00095 GammaPrt( iplow , iphi , ipop , ioQQQ, rate , rate*0.1 );
00096 }
00097 }
00098
00099
00100 ion_collis(ipIRON);
00101
00102
00103 ion_recomb(false,(const double*)dicoef,(const double*)dite,ditcrt,aa,bb,cc,dd,ff,ipIRON);
00104
00105
00106
00107
00108 if( !co.lgNoCOMole )
00109 {
00110 ionbal.PhotoRate_Shell[ipIRON][0][6][0] +=
00111 CO_findrk("S+,Fe=>S,Fe+")*dense.xIonDense[ipSULPHUR][1] +
00112 CO_findrk("Si+,Fe=>Si,Fe+")*dense.xIonDense[ipSILICON][1] +
00113 CO_findrk("C+,Fe=>C,Fe+")*dense.xIonDense[ipCARBON][1];
00114
00115
00116 }
00117
00118
00119 ion_solver(ipIRON,false);
00120
00121
00122
00123 fe.fekcld = 0.;
00124 limit = MIN2(18,dense.IonHigh[ipIRON]);
00125
00126 for( i=dense.IonLow[ipIRON]; i < limit; i++ )
00127 {
00128 ASSERT( i < NDIM + 1 );
00129 fe.fekcld +=
00130 (realnum)(ionbal.PhotoRate_Shell[ipIRON][i][0][0]*dense.xIonDense[ipIRON][i]*
00131 fyield[i]);
00132 }
00133
00134
00135 fe.fekhot = 0.;
00136 limit = MAX2(18,dense.IonLow[ipIRON]);
00137
00138 limit2 = MIN2(ipIRON+1,dense.IonHigh[ipIRON]);
00139 ASSERT( limit2 <= LIMELM + 1 );
00140
00141 for( i=limit; i < limit2; i++ )
00142 {
00143 ASSERT( i < NDIM + 1 );
00144 fe.fekhot +=
00145 (realnum)(ionbal.PhotoRate_Shell[ipIRON][i][0][0]*dense.xIonDense[ipIRON][i]*
00146 fyield[i]);
00147 }
00148
00149
00150
00151 i = 0;
00152
00153 fe.fegrain = ( gv.lgWD01 ) ? 0.f : (realnum)(ionbal.PhotoRate_Shell[ipIRON][i][0][0]*fyield[i]*
00154 gv.elmSumAbund[ipIRON]);
00155
00156 if( trace.lgTrace && trace.lgHeavyBug )
00157 {
00158
00159 fprintf( ioQQQ, " Fe" );
00160 for( i=0; i < 15; i++ )
00161 {
00162 fprintf( ioQQQ, "\t%.1e", dense.xIonDense[ipIRON][i]/dense.gas_phase[ipIRON] );
00163 }
00164 fprintf( ioQQQ, "\n" );
00165 }
00166
00167 if( trace.lgTrace && trace.lgFeBug )
00168 {
00169 fprintf( ioQQQ, " IonIron-Abund:" );
00170 for( i=0; i < 27; i++ )
00171 {
00172 fprintf( ioQQQ, "%10.2e", dense.xIonDense[ipIRON][i] );
00173 }
00174 fprintf( ioQQQ, "\n" );
00175
00176 fprintf( ioQQQ, " IonIron - Ka(Auger)%10.2e\n", fe.fekcld +
00177 fe.fekhot );
00178 }
00179 return;
00180 }