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