Changeset 1853 for trunk/source/hydro_vs_rates.cpp

Timestamp:

03/17/08 07:39:48 (10 months ago)

Author:

rporter

Message:

A set of changes to iso sequence collisions. The main effects are greater collisional excitation from ground states at kT much less then deltaE and a smoothing of excitation onto "collapsed" levels.

Files:

• trunk/source/hydro_vs_rates.cpp (modified) (9 diffs)

trunk/source/hydro_vs_rates.cpp

@@ -14 +14 @@
 #include "iso.h"
 #include "hydro_vs_rates.h"
+#include "lines_service.h"
 
 STATIC double hydro_vs_coll_str( double energy );
@@ -19 +20 @@
 
 static long global_ipISO, global_nelem, global_ipHi, global_ipLo, global_Collider;
-static double global_temp;
+static double global_temp, global_Aul;
 /* These are masses relative to the proton mass of the electron, proton, and alpha particle. */
 static double ColliderMass[4] = {ELECTRON_MASS/PROTON_MASS, 1.0, 4.0, 4.0};
@@ -33 +34 @@
 /* VS80 stands for Vriens and Smeets 1980 */
 /* This routine calculates thermally-averaged collision strengths. */
-double CS_VS80( long int ipISO, long int nelem, long int ipHi, long int ipLo, double temp, long int Collider )
+double CS_VS80( long int ipISO, long int nelem, long int ipHi, long int ipLo, double Aul, double temp, long int Collider )
 {
         double coll_str;
@@ -43 +44 @@
         global_temp = temp;
         global_Collider = Collider;
+        global_Aul = Aul;
 
         /* evaluate collision strength, two options,
@@ -84 +86 @@
         double gamma, p, n;
         double ryd, s;
-        double cross_section, coll_str, stat_weight;
+        double cross_section, coll_str, gLo, gHi, abs_osc_str, Aul;
         long ipISO, nelem, ipHi, ipLo, Collider;
 
@@ -94 +96 @@
         ipLo = global_ipLo;
         Collider = global_Collider;
-
-        stat_weight = StatesElem[ipISO][nelem][ipLo].g;
+        Aul = global_Aul;
+
+        gLo = StatesElem[ipISO][nelem][ipLo].g;
+        gHi = StatesElem[ipISO][nelem][ipHi].g;
 
         /* This comes from equations 14,15, and 16 of Vriens and Smeets. */
@@ -114 +118 @@
 
         /* This is an absorption oscillator strength. */
-        Apn = 2.*ryd/Epn*HydroOscilStr(p,n);
+        abs_osc_str = GetGF( Aul, Epn*RYD_INF/EVRYD, gHi)/gLo;
+        Apn = 2.*ryd/Epn*abs_osc_str;
 
         bp = 1.4*log(p)/p - .7/p - .51/p/p + 1.16/p/p/p - 0.55/p/p/p/p;
@@ -140 +145 @@
 
         /* convert to collision strength */
-        coll_str = cross_section * stat_weight * (energy/EVRYD);
+        coll_str = cross_section * gLo * (energy/EVRYD);
 
         if( nelem==ipCARBON )
@@ -300 +305 @@
         return( HydColIon_v );
 }
+
+/*hydro_vs_deexcit compute collisional deexcitation rates for hydrogen atom, 
+ * from Vriens and Smeets 1980 */
+double hydro_vs_deexcit( long ipISO, long nelem, long ipHi, long ipLo, double Aul )
+{
+        double Anp, bn, Bnp, delta_np;
+        double Eni, Enp;
+        double Gamma_np, p, n, g_p, g_n;
+        double ryd, s, kT_eV, rate, col_str, abs_osc_str;
+        long Collider;
+
+        DEBUG_ENTRY( "hydro_vs_coll_str()" );
+
+        Collider = ipELECTRON;
+
+        kT_eV = EVRYD * phycon.te / TE1RYD;
+
+        /* This comes from equations 24 of Vriens and Smeets. */
+        /* >>refer he-like cs Vriens, L. \& Smeets, A. H. M. Phys. Rev. A, 22, 940 */ 
+        /* Computes the Vriens and Smeets
+         * rate coeff. for collisional dexcitation between any two levels of H.
+         * valid for all nhigh, nlow
+         * at end converts this excitation rate to collision strength   */
+
+        n = (double)StatesElem[ipISO][nelem][ipLo].n;
+        p = (double)StatesElem[ipISO][nelem][ipHi].n;
+
+        ASSERT( n!=p );
+
+        g_n = StatesElem[ipISO][nelem][ipLo].g;
+        g_p = StatesElem[ipISO][nelem][ipHi].g;
+
+        ryd = EVRYD;
+        s = fabs(n-p);
+
+        Enp = EVRYD*(iso.xIsoLevNIonRyd[ipISO][nelem][ipLo] - iso.xIsoLevNIonRyd[ipISO][nelem][ipHi]);
+        Eni = EVRYD*iso.xIsoLevNIonRyd[ipISO][nelem][ipHi];
+
+        ASSERT( Enp > 0. );
+
+        /* This is an absorption oscillator strength. */
+        abs_osc_str = GetGF( Aul, Enp*RYD_INF/EVRYD, g_p)/g_n;
+        Anp = 2.*ryd/Enp*abs_osc_str;
+
+        bn = 1.4*log(n)/n - .7/n - .51/n/n + 1.16/n/n/n - 0.55/n/n/n/n;
+
+        Bnp = 4.*ryd*ryd/p/p/p*(1./Enp/Enp + 4./3.*Eni/POW3(Enp) + bn*Eni*Eni/powi(Enp,4));
+
+        delta_np = exp(-Bnp/Anp) + 0.1*Enp/ryd;
+
+        Gamma_np = ryd*log(1. + n*n*n*kT_eV/ryd) * (3. + 11.*s*s/n/n) /
+                ( 6. + 1.6*p*s + 0.3/s/s + 0.8*(pow(p,1.5))/(pow(s,0.5))*fabs(s-0.6) );
+
+        if( 0.3*kT_eV/ryd+delta_np <= 0 )
+        {
+                rate = 0.;
+        }
+        else
+        {
+                rate = 1.6E-7 * pow(kT_eV, 0.5) * g_n / g_p / ( kT_eV + Gamma_np ) *
+                        ( Anp * log(0.3*kT_eV/ryd + delta_np) + Bnp );
+        }
+
+        col_str = rate / COLL_CONST * phycon.sqrte * StatesElem[ipISO][nelem][ipHi].g;
+
+        return( col_str );
+}