00001
00002
00003
00004 #include "cddefines.h"
00005 #include "taulines.h"
00006 #include "iso.h"
00007 #include "dense.h"
00008 #include "hydrogenic.h"
00009 #include "phycon.h"
00010 #include "rt.h"
00011
00012 void RT_stark(void)
00013 {
00014 long int ipLo,
00015 ipHi,
00016 nelem,
00017 ipISO;
00018
00019 double aa , ah,
00020 stark,
00021 strkla;
00022
00023 DEBUG_ENTRY( "RT_stark()" );
00024
00025
00026 static long int nZoneEval=-1;
00027 if( nzone==nZoneEval )
00028 return;
00029 nZoneEval = nzone;
00030
00031 for( ipISO=ipH_LIKE; ipISO<NISO; ++ipISO )
00032 {
00033
00034 for( nelem=ipISO; nelem<LIMELM; ++nelem )
00035 {
00036 if( nelem >= 2 && !dense.lgElmtOn[nelem] )
00037 continue;
00038
00039 if( !rt.lgStarkON || dense.eden < 1e8 )
00040 {
00041 for( ipHi=0; ipHi < iso.numLevels_max[ipISO][nelem]; ipHi++ )
00042 {
00043 for( ipLo=0; ipLo < iso.numLevels_max[ipISO][nelem]; ipLo++ )
00044 {
00045 iso.pestrk[ipISO][nelem][ipHi][ipLo] = 0.;
00046 iso.pestrk[ipISO][nelem][ipLo][ipHi] = 0.;
00047 }
00048 }
00049 continue;
00050 }
00051
00052
00053
00054
00055
00056
00057 ah = 6.9e-6*1000./1e12/(phycon.sqrte*phycon.te10*phycon.te10*
00058 phycon.te03*phycon.te01*phycon.te01)*dense.eden;
00059
00060
00061 ah *= pow( (double)(nelem+1), -4.5 );
00062
00063
00064
00065 stark = 0.264*pow(ah,0.4);
00066
00067
00068
00069 strkla = 0.538*ah*4.*9.875*(phycon.sqrte/phycon.te10/phycon.te03);
00070
00071
00072
00073
00074
00075
00076 aa = (realnum)SDIV(Transitions[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Emis->TauIn);
00077 aa = pow( aa ,-0.75);
00078 iso.pestrk[ipISO][nelem][0][iso.nLyaLevel[ipISO]] = strkla/2.*MAX2(1.,aa);
00079
00084 iso.pestrk[ipISO][nelem][0][iso.nLyaLevel[ipISO]] = MIN2(.01,iso.pestrk[ipISO][nelem][0][iso.nLyaLevel[ipISO]]);
00085 iso.pestrk[ipISO][nelem][0][iso.nLyaLevel[ipISO]] = 0.;
00086 iso.pestrk[ipISO][nelem][iso.nLyaLevel[ipISO]][0] =
00087 iso.pestrk[ipISO][nelem][0][iso.nLyaLevel[ipISO]]*Transitions[ipISO][nelem][iso.nLyaLevel[ipISO]][0].Emis->Aul;
00088
00089
00090
00091 for( ipHi=3; ipHi < iso.numLevels_local[ipISO][nelem]; ipHi++ )
00092 {
00093 if( Transitions[ipISO][nelem][ipHi][ipH1s].ipCont <= 0 )
00094 continue;
00095
00096 iso.pestrk[ipISO][nelem][0][ipHi] = stark*iso.strkar[ipISO][nelem][0][ipHi]/2.*pow(MAX2(1.,
00097 Transitions[ipISO][nelem][ipHi][ipH1s].Emis->TauIn),-0.75);
00098
00099 iso.pestrk[ipISO][nelem][0][ipHi] = MIN2(.01,iso.pestrk[ipISO][nelem][0][ipHi]);
00100 iso.pestrk[ipISO][nelem][ipHi][0] = Transitions[ipISO][nelem][ipHi][ipH1s].Emis->Aul*
00101 iso.pestrk[ipISO][nelem][0][ipHi];
00102 }
00103
00104
00105 for( ipHi=iso.numLevels_local[ipISO][nelem]; ipHi < iso.numLevels_max[ipISO][nelem]; ipHi++ )
00106 {
00107 iso.pestrk[ipISO][nelem][0][ipHi] = 0.;
00108 iso.pestrk[ipISO][nelem][ipHi][0] = 0.;
00109 }
00110
00111
00112 for( ipLo=ipH2s; ipLo < (iso.numLevels_local[ipISO][nelem] - 1); ipLo++ )
00113 {
00114 for( ipHi=ipLo + 1; ipHi < iso.numLevels_local[ipISO][nelem]; ipHi++ )
00115 {
00116 if( Transitions[ipISO][nelem][ipHi][ipLo].ipCont <= 0 )
00117 continue;
00118
00119 aa = stark*iso.strkar[ipISO][nelem][ipLo][ipHi]*
00120 pow(MAX2(1.,Transitions[ipISO][nelem][ipHi][ipLo].Emis->TauIn),-0.75);
00121 iso.pestrk[ipISO][nelem][ipLo][ipHi] =
00122 (realnum)MIN2(.01,aa);
00123
00124 iso.pestrk[ipISO][nelem][ipHi][ipLo] = Transitions[ipISO][nelem][ipHi][ipLo].Emis->Aul*
00125 iso.pestrk[ipISO][nelem][ipLo][ipHi];
00126 }
00127 }
00128
00129
00130 for( ipLo=(iso.numLevels_local[ipISO][nelem] - 1); ipLo<(iso.numLevels_max[ipISO][nelem] - 1); ipLo++ )
00131 {
00132 for( ipHi=ipLo + 1; ipHi < iso.numLevels_max[ipISO][nelem]; ipHi++ )
00133 {
00134 iso.pestrk[ipISO][nelem][ipLo][ipHi] = 0.;
00135 iso.pestrk[ipISO][nelem][ipHi][ipLo] = 0.;
00136 }
00137 }
00138 }
00139 }
00140
00141 return;
00142 }