cloudy trunk
|
00001 /* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and 00002 * others. For conditions of distribution and use see copyright notice in license.txt */ 00003 /*IonSulph compute ionization balance for sulphur */ 00004 #include "cddefines.h" 00005 #include "dense.h" 00006 #include "ionbal.h" 00007 00008 void IonSulph(void) 00009 { 00010 const int NDIM = ipSULPHUR+1; 00011 00012 static const double dicoef[2][NDIM] = { 00013 {1.62e-3,1.09e-2,3.35e-2,3.14e-2,1.27e-2,1.47e-2,1.34e-2, 00014 2.38e-2,3.19e-2,7.13e-2,.08,.0796,.0134,.402,.241,0.}, 00015 {0.,1.20e-2,6.59e-2,6.89e-2,.187,.129,1.04,1.12,1.40,1.00, 00016 .555,1.63,.304,.298,.281,0.} 00017 }; 00018 static const double dite[2][NDIM] = { 00019 {1.25e5,1.92e5,1.89e5,1.68e5,1.38e5,1.80e6,6.90e5,5.84e5, 00020 5.17e5,6.66e5,6.00e5,5.09e5,2.91e5,2.41e7,2.54e7,0.}, 00021 {0.,1.80e4,1.59e5,8.04e4,1.71e5,1.75e6,2.15e6,2.59e6,2.91e6, 00022 2.32e6,2.41e6,6.37e6,1.04e6,4.67e6,5.30e6,0.} 00023 }; 00024 static const double aa[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0., 00025 0.,0.,0.,0.}; 00026 static const double bb[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0., 00027 0.,0.,0.,0.}; 00028 static const double cc[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0., 00029 0.,0.,0.,0.}; 00030 static const double dd[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0., 00031 0.,0.,0.,0.}; 00032 static const double ff[NDIM] = {0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0.,0., 00033 0.,0.,0.,0.}; 00034 static const double ditcrt[NDIM] = {2.2e4,1.2e4,1.4e4,1.5e4,1.4e4,2.9e5, 00035 1.3e5,1.1e5,9.0e4,9.0e4,9.0e4,8.3e4,6.0e4,5.0e6,9.0e6,1e20}; 00036 00037 bool lgPrtDebug=false; 00038 double save_rec; 00039 00040 DEBUG_ENTRY( "IonSulph()" ); 00041 00042 /* sulphur, atomic number 16 */ 00043 if( !dense.lgElmtOn[ipSULPHUR] ) 00044 { 00045 return; 00046 } 00047 00048 ion_zero(ipSULPHUR); 00049 00050 ion_photo(ipSULPHUR,lgPrtDebug); 00051 00052 /* find collisional ionization rates */ 00053 ion_collis(ipSULPHUR); 00054 00055 /* get recombination coefficients */ 00056 /*lint -e64 type mismatch */ 00057 ion_recomb(false,(const double*)dicoef,(const double*)dite,ditcrt,aa,bb,cc,dd,ff,ipSULPHUR); 00058 /*lint +e64 type mismatch */ 00059 00060 /* >>chng 03 sep 29, synch up ion and co solvers. 00061 * the co solver will have a different S/S+ balance that 00062 * is strongly affected by the chemistry if we are in neutral gas 00063 * (hence the test that the highest stage of ionization is <=2 on 00064 * physics scale) 00065 * in this case use co solver's S/S+ balance */ 00066 save_rec = ionbal.RateRecomTot[ipSULPHUR][0]; 00067 /*>>chng 04 apr 27, do not test for ionhigh being 1, 00068 * no reason for test on upper stage of ionization, 00069 * if codrive called but not evaluted then hevmol is all zero */ 00070 /* >>chng 04 sep 10, rm check on search phase, no reason for it */ 00071 # if 0 00072 if( dense.IonLow[ipSULPHUR]==0 && 00073 co.hevmol[ipSP] > SMALLFLOAT && 00074 ionbal.RateIonizTot[ipSULPHUR][0]*co.hevmol[ipATS]> 0. ) 00075 { 00076 ionbal.RateRecomTot[ipSULPHUR][0] = 00077 ionbal.RateIonizTot[ipSULPHUR][0]* 00078 co.hevmol[ipATS]/co.hevmol[ipSP]; 00079 } 00080 # endif 00081 00082 /* solve for ionization balance */ 00083 /*if( nzone>700 ) lgPrtDebug = true;*/ 00084 ion_solver(ipSULPHUR,lgPrtDebug); 00085 00086 /* reset the var we just hosed */ 00087 /* if ion just trimed down, old rec may be zero, must not use it*/ 00088 if( save_rec > 0. ) 00089 ionbal.RateRecomTot[ipSULPHUR][0] = save_rec; 00090 return; 00091 }