00001
00002
00003
00004 #include "cddefines.h"
00005 #include "opacity.h"
00006 #include "oxy.h"
00007 #include "thermal.h"
00008 #include "dense.h"
00009 #include "iso.h"
00010 #include "trace.h"
00011 #include "rfield.h"
00012 #include "atmdat.h"
00013 #include "atoms.h"
00014 #include "gammas.h"
00015 #include "ionbal.h"
00016
00017 void IonOxyge(void)
00018 {
00019 const int NDIM = ipOXYGEN+1;
00020
00021 static const double dicoef[2][NDIM] = {
00022 {1.11e-3,5.07e-3,1.48e-2,1.84e-2,4.13e-3,1.06e-1,6.23e-2,0.},
00023 {.0925,.181,.305,.1,.162,.34,.304,0.}
00024 };
00025 static const double dite[2][NDIM] = {
00026 {1.75e5,1.98e5,2.41e5,2.12e5,1.25e5,6.25e6,7.01e6,0.},
00027 {1.45e5,3.35e5,2.83e5,2.83e5,2.27e5,1.12e6,1.47e6,0.}
00028 };
00029 static const double ditcrt[NDIM] = {2.7e4,2.2e4,2.4e4,2.5e4,1.6e4,1.0e6,1.5e6,1e20};
00030 static const double aa[NDIM] = {0.,-0.0036,0.,0.0061,-2.8425,0.,0.,0.};
00031 static const double bb[NDIM] = {0.0238,0.7519,21.8790,0.2269,0.2283,0.,0.,0.};
00032 static const double cc[NDIM] = {0.0659,1.5252,16.2730,32.1419,40.4072,0.,0.,0.};
00033 static const double dd[NDIM] = {0.0349,-0.0838,-0.7020,1.9939,-3.4956,0.,0.,0.};
00034 static const double ff[NDIM] = {0.5334,0.2769,1.1899,-0.0646,1.7558,0.,0.,0.};
00035
00036 bool lgDebug = false;
00037 long int iup;
00038 double aeff;
00039
00040 DEBUG_ENTRY( "IonOxyge()" );
00041
00042
00043 if( !dense.lgElmtOn[ipOXYGEN] )
00044 {
00045 oxy.poiii2Max = 0.;
00046 oxy.poiii3Max = 0.;
00047 oxy.r4363Max = 0.;
00048 oxy.r5007Max = 0.;
00049 oxy.poiii2 = 0.;
00050 oxy.p1666 = 0.;
00051 oxy.AugerO3 = 0.;
00052 oxy.p1401 = 0.;
00053 oxy.s3727 = 0.;
00054 oxy.s7325 = 0.;
00055 thermal.heating[7][9] = 0.;
00056 oxy.poimax = 0.;
00057 return;
00058 }
00059
00060 ion_zero(ipOXYGEN);
00061
00062 ion_photo(ipOXYGEN,false);
00063
00064
00065 ion_collis(ipOXYGEN);
00066
00067
00068
00069 ion_recomb(false,(const double*)dicoef,(const double*)dite,ditcrt,aa,bb,cc,dd,ff,ipOXYGEN);
00070
00071
00072
00075 oxy.p1666 = ionbal.PhotoRate_Shell[ipOXYGEN][3][1][0];
00076
00077 oxy.p1401 = ionbal.PhotoRate_Shell[ipOXYGEN][2][1][0];
00078
00079 t_phoHeat dummy;
00080
00081
00082
00083
00084
00085
00086 oxy.d5007r = (realnum)(GammaK(opac.ipo3exc[0],opac.ipo3exc[1],
00087 opac.ipo3exc[2] , 1., &dummy ));
00088
00089
00090 oxy.d4363 = (realnum)(GammaK(opac.ipo3exc3[0],opac.ipo3exc3[1],
00091 opac.ipo3exc3[2] , 1., &dummy ));
00092
00093
00094 oxy.d6300 = (realnum)(GammaK(opac.ipo1exc[0],opac.ipo1exc[1],
00095 opac.ipo1exc[2] , 1., &dummy ));
00096
00097
00098 aeff = 0.0263 + oxy.d5007r;
00099
00100
00101 oxy.poiii2 = (realnum)(atom_pop2(2.5,9.,5.,aeff,2.88e4,1.)/aeff);
00102 {
00103
00104 enum {DEBUG_LOC=false};
00105
00106 if( DEBUG_LOC )
00107 {
00108 fprintf(ioQQQ,"pop rel %.1e rate %.1e grnd rate %.1e\n",
00109 oxy.poiii2 , oxy.d5007r ,ionbal.PhotoRate_Shell[ipOXYGEN][2][2][0] );
00110 }
00111 }
00112
00113
00114 if( nzone > 0 )
00115 {
00116
00117 ionbal.PhotoRate_Shell[ipOXYGEN][0][2][0] = ionbal.PhotoRate_Shell[ipOXYGEN][0][2][0]*
00118 (1. - oxy.poiexc) + oxy.d6300*oxy.poiexc;
00119
00120
00121 ionbal.PhotoRate_Shell[ipOXYGEN][2][2][0] = ionbal.PhotoRate_Shell[ipOXYGEN][2][2][0]*
00122 (1. - oxy.poiii2 - oxy.poiii3) + oxy.d5007r*oxy.poiii2 +
00123 oxy.d4363*oxy.poiii3;
00124
00125 if( ionbal.PhotoRate_Shell[ipOXYGEN][2][2][0] > 1e-30 && dense.IonLow[ipOXYGEN] <= 2 )
00126 {
00127 if( (oxy.d5007r*oxy.poiii2 + oxy.d4363*oxy.poiii3)/
00128 ionbal.PhotoRate_Shell[ipOXYGEN][2][2][0] > (oxy.r4363Max +
00129 oxy.r5007Max) )
00130 {
00131 oxy.poiii2Max = (realnum)(oxy.d5007r*oxy.poiii2/ionbal.PhotoRate_Shell[ipOXYGEN][2][2][0]);
00132 oxy.poiii3Max = (realnum)(oxy.d4363*oxy.poiii3/ionbal.PhotoRate_Shell[ipOXYGEN][2][2][0]);
00133 }
00134 oxy.r4363Max = (realnum)(MAX2(oxy.r4363Max,oxy.d4363));
00135 oxy.r5007Max = (realnum)(MAX2(oxy.r5007Max,oxy.d5007r));
00136 }
00137
00138
00139 if( dense.IonLow[ipOXYGEN] <= 0 && (ionbal.PhotoRate_Shell[ipOXYGEN][0][2][0] +
00140 atmdat.HCharExcIonOf[ipOXYGEN][0]*dense.xIonDense[ipHYDROGEN][1]) > 1e-30 )
00141 {
00142 oxy.poimax = (realnum)(MAX2(oxy.poimax,oxy.d6300*oxy.poiexc/
00143 (ionbal.PhotoRate_Shell[ipOXYGEN][0][2][0]+
00144 atmdat.HCharExcIonOf[ipOXYGEN][0]* dense.xIonDense[ipHYDROGEN][1])));
00145 }
00146 }
00147 else
00148 {
00149 oxy.poiii2Max = 0.;
00150 oxy.poiii3Max = 0.;
00151 oxy.r4363Max = 0.;
00152 oxy.r5007Max = 0.;
00153 oxy.poimax = 0.;
00154 }
00155
00156
00157 if( dense.IonLow[ipOXYGEN] == 0 && oxy.i2d < rfield.nflux )
00158 {
00159 oxy.s3727 = (realnum)(GammaK(oxy.i2d,oxy.i2p,opac.iopo2d , 1., &dummy ));
00160
00161 iup = MIN2(iso.ipIsoLevNIonCon[ipH_LIKE][1][0],rfield.nflux);
00162 oxy.s7325 = (realnum)(GammaK(oxy.i2d,iup,opac.iopo2d , 1., &dummy ));
00163
00164 oxy.s7325 -= oxy.s3727;
00165 oxy.s3727 = oxy.s3727 + oxy.s7325;
00166
00167
00168 oxy.s7325 *= 0.66f;
00169 }
00170 else
00171 {
00172 oxy.s3727 = 0.;
00173 oxy.s7325 = 0.;
00174 }
00175
00176 oxy.AugerO3 = (realnum)ionbal.PhotoRate_Shell[ipOXYGEN][0][0][0];
00177
00178
00179 if(0 && nzone > 100 )
00180 lgDebug = true;
00181 else
00182 lgDebug = false;
00183 ion_solver(ipOXYGEN,lgDebug);
00184 if( lgDebug )
00185 fprintf(ioQQQ,"DEBUG O\t%.3e\t%.3e\tH\t%.3e\t%.3e\n",
00186 dense.xIonDense[ipOXYGEN][0],
00187 dense.xIonDense[ipOXYGEN][1],
00188 dense.xIonDense[ipHYDROGEN][0],
00189 dense.xIonDense[ipHYDROGEN][1]);
00190
00191
00192 oxy.p1666 *= dense.xIonDense[ipOXYGEN][1]*0.3;
00193 oxy.p1401 *= dense.xIonDense[ipOXYGEN][2]*0.43;
00194 oxy.s3727 *= dense.xIonDense[ipOXYGEN][0];
00195 oxy.s7325 *= dense.xIonDense[ipOXYGEN][0];
00196 oxy.AugerO3 *= dense.xIonDense[ipOXYGEN][0];
00197
00198 if( trace.lgTrace )
00199 {
00200 fprintf( ioQQQ, " IonOxyge returns; frac=" );
00201 for( int i=1; i <= 9; i++ )
00202 {
00203 fprintf( ioQQQ, " %10.3e", dense.xIonDense[ipOXYGEN][i-1]/
00204 dense.gas_phase[ipOXYGEN] );
00205 }
00206 fprintf( ioQQQ, "\n" );
00207 }
00208 return;
00209 }