00001
00002
00003
00004 #include "cddefines.h"
00005 #include "coolheavy.h"
00006 #include "dense.h"
00007 #include "taulines.h"
00008 #include "h2.h"
00009 #include "phycon.h"
00010 #include "embesq.h"
00011 #include "hmi.h"
00012 #include "oxy.h"
00013 #include "colden.h"
00014 #include "mole.h"
00015 #include "ligbar.h"
00016 #include "thermal.h"
00017 #include "lines_service.h"
00018 #include "atoms.h"
00019 #include "cooling.h"
00020
00021 #if defined (__ICC) && defined(__ia64) && __INTEL_COMPILER < 910
00022 #pragma optimization_level 1
00023 #endif
00024 void CoolOxyg(void)
00025 {
00026 double a21,
00027 a31,
00028 a32,
00029 a41,
00030 a42,
00031 a43,
00032 a51,
00033 a52,
00034 a53,
00035 a54,
00036 a6300,
00037 a6363,
00038 aeff,
00039 cs2s2p,
00040 cs2s3p,
00041 p[5],
00042 r12 ,
00043 r13;
00044
00045 static double go2[5]={4.,6.,4.,4.,2.};
00046 static double exo2[4]={26808.4,21.0,13637.5,1.5};
00047
00048
00049 static double oi_cs_tsave=-1. , oi_te_tsave=-1. , oi_dcdt_tsave=-1.;
00050 long int i;
00051 double rate_OH_dissoc;
00052
00053 DEBUG_ENTRY( "CoolOxyg()" );
00054
00055 {
00056 enum{DEBUG_LOC=false};
00057 if( DEBUG_LOC )
00058 {
00059
00060 double TeTest = phycon.TEMP_LIMIT_LOW;
00061 double oi_a, oi_b, oi_c, oi_d, oi_e, oi_f;
00062 double oii_a,
00063 oii_b,
00064 oii_c,
00065 oii_d,
00066 oii_e,
00067 oii_f,
00068 oii_g,
00069 oii_h,
00070 oii_i,
00071 oii_j,
00072 oii_k;
00073 double oiii_a,
00074 oiii_b,
00075 oiii_c,
00076 oiii_d,
00077 oiii_e,
00078 oiii_f,
00079 oiii_g,
00080 oiii_h,
00081 oiii_i,
00082 oiii_j;
00083 double oiv_a,oiv_b;
00084 double ov_a,ov_b;
00085 double sii_a,
00086 sii_b,
00087 sii_c,
00088 sii_d,
00089 sii_e,
00090 sii_f,
00091 sii_g,
00092 sii_h,
00093 sii_i,
00094 sii_j,
00095 sii_k;
00096 double siii_a,
00097 siii_b,
00098 siii_c,
00099 siii_d,
00100 siii_e,
00101 siii_f,
00102 siii_g,
00103 siii_h,
00104 siii_i,
00105 siii_j,
00106 siii_k,
00107 siii_l;
00108 double siv_a;
00109 double sviii_a;
00110 double neiii_a,neiii_b,neiii_c,neiii_d,neiii_e;
00111 TempChange( TeTest , true );
00112 oi_cs(oi_a, oi_b, oi_c, oi_d, oi_e, oi_f);
00113 oii_cs(oii_a, oii_b, oii_c, oii_d, oii_e, oii_f, oii_g,
00114 oii_h, oii_i, oii_j, oii_k);
00115 oiii_cs(oiii_a, oiii_b, oiii_c, oiii_d, oiii_e, oiii_f,
00116 oiii_g, oiii_h, oiii_i, oiii_j);
00117 oiv_cs(oiv_a, oiv_b);
00118 ov_cs(ov_a, ov_b);
00119 sii_cs(sii_a, sii_b, sii_c, sii_d, sii_e, sii_f, sii_g,
00120 sii_h, sii_i, sii_j, sii_k);
00121 siii_cs(siii_a, siii_b, siii_c, siii_d, siii_e, siii_f,
00122 siii_g, siii_h, siii_i, siii_j, siii_k,siii_l);
00123 siv_cs(siv_a);
00124 sviii_cs(sviii_a);
00125 neiii_cs(neiii_a,neiii_b,neiii_c,neiii_d,neiii_e);
00126 double tinc = 0.05;
00127 for( TeTest=phycon.TEMP_LIMIT_LOW; TeTest<phycon.TEMP_LIMIT_HIGH;
00128 TeTest *=(1+tinc) )
00129 {
00130 double oi_aold=oi_a,
00131 oi_bold=oi_b,
00132 oi_cold=oi_c,
00133 oi_dold=oi_d,
00134 oi_eold=oi_e,
00135 oi_fold=oi_f;
00136 double oii_aold=oii_a,
00137 oii_bold=oii_b,
00138 oii_cold=oii_c,
00139 oii_dold=oii_d,
00140 oii_eold=oii_e,
00141 oii_fold=oii_f,
00142 oii_gold=oii_g,
00143 oii_hold=oii_h,
00144 oii_iold=oii_i,
00145 oii_jold=oii_j,
00146 oii_kold=oii_k;
00147 double oiii_aold=oiii_a,
00148 oiii_bold=oiii_b,
00149 oiii_cold=oiii_c,
00150 oiii_dold=oiii_d,
00151 oiii_eold=oiii_e,
00152 oiii_fold=oiii_f,
00153 oiii_gold=oiii_g,
00154 oiii_hold=oiii_h,
00155 oiii_iold=oiii_i,
00156 oiii_jold=oiii_j;
00157 double oiv_aold=oiv_a,oiv_bold=oiv_b;
00158 double ov_aold=ov_a,ov_bold=ov_b;
00159 double sii_aold=sii_a,
00160 sii_bold=sii_b,
00161 sii_cold=sii_c,
00162 sii_dold=sii_d,
00163 sii_eold=sii_e,
00164 sii_fold=sii_f,
00165 sii_gold=sii_g,
00166 sii_hold=sii_h,
00167 sii_iold=sii_i,
00168 sii_jold=sii_j,
00169 sii_kold=sii_k;
00170 double siii_aold=siii_a,
00171 siii_bold=siii_b,
00172 siii_cold=siii_c,
00173 siii_dold=siii_d,
00174 siii_eold=siii_e,
00175 siii_fold=siii_f,
00176 siii_gold=siii_g,
00177 siii_hold=siii_h,
00178 siii_iold=siii_i,
00179 siii_jold=siii_j,
00180 siii_kold=siii_k,
00181 siii_lold=siii_l;
00182 double siv_aold=siv_a;
00183 double sviii_aold=sviii_a;
00184 double neiii_aold=neiii_a,
00185 neiii_bold=neiii_b,
00186 neiii_cold=neiii_c,
00187 neiii_dold=neiii_d,
00188 neiii_eold=neiii_e;
00189 TempChange( TeTest , true );
00190 oi_cs(oi_a, oi_b, oi_c, oi_d, oi_e, oi_f);
00191 ASSERT( fabs(oi_a-oi_aold) <= oi_a*3.*tinc);
00192 ASSERT( fabs(oi_b-oi_bold) <= oi_b*3.*tinc);
00193 ASSERT( fabs(oi_c-oi_cold) <= oi_c*3.*tinc);
00194 ASSERT( fabs(oi_d-oi_dold) <= oi_d*3.*tinc);
00195 ASSERT( fabs(oi_e-oi_eold) <= oi_e*3.*tinc);
00196 ASSERT( fabs(oi_f-oi_fold) <= oi_f*3.*tinc);
00197 oii_cs(oii_a, oii_b, oii_c, oii_d, oii_e, oii_f,
00198 oii_g, oii_h, oii_i, oii_j, oii_k);
00199 ASSERT( fabs(oii_a-oii_aold) <= oii_a*3.*tinc);
00200 ASSERT( fabs(oii_b-oii_bold) <= oii_b*3.*tinc);
00201 ASSERT( fabs(oii_c-oii_cold) <= oii_c*3.*tinc);
00202 ASSERT( fabs(oii_d-oii_dold) <= oii_d*3.*tinc);
00203 ASSERT( fabs(oii_e-oii_eold) <= oii_e*3.*tinc);
00204 ASSERT( fabs(oii_f-oii_fold) <= oii_f*3.*tinc);
00205 ASSERT( fabs(oii_g-oii_gold) <= oii_g*3.*tinc);
00206 ASSERT( fabs(oii_h-oii_hold) <= oii_h*3.*tinc);
00207 ASSERT( fabs(oii_i-oii_iold) <= oii_i*3.*tinc);
00208 ASSERT( fabs(oii_j-oii_jold) <= oii_j*3.*tinc);
00209 ASSERT( fabs(oii_k-oii_kold) <= oii_k*3.*tinc);
00210 oiii_cs(oiii_a, oiii_b, oiii_c, oiii_d, oiii_e,
00211 oiii_f, oiii_g, oiii_h, oiii_i, oiii_j);
00212 ASSERT( fabs(oiii_a-oiii_aold) <= oiii_a*3.*tinc);
00213 ASSERT( fabs(oiii_b-oiii_bold) <= oiii_b*3.*tinc);
00214 ASSERT( fabs(oiii_c-oiii_cold) <= oiii_c*3.*tinc);
00215 ASSERT( fabs(oiii_d-oiii_dold) <= oiii_d*3.*tinc);
00216 ASSERT( fabs(oiii_e-oiii_eold) <= oiii_e*3.*tinc);
00217 ASSERT( fabs(oiii_f-oiii_fold) <= oiii_f*3.*tinc);
00218 ASSERT( fabs(oiii_g-oiii_gold) <= oiii_g*3.*tinc);
00219 ASSERT( fabs(oiii_h-oiii_hold) <= oiii_h*3.*tinc);
00220 ASSERT( fabs(oiii_i-oiii_iold) <= oiii_i*3.*tinc);
00221 ASSERT( fabs(oiii_j-oiii_jold) <= oiii_j*3.*tinc);
00222 oiv_cs(oiv_a, oiv_b);
00223 ASSERT( fabs(oiv_a-oiv_aold) <= oiv_a*3.*tinc);
00224 ASSERT( fabs(oiv_b-oiv_bold) <= oiv_b*3.*tinc);
00225 ov_cs(ov_a, ov_b);
00226 ASSERT( fabs(ov_a-ov_aold) <= ov_a*3.*tinc);
00227 ASSERT( fabs(ov_b-ov_bold) <= ov_b*3.*tinc);
00228 sii_cs(sii_a, sii_b, sii_c, sii_d, sii_e, sii_f,
00229 sii_g, sii_h, sii_i, sii_j, sii_k);
00230 ASSERT( fabs(sii_a-sii_aold) <= sii_a*3.*tinc);
00231 ASSERT( fabs(sii_b-sii_bold) <= sii_b*3.*tinc);
00232 ASSERT( fabs(sii_c-sii_cold) <= sii_c*3.*tinc);
00233 ASSERT( fabs(sii_d-sii_dold) <= sii_d*3.*tinc);
00234 ASSERT( fabs(sii_e-sii_eold) <= sii_e*3.*tinc);
00235 ASSERT( fabs(sii_f-sii_fold) <= sii_f*3.*tinc);
00236 ASSERT( fabs(sii_g-sii_gold) <= sii_g*3.*tinc);
00237 ASSERT( fabs(sii_h-sii_hold) <= sii_h*3.*tinc);
00238 ASSERT( fabs(sii_i-sii_iold) <= sii_i*3.*tinc);
00239 ASSERT( fabs(sii_j-sii_jold) <= sii_j*3.*tinc);
00240 ASSERT( fabs(sii_k-sii_kold) <= sii_k*3.*tinc);
00241 siii_cs(siii_a,
00242 siii_b,
00243 siii_c,
00244 siii_d,
00245 siii_e,
00246 siii_f,
00247 siii_g,
00248 siii_h,
00249 siii_i,
00250 siii_j,
00251 siii_k,
00252 siii_l);
00253 ASSERT( fabs(siii_a-siii_aold) <= siii_a*3.*tinc);
00254 ASSERT( fabs(siii_b-siii_bold) <= siii_b*3.*tinc);
00255 ASSERT( fabs(siii_c-siii_cold) <= siii_c*3.*tinc);
00256 ASSERT( fabs(siii_d-siii_dold) <= siii_d*3.*tinc);
00257 ASSERT( fabs(siii_e-siii_eold) <= siii_e*3.*tinc);
00258 ASSERT( fabs(siii_f-siii_fold) <= siii_f*3.*tinc);
00259 ASSERT( fabs(siii_g-siii_gold) <= siii_g*3.*tinc);
00260 ASSERT( fabs(siii_h-siii_hold) <= siii_h*3.*tinc);
00261 ASSERT( fabs(siii_i-siii_iold) <= siii_i*3.*tinc);
00262 ASSERT( fabs(siii_j-siii_jold) <= siii_j*3.*tinc);
00263 ASSERT( fabs(siii_k-siii_kold) <= siii_k*3.*tinc);
00264 ASSERT( fabs(siii_l-siii_lold) <= siii_l*3.*tinc);
00265 siv_cs(siv_a);
00266 ASSERT( fabs(siv_a-siv_aold) <= siv_a*3.*tinc);
00267 sviii_cs(sviii_a);
00268 ASSERT( fabs(sviii_a-sviii_aold) <= sviii_a*3.*tinc);
00269 neiii_cs(neiii_a, neiii_b, neiii_c, neiii_d, neiii_e);
00270 ASSERT( fabs(neiii_a-neiii_aold) <= neiii_a*3.*tinc);
00271 ASSERT( fabs(neiii_b-neiii_bold) <= neiii_b*3.*tinc);
00272 ASSERT( fabs(neiii_c-neiii_cold) <= neiii_c*3.*tinc);
00273 ASSERT( fabs(neiii_d-neiii_dold) <= neiii_d*3.*tinc);
00274 ASSERT( fabs(neiii_e-neiii_eold) <= neiii_e*3.*tinc);
00275 }
00276 cdEXIT(EXIT_SUCCESS);
00277 }
00278 }
00279
00280
00281
00282
00283
00284 atom_oi_calc(&CoolHeavy.coolOi);
00285 CoolAdd("o 1",8446,CoolHeavy.coolOi);
00286 double oi_cs3P23P1,
00287 oi_cs3P23P0,
00288 oi_cs3P13P0,
00289 oi_cs3P1D2,
00290 oi_cs1D21S0,
00291 oi_cs3P1S0;
00292
00293
00294
00295
00296
00297
00298
00299
00300
00301
00302
00303 oi_cs(oi_cs3P23P1,
00304 oi_cs3P23P0,
00305 oi_cs3P13P0,
00306 oi_cs3P1D2,
00307 oi_cs1D21S0,
00308 oi_cs3P1S0);
00309 double csh01=-1.,
00310 cshe01=-1.,
00311 csh201=-1.,
00312 csh12=-1.,
00313 cshe12=-1.,
00314 csh212=-1.,
00315 csh02=-1.,
00316 cshe02=-1.,
00317 csh202 =-1.,
00318 csh2o01=-1.,
00319 csh2o02=-1.,
00320 csh2o12=-1.,
00321 csh2p01=-1.,
00322 csh2p02=-1.,
00323 csh2p12=-1.,
00324 csp01=-1.,
00325 csp02=-1.,
00326 csp12=-1.;
00327
00328
00329
00330
00331
00332
00333
00334
00335
00336
00337 oi_othercs(csh01,cshe01,csh201,csh12,cshe12,csh212,csh02,cshe02,csh202,
00338 csh2o01,csh2o02,csh2o12,csh2p01,csh2p02,csh2p12,csp01,csp02,csp12);
00339
00340 oi_cs3P23P1 = oi_cs3P23P1+csp01+3.*(csh01*dense.xIonDense[ipHYDROGEN][0]
00341 + cshe01*dense.xIonDense[ipHELIUM][0] + csh201*hmi.H2_total)/
00342 dense.cdsqte;
00343 oi_cs3P13P0 = oi_cs3P13P0+csp12+(csh12*dense.xIonDense[ipHYDROGEN][0] +
00344 cshe12*dense.xIonDense[ipHELIUM][0] + csh212*hmi.H2_total)/
00345 dense.cdsqte;
00346 oi_cs3P23P0 = oi_cs3P23P0+csp02+(csh02*dense.xIonDense[ipHYDROGEN][0] +
00347 cshe02*dense.xIonDense[ipHELIUM][0] + csh202*hmi.H2_total)/
00348 dense.cdsqte;
00349
00350
00351
00352
00353 a6300 = TauLines[ipT6300].Emis->Aul*TauLines[ipT6300].Emis->Pesc;
00354 TauLines[ipT6300].Emis->PopOpc = (dense.xIonDense[ipOXYGEN][0]*5./5.);
00355 TauLines[ipT6300].Lo->Pop = (dense.xIonDense[ipOXYGEN][0]*5./5.);
00356 TauLines[ipT6300].Hi->Pop = 0.;
00357 TauLines[ipT6300].Coll.col_str = (realnum)(oi_cs3P1D2*5./9.);
00358 TauLines[ipT6363].Emis->PopOpc = (dense.xIonDense[ipOXYGEN][0]*5./3.);
00359 TauLines[ipT6363].Lo->Pop = (dense.xIonDense[ipOXYGEN][0]*5./3.);
00360 TauLines[ipT6363].Hi->Pop = 0.;
00361 TauLines[ipT6363].Coll.col_str = (realnum)(oi_cs3P1D2*3./9.);
00362 a6363 = TauLines[ipT6363].Emis->Aul*TauLines[ipT6363].Emis->Pesc;
00363 a21 = a6300 + a6363 + oxy.d6300;
00364 a32 = TauLines[ipT5577].Emis->Aul*TauLines[ipT5577].Emis->Pesc;
00365 PutCS(oi_cs1D21S0,&TauLines[ipT5577]);
00366
00367
00368
00369 rate_OH_dissoc = CO_findrate("PHOTON,OH=>O,H");
00370 r12 = rate_OH_dissoc*0.55/SDIV( dense.xIonDense[ipOXYGEN][0] );
00371 r13 = rate_OH_dissoc*0.05/SDIV( dense.xIonDense[ipOXYGEN][0] );
00372
00373 CoolHeavy.c6300_frac_emit = (a6300+a6363)/(a6300+a6363+oi_cs3P1D2*dense.cdsqte/5.);
00374 CoolHeavy.c5577_frac_emit = (a32)/(a32+oi_cs1D21S0*dense.cdsqte/3.);
00375
00376
00377 CoolHeavy.c5577 = atom_pop3(9.,5.,1.,oi_cs3P1D2,oi_cs3P1S0,oi_cs1D21S0,a21,7.56e-2,a32,
00378 22590.,25775.7,&oxy.poiexc,dense.xIonDense[ipOXYGEN][0],0.,r12,r13)*a32*
00379 3.57e-12;
00380 TauLines[ipT5577].Emis->PopOpc = oxy.poiexc;
00381 TauLines[ipT5577].Lo->Pop = oxy.poiexc;
00382 TauLines[ipT5577].Hi->Pop = 0.;
00383
00384
00385 CoolHeavy.c5577 *= (1.-r13/(SDIV(atoms.c13)));
00386 CoolHeavy.c6300 = oxy.poiexc*a6300*TauLines[ipT6300].EnergyErg *
00387 (1.-r12/(SDIV(atoms.c12)));
00388 CoolHeavy.c6363 = oxy.poiexc*a6363*TauLines[ipT6363].EnergyErg *
00389 (1.-r12/(SDIV(atoms.c12)));
00390
00391 thermal.dCooldT += CoolHeavy.c6300*(2.28e4*thermal.tsq1 + thermal.halfte) *
00392
00393 (1.-r12/(SDIV(atoms.c12)));
00394 oxy.poiexc /= (realnum)MAX2(1e-20,dense.xIonDense[ipOXYGEN][0]);
00395 CoolAdd("o 1",5577,CoolHeavy.c5577);
00396 CoolAdd("o 1",6300,CoolHeavy.c6300);
00397 CoolAdd("o 1",6363,CoolHeavy.c6363);
00398 PutCS(oi_cs3P23P1,&TauLines[ipT63]);
00399 PutCS(oi_cs3P13P0,&TauLines[ipT146]);
00400 PutCS(oi_cs3P23P0,&TauDummy);
00401 atom_level3(&TauLines[ipT63],&TauLines[ipT146],&TauDummy);
00402
00403 for( i=0; i<3; ++i)
00404 colden.O1Pops[i] = (realnum)atoms.PopLevels[i];
00405
00406 if( !fp_equal(phycon.te,oi_te_tsave) )
00407 {
00408
00409 if(oi_te_tsave>0. )
00410 oi_dcdt_tsave = (oi_cs3P23P1-oi_cs_tsave) /
00411 (phycon.te-oi_te_tsave);
00412 else
00413 oi_dcdt_tsave = 0.;
00414 oi_cs_tsave = oi_cs3P23P1;
00415 oi_te_tsave = phycon.te;
00416
00417
00418
00419
00420 oi_dcdt_tsave = MAX2( 0. , oi_dcdt_tsave);
00421
00422 oi_dcdt_tsave = MIN2( TauLines[ipT63].EnergyK*thermal.tsq1*4.,oi_dcdt_tsave);
00423 }
00424
00425
00426 thermal.dCooldT += TauLines[ipT63].Coll.cool*oi_dcdt_tsave;
00427
00428 {
00429 enum{DEBUG_LOC=false};
00430 if( DEBUG_LOC )
00431 {
00432 fprintf(ioQQQ,"DEBUG OI\t%.2f\tte\t%.5e\tcool\t%.5e\tcs\t%.4e\told\t%.4e\tnew\t%.4e\n",
00433 fnzone,
00434 phycon.te,
00435 TauLines[ipT63].Coll.cool ,
00436 TauLines[ipT63].Coll.col_str ,
00437 TauLines[ipT63].Coll.cool*TauLines[ipT63].EnergyK*thermal.tsq1,
00438 TauLines[ipT63].Coll.cool*oi_dcdt_tsave );
00439 }
00440 }
00441
00442
00443
00444
00445 double oii_cs4S32D5,
00446 oii_cs4S32D3,
00447 oii_cs2D52D3,
00448 oii_cs4S32P3,
00449 oii_cs2D52P3,
00450 oii_cs2D32P3,
00451 oii_cs4S32P1,
00452 oii_cs2D52P1,
00453 oii_cs2D32P1,
00454 oii_cs2P32P1,
00455 oii_cs4S34P;
00456
00457
00458
00459
00460
00461
00462 a21 = 4.12e-5;
00463 a31 = 1.63e-4;
00464 a32 = 1.24e-7;
00465 a41 = 5.65e-2;
00466 a42 = 1.11e-1;
00467 a43 = 5.87e-2;
00468 a51 = 2.27e-2;
00469 a52 = 5.82e-2;
00470 a53 = 9.67e-2;
00471 a54 = 3.15e-10;
00472
00473
00474
00475
00476
00477
00478
00479 oii_cs(oii_cs4S32D5,
00480 oii_cs4S32D3,
00481 oii_cs2D52D3,
00482 oii_cs4S32P3,
00483 oii_cs2D52P3,
00484 oii_cs2D32P3,
00485 oii_cs4S32P1,
00486 oii_cs2D52P1,
00487 oii_cs2D32P1,
00488 oii_cs2P32P1,
00489 oii_cs4S34P);
00490
00491 double Cooling , CoolingDeriv;
00492 atom_pop5(go2,exo2,oii_cs4S32D5,oii_cs4S32D3,oii_cs4S32P3,oii_cs4S32P1,
00493 oii_cs2D52D3,oii_cs2D52P3,oii_cs2D52P1,oii_cs2D32P3,oii_cs2D32P1,
00494 oii_cs2P32P1,a21,a31,a41,a51,a32,a42,a52,a43,a53,a54,p,
00495 dense.xIonDense[ipOXYGEN][1], &Cooling , &CoolingDeriv, 0.,0.,0.,0.);
00496
00497 CoolHeavy.O3730 = (realnum)(p[1]*a21*5.34e-12);
00498 CoolHeavy.O3726 = (realnum)(p[2]*a31*5.34e-12);
00499 CoolHeavy.O2471 = (realnum)((p[3]*a41 + p[4]*a51)*8.05e-12);
00500 CoolHeavy.O7323 = (realnum)((p[3]*a42 + p[4]*a52)*2.72e-12);
00501 CoolHeavy.O7332 = (realnum)((p[3]*a43 + p[4]*a53)*2.71e-12);
00502 CoolHeavy.c3727 = CoolHeavy.O3730 + CoolHeavy.O3726;
00503 CoolHeavy.c7325 = CoolHeavy.O7323 + CoolHeavy.O7332;
00504
00505
00506 CoolAdd("O 2",3727,Cooling );
00507 thermal.dCooldT += CoolingDeriv;
00508
00509
00510
00511 if( (p[3] + p[4]) > SMALLFLOAT )
00512 CoolHeavy.O2_A3_tot = (p[3]*(a41+a42+a43) + p[4]*(a51+a52+a53) ) /
00513 ( (p[3]*(a41+a42+a43) + p[4]*(a51+a52+a53) ) +
00514 ( p[3]*(oii_cs4S32P3+oii_cs2D52P3+oii_cs2D32P3)/go2[3] +
00515 p[4]*(oii_cs4S32P1+oii_cs2D52P1+oii_cs2D32P1)/go2[4]) *
00516 dense.cdsqte );
00517 else
00518 CoolHeavy.O2_A3_tot = 0.;
00519
00520 if( (p[1] + p[2]) > SMALLFLOAT )
00521 CoolHeavy.O2_A2_tot = (p[1]*a21 + p[2]*a31 ) /
00522 ( (p[1]*a21 + p[2]*a31 ) +
00523 ( p[1]*oii_cs4S32D5/go2[1] + p[2]*oii_cs4S32D3/go2[2]) *
00524 dense.cdsqte );
00525 else
00526 CoolHeavy.O2_A2_tot = 0.;
00527
00528
00529
00530
00531 PutCS(oii_cs4S34P,&TauLines[ipT834]);
00532 atom_level2(&TauLines[ipT834]);
00533
00534
00535
00536
00537 double oiii_cs3P25S2,
00538 oiii_cs3P15S2,
00539 oiii_cs3P05S2,
00540 oiii_cs3P1D2,
00541 oiii_cs1D21S0,
00542 oiii_cs3P1S0,
00543 oiii_cs3P03P1,
00544 oiii_cs3P13P2,
00545 oiii_cs3P03P2,
00546 oiii_cs3P3D;
00547
00548
00549
00550
00551
00552 oiii_cs(oiii_cs3P25S2,
00553 oiii_cs3P15S2,
00554 oiii_cs3P05S2,
00555 oiii_cs3P1D2,
00556 oiii_cs1D21S0,
00557 oiii_cs3P1S0,
00558 oiii_cs3P03P1,
00559 oiii_cs3P13P2,
00560 oiii_cs3P03P2,
00561 oiii_cs3P3D);
00562
00563 PutCS(oiii_cs3P25S2,&TauLines[ipT1666]);
00564 PutCS(oiii_cs3P15S2,&TauLines[ipT1661]);
00565 PutCS(oiii_cs3P05S2,&TauDummy);
00566 atom_level3(&TauDummy,&TauLines[ipT1666],&TauLines[ipT1661]);
00567 TauLines[ipT304].Emis->PopOpc = dense.xIonDense[ipOXYGEN][2];
00568 TauLines[ipT304].Lo->Pop = dense.xIonDense[ipOXYGEN][2];
00569 TauLines[ipT304].Hi->Pop = 0.;
00570
00571
00572
00573
00574
00575 aeff = 0.027242 + oxy.d5007r;
00576 a21 = aeff - oxy.d5007r;
00577 a31 = 0.2262;
00578 a32 = 1.685;
00579 oxy.o3ex23 = 32947.;
00580 oxy.o3br32 = (realnum)(a32/(a31 + a32));
00581 oxy.o3enro = (realnum)(4.56e-12/3.98e-12);
00582
00583
00584 oxy.poiii3 = (realnum)(atom_pop3(9.,5.,1.,oiii_cs3P1D2,oiii_cs3P1S0,
00585 oiii_cs1D21S0,a21,a31,a32,28990.,oxy.o3ex23,&oxy.poiii2,
00586 dense.xIonDense[ipOXYGEN][2],oxy.d5007r,0.,0.));
00587 CoolHeavy.c4363 = oxy.poiii3*a32*4.56e-12;
00588 CoolHeavy.c5007 = oxy.poiii2*a21*3.98e-12;
00589 oxy.d5007t = (realnum)(CoolHeavy.c5007*oxy.d5007r/aeff);
00590 thermal.dCooldT += CoolHeavy.c5007*(2.88e4*thermal.tsq1 - thermal.halfte);
00591 thermal.dCooldT += CoolHeavy.c4363*6.20e4*thermal.tsq1;
00592 CoolAdd("O 3",5007,CoolHeavy.c5007);
00593 CoolAdd("O 3",4363,CoolHeavy.c4363);
00594 CoolAdd("O 3",2331,CoolHeavy.c4363*0.236);
00595
00596 if( MAX2(oxy.poiii2,oxy.poiii3) > 0.f && dense.xIonDense[ipOXYGEN][2]>SMALLFLOAT)
00597 {
00598 oxy.poiii3 /= dense.xIonDense[ipOXYGEN][2];
00599 oxy.poiii2 /= dense.xIonDense[ipOXYGEN][2];
00600 }
00601
00602
00603 PutCS(oiii_cs3P03P1,&TauLines[ipTO88]);
00604 PutCS(oiii_cs3P13P2,&TauLines[ipT52]);
00605 PutCS(oiii_cs3P03P2,&TauDummy);
00606 atom_level3(&TauLines[ipTO88],&TauLines[ipT52],&TauDummy);
00607
00608
00609 PutCS(oiii_cs3P3D,&TauLines[ipT835]);
00610 atom_level2(&TauLines[ipT835]);
00611
00612
00613
00614
00615 double oiv_cs2P2D,oiv_cs2P12P3;
00616
00617
00618
00619
00620
00621
00622 oiv_cs(oiv_cs2P2D,oiv_cs2P12P3);
00623
00624
00625 PutCS(oiv_cs2P2D,&TauLines[ipT789]);
00626 atom_level2(&TauLines[ipT789]);
00627
00628
00629
00630
00631
00632
00633
00634
00635 PutCS(oiv_cs2P12P3,&TauLines[ipT26]);
00636 static vector< pair<transition*,double> > O4Pump;
00637 O4Pump.reserve(48);
00638
00639 if( O4Pump.empty() )
00640 {
00641
00642 pair<transition*,double> pp( &TauLines[ipT789], 1./6. );
00643 O4Pump.push_back( pp );
00644
00645 for( i=0; i < nWindLine; ++i )
00646 {
00647
00648 if( TauLine2[i].Hi->nelem == 8 && TauLine2[i].Hi->IonStg == 4 )
00649 {
00650 # if 0
00651 DumpLine( &TauLine2[i] );
00652 # endif
00653 double branch_ratio;
00654
00655
00656 if( fp_equal( TauLine2[i].Hi->g, realnum(2.) ) )
00657 branch_ratio = 2./3.;
00658 else if( fp_equal( TauLine2[i].Hi->g, realnum(6.) ) )
00659 branch_ratio = 1./2.;
00660 else if( fp_equal( TauLine2[i].Hi->g, realnum(10.) ) )
00661 branch_ratio = 1./6.;
00662 else
00663 TotalInsanity();
00664 pair<transition*,double> pp2( &TauLine2[i], branch_ratio );
00665 O4Pump.push_back( pp2 );
00666 }
00667 }
00668 }
00669
00670
00671 double pump_rate = 0.;
00672 vector< pair<transition*,double> >::const_iterator o4p;
00673 for( o4p=O4Pump.begin(); o4p != O4Pump.end(); ++o4p )
00674 {
00675 const transition* t = o4p->first;
00676 double branch_ratio = o4p->second;
00677 pump_rate += t->Emis->pump*branch_ratio;
00678 # if 0
00679 dprintf( ioQQQ, "O IV %.3e %.3e\n",
00680 t->WLAng , t->Emis->pump*branch_ratio );
00681 # endif
00682 }
00683
00684
00685
00686
00687
00688 AtomSeqBoron(&TauLines[ipT26],
00689 &TauLines[ipO4_1400],
00690 &TauLines[ipO4_1397],
00691 &TauLines[ipO4_1407],
00692 &TauLines[ipO4_1405],
00693 &TauLines[ipO4_1401],
00694 0.1367 , 0.1560 , 1.1564 , 0.0457 , pump_rate,"O 4");
00695
00696
00697
00698
00699 double ov_cs1S01P1,ov_cs1S03P;
00700
00701
00702
00703
00704
00705 ov_cs(ov_cs1S01P1,ov_cs1S03P);
00706 PutCS(ov_cs1S01P1,&TauLines[ipT630]);
00707 atom_level2(&TauLines[ipT630]);
00708
00709 PutCS(ov_cs1S03P,&TauLines[ipT1214]);
00710
00711
00712
00713
00714
00715 AtomSeqBeryllium(.87,0.602,2.86,&TauLines[ipT1214],.02198);
00716 embesq.em1218 = (realnum)(atoms.PopLevels[3]*0.0216*1.64e-11);
00717
00718
00719
00720
00721 ligbar(8,&TauLines[ipT1032],&TauLines[ipT150],&cs2s2p,&cs2s3p);
00722 PutCS(cs2s2p,&TauLines[ipT1032]);
00723 PutCS(cs2s2p*0.5,&TauLines[ipT1037]);
00724
00725 PutCS(1.0,&TauDummy);
00726
00727 atom_level3(&TauLines[ipT1037],&TauDummy,&TauLines[ipT1032]);
00728
00729 PutCS(cs2s3p,&TauLines[ipT150]);
00730 atom_level2(&TauLines[ipT150]);
00731 return;
00732 }
00733
00734
00735
00736
00737
00738
00739
00740
00741
00742
00743
00744
00745
00746
00747 void oi_cs(double& oi_cs3P23P1,double& oi_cs3P23P0,double& oi_cs3P13P0,
00748 double& oi_cs3P1D2,double& oi_cs1D21S0,double& oi_cs3P1S0)
00749 {
00750 double a,
00751 b,
00752 c,
00753 d;
00754
00755 DEBUG_ENTRY( "oi_cs()" );
00756
00757
00758
00759 if( phycon.te <= 3.0e3 )
00760 {
00761 oi_cs3P23P1 = 1.49e-4*phycon.sqrte/phycon.te02/phycon.te02;
00762 }
00763 else if( phycon.te <= 1.0e4 )
00764 {
00765 a = -5.5634127e-04;
00766 b = 8.3458102e-08;
00767 c = 2.3068232e-04;
00768 oi_cs3P23P1 = 0.228f*(a + b*phycon.te32 + c*phycon.sqrte);
00769 }
00770 else
00771 {
00772 oi_cs3P23P1 = 0.228*MIN2(0.222,5.563e-06*phycon.te*phycon.te05*
00773 phycon.te02);
00774 }
00775
00776
00777 if( phycon.te <= 3.0e3 )
00778 {
00779 oi_cs3P23P0 = 4.98e-5*phycon.sqrte;
00780 }
00781 else if( phycon.te <= 1.0e4 )
00782 {
00783 a = -3.7178028e-04;
00784 b = 2.0569267e-08;
00785 c = 1.1898539e-04;
00786 oi_cs3P23P0 = 0.288*(a + b*phycon.te32 + c*phycon.sqrte);
00787 }
00788 else
00789 {
00790 oi_cs3P23P0 = 0.288*MIN2(0.0589,1.015e-05*phycon.te/phycon.te10/
00791 phycon.te02/phycon.te005);
00792 }
00793
00794
00795 if( phycon.te <= 3.0e3 )
00796 {
00797 oi_cs3P13P0 = 1.83e-9*phycon.te*phycon.te30*phycon.te05;
00798 }
00799 else if( phycon.te <= 1.0e4 )
00800 {
00801 a = -5.9364373e-04;
00802 b = 0.02946867;
00803 c = 10768.675;
00804 d = 3871.9826;
00805 oi_cs3P13P0 = 0.0269*(a + b*exp(-0.5*POW2((phycon.te-c)/d)));
00806 }
00807 else
00808 {
00809 oi_cs3P13P0= 0.0269*MIN2(0.074,7.794e-08*phycon.te32/phycon.te10/
00810 phycon.te01);
00811 }
00812
00813
00814
00815
00816
00817
00818
00819
00820
00821
00822
00823 if(phycon.te < 8E3)
00824 oi_cs3P1D2 = (4E-08)*(phycon.te*phycon.te70*phycon.te05);
00825 else if(phycon.te < 2E4)
00826 oi_cs3P1D2 = (4.630155E-05)*((phycon.te/phycon.te04)*phycon.te005*phycon.te0001);
00827 else
00828 oi_cs3P1D2 = (1.5394E-03)*(phycon.sqrte*phycon.te10*phycon.te01*phycon.te001*phycon.te0003);
00829
00830
00831
00832
00833
00834
00835 double te_scale = phycon.te / 6000.;
00836 double rate_H0 = (1.74*te_scale + 0.6)*1e-12*sexp(0.47*te_scale) / sqrt(te_scale ) *
00837 dense.xIonDense[ipHYDROGEN][0];
00838 oi_cs3P1D2 += ConvRate2CS( 5.f , (realnum)rate_H0 );
00839
00840 if(phycon.te < 5E3)
00841 oi_cs1D21S0 = (7E-08)*(phycon.te*phycon.sqrte*phycon.te10*phycon.te007*phycon.te0001);
00842 else if(phycon.te<2E4)
00843 oi_cs1D21S0 = (1.98479e-04)*((phycon.te70/phycon.te03)*phycon.te003*phycon.te0007);
00844 else
00845 oi_cs1D21S0 = (9.31E-04)*(phycon.sqrte*phycon.te01*phycon.te007*phycon.te0005*phycon.te0001);
00846
00847
00848 if(phycon.te < 2E4)
00849 oi_cs3P1S0 = (2E-05)*(phycon.sqrte*phycon.te30*phycon.te05*phycon.te01*(phycon.te004/phycon.te0002));
00850 else
00851 oi_cs3P1S0 = (1.054E-03)*(phycon.sqrte/phycon.te04)*phycon.te003*phycon.te0005;
00852
00853
00854 return;
00855 }
00856
00857
00858
00859
00860
00861
00862
00863
00864
00865 void oi_othercs(double& csh01,double& cshe01,double& csh201,double& csh12,double& cshe12,
00866 double& csh212,double& csh02,double& cshe02,double& csh202,double& csh2o01,
00867 double& csh2o02,double& csh2o12,double& csh2p01,double& csh2p02,double& csh2p12,
00868 double& csp01,double& csp02,double& csp12)
00869 {
00870 DEBUG_ENTRY( "oi_othercs()" );
00871
00872
00873
00874
00875
00876
00877
00878
00879
00880
00881
00882
00883
00884
00885
00886
00887
00888
00889
00890
00891
00892
00893
00894
00895
00896
00897
00898
00899
00900
00901
00902 csh12 = MAX2(2.00e-11*phycon.te30*phycon.te05*phycon.te02,
00903 7.67e-12*phycon.sqrte*phycon.te03);
00904
00905
00906 csh01 = MIN2(3.12e-12*phycon.te70*phycon.te02*phycon.te02,
00907 7.51e-12*phycon.sqrte*phycon.te05*phycon.te03);
00908
00909
00910 csh02 = MIN2(6.96e-13*phycon.te/phycon.te10/phycon.te02,
00911 1.49e-12*phycon.te70*phycon.te05);
00912
00913
00914
00915
00916
00917
00918 cshe12 = MIN2(8.09e-16*phycon.te32*phycon.te10*phycon.te03,
00919 9.72e-15*phycon.te*phycon.te20);
00920
00921 cshe01 = 1.57e-12*phycon.te70/phycon.te01;
00922
00923 cshe02 = MIN2(1.80e-12*phycon.te70*phycon.te05,
00924 4.45e-12*phycon.te70/phycon.te10);
00925
00926 if( phycon.te<=1.5e3 )
00927 {
00928
00929 double ortho_frac = h2.ortho_density/SDIV(hmi.H2_total);
00930
00931
00932
00933
00934 csh2o12 = 2.21e-14*phycon.te*phycon.te10/phycon.te01;
00935 csh2p12 = 9.45e-15*phycon.te*phycon.te20/phycon.te01;
00936 csh212 = ortho_frac*csh2o12 + (1.-ortho_frac)*csh2p12;
00937
00938 csh2o01 = 2.37e-11*phycon.te30*phycon.te10/phycon.te02;
00939 csh2p01 = 3.40e-11*phycon.te30*phycon.te02;
00940 csh201 = ortho_frac*csh2o01 + (1.-ortho_frac)*csh2p01;
00941
00942 csh2o02 = 4.86e-11*phycon.te30*phycon.te02*phycon.te02;
00943 csh2p02 = 6.82e-11*phycon.te30/phycon.te02;
00944 csh202 = ortho_frac*csh2o02 + (1.-ortho_frac)*csh2p02;
00945 }
00946 else
00947 {
00948 csh212 = 0.;
00949 csh201 = 0.;
00950 csh202 = 0.;
00951 }
00952
00953
00954
00955
00956
00957
00958 if( phycon.te<=1000. )
00959 {
00960 csp01 = MAX2(2.22e-5*phycon.te/phycon.te10,
00961
00962
00963 2.69e-6*phycon.te*phycon.te30)*dense.xIonDense[ipHYDROGEN][1]/dense.EdenHCorr;
00964
00965 csp02 = MIN2(7.07e-8*phycon.te32*phycon.te10,2.46e-7*
00966
00967
00968 phycon.te32/phycon.te10)*dense.xIonDense[ipHYDROGEN][1]/dense.EdenHCorr;
00969 }
00970 else
00971 {
00972 csp01 = MIN2(2.69e-6*phycon.te*phycon.te30,9.21e-5*phycon.te/phycon.te10/
00973
00974
00975 phycon.te03)*dense.xIonDense[ipHYDROGEN][1]/dense.EdenHCorr;
00976
00977 csp02 = MIN2(2.46e-7*phycon.te32/phycon.te10,5.21e-5*phycon.te/
00978
00979
00980 phycon.te20)*dense.xIonDense[ipHYDROGEN][1]/dense.EdenHCorr;
00981 }
00982
00983 csp12 = MIN2(2.35e-6*phycon.te*phycon.te05*phycon.te01,3.98e-5*
00984
00985
00986
00987 phycon.te70/phycon.te01)*dense.xIonDense[ipHYDROGEN][1]/dense.EdenHCorr;
00988
00989 return;
00990 }
00991
00992
00993
00994
00995
00996
00997
00998
00999
01000 void oii_cs(double& oii_cs4S32D5,
01001 double& oii_cs4S32D3,
01002 double& oii_cs2D52D3,
01003 double& oii_cs4S32P3,
01004 double& oii_cs2D52P3,
01005 double& oii_cs2D32P3,
01006 double& oii_cs4S32P1,
01007 double& oii_cs2D52P1,
01008 double& oii_cs2D32P1,
01009 double& oii_cs2P32P1,
01010 double& oii_cs4S34P)
01011 {
01012 DEBUG_ENTRY( "oii_cs()" );
01013
01014
01015 oii_cs4S32D5 = (realnum)(0.7776*(phycon.te007*phycon.te0005*phycon.te0001));
01016 oii_cs4S32D3 = (realnum)(0.5643/phycon.te002);
01017 oii_cs2D52D3 = (realnum)(2.2277/(phycon.te07/(phycon.te003*phycon.te0001)));
01018 oii_cs4S32P3 = (realnum)(0.2004*phycon.te02*(phycon.te007/phycon.te0004));
01019 oii_cs2D52P3 = (realnum)(0.6211*phycon.te03*phycon.te004*phycon.te0002);
01020 oii_cs2D32P3 = (realnum)(0.3159*phycon.te04*(phycon.te004/phycon.te0004));
01021 oii_cs4S32P1 = (realnum)(0.1112*(phycon.te02/(phycon.te001*phycon.te0004)));
01022 oii_cs2D52P1 = (realnum)(0.2337*phycon.te04*phycon.te0004);
01023 oii_cs2D32P1 = (realnum)(0.2226*phycon.te04*phycon.te003*phycon.te0001);
01024 oii_cs2P32P1 = (realnum)(0.1943*phycon.te04*(phycon.te002/phycon.te0004));
01025
01026
01027
01028
01029 oii_cs4S34P = 0.355*phycon.te10*phycon.te10*phycon.te003*phycon.te001*phycon.te001;
01030
01031 return;
01032 }
01033
01034
01035
01036
01037
01038
01039 void oiii_cs(double& oiii_cs3P25S2,
01040 double& oiii_cs3P15S2,
01041 double& oiii_cs3P05S2,
01042 double& oiii_cs3P1D2,
01043 double& oiii_cs1D21S0,
01044 double& oiii_cs3P1S0,
01045 double& oiii_cs3P03P1,
01046 double& oiii_cs3P13P2,
01047 double &oiii_cs3P03P2,
01048 double& oiii_cs3P3D)
01049 {
01050 DEBUG_ENTRY( "oiii_cs()" );
01051
01052
01053
01054
01055
01056
01057
01058
01059
01060
01061 if(phycon.te < 30000)
01062 oiii_cs3P25S2 = (realnum)(0.2519*(phycon.te07*phycon.te02*phycon.te007*phycon.te001*phycon.te0002));
01063 else
01064 oiii_cs3P25S2 = (realnum)(6.166388*(1/(phycon.te20*phycon.te01*phycon.te002)));
01065
01066
01067
01068
01069
01070
01071 if(phycon.te < 30000)
01072 oiii_cs3P15S2 = (realnum)((0.1511)*(phycon.te07*phycon.te02*phycon.te007*phycon.te001*phycon.te0002));
01073 else
01074 oiii_cs3P15S2 = (realnum)((3.668474)*(1/(phycon.te20*phycon.te01*phycon.te001*phycon.te0002)));
01075
01076
01077
01078
01079
01080
01081 if(phycon.te < 30000)
01082 oiii_cs3P05S2 = (realnum)(0.0504*((phycon.te10/phycon.te01)*phycon.te005*phycon.te003*phycon.te0002));
01083 else
01084 oiii_cs3P05S2 = (realnum)(1.223633/(phycon.te20*phycon.te01*phycon.te001*phycon.te0002));
01085
01086
01087
01088
01089
01090
01091
01092
01093
01094
01095
01096
01097
01098
01099
01100
01101
01102
01103
01104
01105
01106
01107
01108
01109
01110
01111
01112
01113 if(phycon.te < 3E4)
01114 {
01115 oiii_cs3P1D2 = (realnum)(0.9062*(phycon.te10/phycon.te002));
01116 oiii_cs1D21S0 = (realnum)(0.0995*phycon.te10*phycon.te07*(phycon.te007/phycon.te0002));
01117 oiii_cs3P1S0 = (realnum)(0.1237*phycon.te07*phycon.te02*phycon.te005*phycon.te0005);
01118
01119 }
01120 else if(phycon.te < 6E4)
01121 {
01122 oiii_cs3P1D2 = (realnum)(1.710073*(phycon.te04/phycon.te004)*phycon.te0004);
01123 oiii_cs3P1S0 = (realnum)(0.1963109*phycon.te05*phycon.te0007);
01124 oiii_cs1D21S0 = (realnum)(0.781266/(phycon.te02*phycon.te003*phycon.te0001));
01125 }
01126 else
01127 {
01128 oiii_cs3P1D2 = (realnum)(6.452638/((phycon.te10/phycon.te02)*phycon.te004*phycon.te0003));
01129 oiii_cs3P1S0 = (realnum)(0.841578/(phycon.te07*phycon.te01*(phycon.te002/phycon.te0004)));
01130 oiii_cs1D21S0 = (realnum)(0.781266/(phycon.te02*phycon.te003*phycon.te0001));
01131 }
01132
01133
01134
01135
01136
01137
01138
01139
01140 if(phycon.te < 3E4)
01141 oiii_cs3P03P1 = (realnum)(0.3993*(phycon.te03/phycon.te001)*phycon.te0001);
01142 else if(phycon.te < 1E5)
01143 oiii_cs3P03P1 = (realnum)(0.245712*phycon.te07*phycon.te005*phycon.te001*phycon.te0002);
01144 else
01145 oiii_cs3P03P1 = (realnum)(1.310467/((phycon.te07/phycon.te001)*phycon.te0002));
01146
01147
01148 if(phycon.te < 3E4)
01149 oiii_cs3P13P2 = (realnum)(0.7812*(phycon.te05/phycon.te0001));
01150 else if(phycon.te < 1.2E5)
01151 oiii_cs3P13P2 = (realnum)(0.516157*phycon.te07*phycon.te02*phycon.te0001);
01152 else
01153 oiii_cs3P13P2 = (realnum)(4.038402/(phycon.te05*phycon.te03*phycon.te005*phycon.te0007*phycon.te0001));
01154
01155
01156 if(phycon.te < 3E4)
01157 oiii_cs3P03P2 = (realnum)(0.1337*phycon.te07*phycon.te002*phycon.te0002);
01158 else if(phycon.te < 1.2E5)
01159 oiii_cs3P03P2 = (realnum)(0.086446*phycon.te10*phycon.te01*phycon.te004*phycon.te0005);
01160 else
01161 oiii_cs3P03P2 = (realnum)(0.82031/(phycon.te07*phycon.te007*phycon.te0007*phycon.te0002));
01162
01163
01164
01165
01166
01167
01168
01169 if(phycon.te < 3E4)
01170 oiii_cs3P3D = (realnum)(4.74*phycon.te04*phycon.te0002);
01171 else
01172 oiii_cs3P3D = (realnum)(0.533*phycon.te20*phycon.te05*phycon.te002*phycon.te0001);
01173 return;
01174 }
01175
01176
01177
01178
01179
01180
01181
01182 void oiv_cs(double& oiv_cs2P2D,double& oiv_cs2P12P3)
01183 {
01184 double Te_bounded,Te_bounded_log;
01185
01186 DEBUG_ENTRY( "oiv_cs()" );
01187
01188 Te_bounded = MAX2(phycon.te,4500.);
01189 Te_bounded = MIN2(Te_bounded,450000.);
01190 Te_bounded_log = log(Te_bounded);
01191
01192
01193 oiv_cs2P2D = -3.0102462 + 109.22973/Te_bounded_log - 18666.357/Te_bounded;
01194
01195
01196
01197
01198
01199
01200
01201
01202 oiv_cs2P12P3 = (realnum)(exp(3.265723 - 0.00014686984*phycon.alnte*phycon.sqrte
01203 -22.035066/phycon.alnte));
01204
01205
01206 oiv_cs2P12P3 = MAX2( oiv_cs2P12P3 , 3.25e-2);
01207 return;
01208 }
01209
01210
01211
01212
01213
01214 void ov_cs(double& ov_cs1S01P1,double& ov_cs1S03P)
01215 {
01216 DEBUG_ENTRY( "ov_cs()" );
01217
01218
01219
01220 ov_cs1S01P1 = MIN2(4.0,1.317*phycon.te10/phycon.te03);
01221
01222
01223
01224 if( phycon.te <= 3.16e4 )
01225 {
01226 ov_cs1S03P = 3.224/(phycon.te10*phycon.te03*phycon.te03*phycon.te003);
01227 }
01228 else
01229 {
01230 ov_cs1S03P = 10.549/(phycon.te10*phycon.te10*phycon.te10/phycon.te03);
01231 }
01232
01233 return;
01234 }
01235