cloudy trunk
Loading...
Searching...
No Matches
hydrocollid.cpp
Go to the documentation of this file.
1/* This file is part of Cloudy and is copyright (C)1978-2013 by Gary J. Ferland and
2 * others. For conditions of distribution and use see copyright notice in license.txt */
3/*HCSAR_interp interpolate on collision strengths */
4/*C6cs123 line collision rates for lower levels of hydrogenic carbon, n=1,2,3 */
5/*Ca20cs123 line collision rates for lower levels of hydrogenic calcium, n=1,2,3 */
6/*Hydcs123 Hydrogenic de-excitation collision rates n=1,2,3 */
7/*H1cs123 hydrogen collision data levels involving 1s,2s,2p,3. */
8/*Ne10cs123 line collision rates for lower levels of hydrogenic neon, n=1,2,3 */
9/*He2cs123 line collision strengths for lower levels of helium ion, n=1,2,3, by K Korista */
10/*Fe26cs123 line collision rates for lower levels of hydrogenic iron, n=1,2,3 */
11#include "cddefines.h"
12#include "atmdat.h"
13#include "dense.h"
14#include "helike_cs.h"
15#include "hydro_vs_rates.h"
16#include "iso.h"
17#include "opacity.h"
18#include "phycon.h"
19#include "physconst.h"
20#include "taulines.h"
21
22STATIC double Fe26cs123(long int i, long int j);
23STATIC double He2cs123(long int i, long int j);
24STATIC double Hydcs123(long int ilow, long int ihigh, long int iz, long int chType);
25STATIC double C6cs123(long int i, long int j);
26STATIC double Ca20cs123(long int i, long int j);
27STATIC double Ne10cs123(long int i, long int j);
28
29STATIC realnum HCSAR_interp( int ipLo , int ipHi );
30STATIC double CS_ThermAve_PR78(long ipISO, long nelem, long nHi, long nLo, double deltaE, double temp );
31STATIC double Therm_ave_coll_str_int_PR78( double EOverKT );
32STATIC double CS_PercivalRichards78( double Ebar );
33
35static double kTRyd, global_deltaE;
36
37static const realnum HCSTE[NHCSTE] = {5802.f,11604.f,34812.f,58020.f,116040.f,174060.f,232080.f,290100.f};
38
39/*HCSAR_interp interpolate on collision strengths */
40STATIC realnum HCSAR_interp( int ipLo , int ipHi )
41{
42
43 static int ip=1;
44 realnum cs;
45
46 DEBUG_ENTRY( "HCSAR_interp()" );
47
48 if( ipLo==1 && ipHi==2 )
49 {
50 fprintf(ioQQQ,"HCSAR_interp was called for the 2s-2p transition, which it cannot do\n");
52 }
53 if( phycon.te <= HCSTE[0] )
54 {
55 cs = t_ADfA::Inst().h_coll_str( ipLo , ipHi , 0 );
56 }
57 else if( phycon.te >= HCSTE[NHCSTE-1] )
58 {
59 cs = t_ADfA::Inst().h_coll_str( ipLo , ipHi , NHCSTE-1 );
60 }
61 else
62 {
63 /* the ip index is most likely correct since it points to the last temperature */
64 if( (HCSTE[ip-1] >= phycon.te ) || ( phycon.te > HCSTE[ip]) )
65 {
66 /* we must find the temperature in the array */
67 for( ip=1; ip<NHCSTE; ++ip )
68 {
69 if( (HCSTE[ip-1] < phycon.te ) && ( phycon.te <= HCSTE[ip]) )
70 break;
71 }
72 }
73 /* we now have the index */
74 cs = t_ADfA::Inst().h_coll_str( ipLo , ipHi , ip-1 ) +
75 (t_ADfA::Inst().h_coll_str( ipLo , ipHi , ip ) - t_ADfA::Inst().h_coll_str( ipLo , ipHi , ip-1 ) ) / (HCSTE[ip]-HCSTE[ip-1] ) *
76 ((realnum)phycon.te - HCSTE[ip-1] );
77 if( cs <= 0.)
78 {
79 fprintf(ioQQQ," insane cs returned by HCSAR_interp, values are\n");
80 fprintf(ioQQQ,"%.3f %.3f \n", t_ADfA::Inst().h_coll_str( ipLo , ipHi , ip-1 ),t_ADfA::Inst().h_coll_str( ipLo , ipHi , ip ) );
81 }
82 }
83 return(cs);
84}
85
86/*Hydcs123 Hydrogenic de-excitation collision strengths between levels n=1,2,3,
87 * for any charge. routine only called by HydroCSInterp to fill in hydroline arrays
88 * with collision strengths */
90 /* lower principal quantum number, 1, 2, or 3, in this routine
91 * 1 is 1s, 2 is 2s, 3 is 2p, and 4 is 3s, 5 is 3p, and 6 is 3d */
92 long int ipLow,
93 /* upper principal quantum nubmer, 2, 3, or 4 */
94 long int ipHi,
95 /* charge, 0 for hydrogen, 1 for helium, etc */
96 long int nelem,
97 /* = 'e' for electron collisions, ='p' for proton */
98 long int chType)
99{
100 long int i,
101 j,
102 k;
103 double C,
104 D,
105 EE,
106 expq ,
107 Hydcs123_v,
108 Ratehigh,
109 Ratelow,
110 TeUse,
111 gLo,
112 gHi,
113 q,
114 rate,
115 slope,
116 temp,
117 temphigh,
118 templow,
119 tev,
120 x,
121 QuanNLo,
122 QuanNUp,
123 Charge,
124 ChargeSquared,
125 zhigh,
126 zlow;
127 static const double ap[5] = {-2113.113,729.0084,1055.397,854.632,938.9912};
128 static const double bp[5] = {-6783.515,-377.7190,724.1936,493.1107,735.7466};
129 static const double cp[5] = {-3049.719,226.2320,637.8630,388.5465,554.6369};
130 static const double dp[5] = {3514.5153,88.60169,-470.4055,-329.4914,-450.8459};
131 static const double ep[5] = {0.005251557,0.009059154,0.008725781,0.009952418,0.01098687};
132 static const double ae[5] = {-767.5859,-643.1189,-461.6836,-429.0543,-406.5285};
133 static const double be[5] = {-1731.9178,-1442.548,-1055.364,-980.3079,-930.9266};
134 static const double ce[5] = {-939.1834,-789.9569,-569.1451,-530.1974,-502.0939};
135 static const double de[5] = {927.4773,773.2008,564.3272,524.2944,497.7763};
136 static const double ee[5] = {-0.002528027,-0.003793665,-0.002122103,-0.002234207,-0.002317720};
137 static const double A[2] = {4.4394,0.0};
138 static const double B[2] = {0.8949,0.8879};
139 static const double C0[2] = {-0.6012,-0.2474};
140 static const double C1[2] = {-3.9710,-3.7562};
141 static const double C2[2] = {-4.2176,2.0491};
142 static const double D0[2] = {2.930,0.0539};
143 static const double D1[2] = {1.7990,3.4009};
144 static const double D2[2] = {4.9347,-1.7770};
145
146 DEBUG_ENTRY( "Hydcs123()" );
147 /* Hydrogenic de-excitation collision rates n=1,2,3
148 * >>refer h1 cs Callaway, J. 1983, Phys Let A, 96, 83
149 * >>refer h1 cs Zygelman, B., & Dalgarno, A. 1987, Phys Rev A, 35, 4085
150 * for 2p-2s only.
151 * The fit from Callaway is in nuclear charge for 1s - 2s,2p only.
152 * For transtions involving level 3, interpolation in Z involving
153 * the functions He2cs123,C6cs123,Ne10cs123,Ca20cs123, Fe26cs123.
154 *
155 * The fits from ZD are for 2p-2s for Z=2,6,12,16,18 only other charges are
156 * interpolated, both electron and proton rates are included,
157 * the variable chType is either 'e' or 'p'.
158 *
159 * ipLow is the lower level and runs from 1 to 3 (1s, 2s, 2p)
160 * ipHi is the upper level and runs from 2 to 6 (2s, 2p, 3s, 3p, 3d) */
161
162 /* for Callaway fit: */
163 /* for Zygelman and Dalgarno: */
164
165 /* first entry is 2p, then 2s */
166
167 /* fit in nuclear charge Z for 2p-2s collisions in hydrogenic species
168 * equation is a+bx+cx^2ln(x)+dexp(x)+eln(x)/x^2, where x=te/Z^2 in a.u.
169 * first are the proton rates: */
170 /* these are electron rates: */
171
172 /* following is for charged species */
173 /* charge is on scale with He+=1, Li++=2, etc */
174 ASSERT( nelem > ipHYDROGEN );
175 ASSERT( nelem < LIMELM );
176
177 /* these are the pointers to upper and lower levels. 1=1s, 2=2s, 3=2p, 4=3 */
178 ASSERT( ipLow > 0);
179 ASSERT( ipLow <= 3);
180 ASSERT( ipHi > 1 );
181 ASSERT( ipHi <=6 );
182
183 /* set quantum numbers and stat. weights of the transitions: */
184 if( ipHi == 6 )
185 {
186 /* upper is n=3 then set level, stat. weight */
187 QuanNUp = 3.;
188 gHi = 10.;
189 /* following will be set here even though it is not used in this case,
190 * to prevent good compilers from falsing on i not set,
191 * there is assert when used to make sure it is ok */
192 i = -1;
193 }
194 else if( ipHi == 5 )
195 {
196 /* upper is n=3 then set level, stat. weight */
197 QuanNUp = 3.;
198 gHi = 6.;
199 /* following will be set here even though it is not used in this case,
200 * to prevent good compilers from falsing on i not set,
201 * there is assert when used to make sure it is ok */
202 i = -1;
203 }
204 else if( ipHi == 4 )
205 {
206 /* upper is n=3 then set level, stat. weight */
207 QuanNUp = 3.;
208 gHi = 2.;
209 /* following will be set here even though it is not used in this case,
210 * to prevent good compilers from falsing on i not set,
211 * there is assert when used to make sure it is ok */
212 i = -1;
213 }
214 else if( ipHi == 3 )
215 {
216 /* upper is nl=2p then set level, stat. weight */
217 QuanNUp = 2.;
218 gHi = 6.;
219 /* used to point within vectors defined above */
220 i = 0;
221 }
222 else if( ipHi == 2 )
223 {
224 /* upper is nl=2s then set level, stat. weight */
225 QuanNUp = 2.;
226 gHi = 2.;
227 /* used to point within vectors defined above */
228 i = 1;
229 }
230 else
231 {
232 fprintf( ioQQQ, " Insane levels in Hydcs123\n" );
234 }
235
236 /* which lower level? */
237 if( ipLow == 1 )
238 {
239 /* lower is n=1 then set level, stat. weight */
240 QuanNLo = 1.;
241 gLo = 2.;
242 }
243 else if( ipLow == 2 )
244 {
245 /* lower is nl=2s then set level, stat. weight */
246 QuanNLo = 2.;
247 gLo = 2.;
248 }
249 else if( ipLow == 3 )
250 {
251 QuanNLo = 2.;
252 gLo = 6.;
253 }
254 else
255 {
256 fprintf( ioQQQ, " Insane levels in Hydcs123\n" );
258 }
259
260 /* following is the physical charge */
261 Charge = (double)(nelem + 1);
262 /* square of charge */
263 ChargeSquared = Charge*Charge;
264
265 if( ipLow == 2 && ipHi == 3 )
266 {
267 /*************************************** this is 2p-2s:
268 * series of if statements determines which entries Charge is between. */
269 if( nelem < 1 )
270 {
271 /* this can't happen since routine returned above when ip=1 and
272 * special atomic hydrogen routine called */
273 fprintf( ioQQQ, " insane charge given to Hydcs123\n" );
275 }
276 else if( nelem == 1 )
277 {
278 zlow = 2.;
279 j = 1;
280 zhigh = 2.;
281 k = 1;
282 }
283 /* Li through C */
284 else if( nelem <= 5 )
285 {
286 zlow = 2.;
287 j = 1;
288 zhigh = 6.;
289 k = 2;
290 }
291 else if( nelem <= 11 )
292 {
293 zlow = 6.;
294 j = 2;
295 zhigh = 12.;
296 k = 3;
297 }
298 else if( nelem <= 15 )
299 {
300 zlow = 12.;
301 j = 3;
302 zhigh = 16.;
303 k = 4;
304 }
305 else if( nelem <= 17 )
306 {
307 zlow = 16.;
308 j = 4;
309 zhigh = 18.;
310 k = 5;
311 }
312 /* following changed to else from else if,
313 * to prevent false comment in good compilers */
314 /*else if( nelem > 18 )*/
315 else
316 {
317 zlow = 18.;
318 j = 5;
319 zhigh = 18.;
320 k = 5;
321 }
322
323 /* convert Te to a.u./Z^2
324 * determine rate at the low Charge */
325 x = EVRYD/TE1RYD*phycon.te/(27.211396*pow2(zlow));
326 TeUse = MIN2(x,0.80);
327 x = MAX2(0.025,TeUse);
328
329 /* what type of collision are we dealing with? */
330 if( chType == 'e' )
331 {
332 /* electron collisions */
333 Ratelow = ae[j-1] + be[j-1]*x + ce[j-1]*pow2(x)*log(x) + de[j-1]*
334 exp(x) + ee[j-1]*log(x)/pow2(x);
335 }
336 else if( chType == 'p' )
337 {
338 Ratelow = ap[j-1] + bp[j-1]*x + cp[j-1]*pow2(x)*log(x) + dp[j-1]*
339 exp(x) + ep[j-1]*log(x)/pow2(x);
340 }
341 else
342 {
343 /* this can't happen */
344 fprintf( ioQQQ, " insane collision species given to Hydcs123\n" );
346 }
347
348 /* determine rate at the high Charge */
349 x = EVRYD/TE1RYD*phycon.te/(27.211396*pow2(zhigh));
350 TeUse = MIN2(x,0.80);
351 x = MAX2(0.025,TeUse);
352 if( chType == 'e' )
353 {
354 Ratehigh = ae[k-1] + be[k-1]*x + ce[k-1]*pow2(x)*log(x) +
355 de[k-1]*exp(x) + ee[k-1]*log(x)/pow2(x);
356 }
357 else
358 {
359 Ratehigh = ap[k-1] + bp[k-1]*x + cp[k-1]*pow2(x)*log(x) +
360 dp[k-1]*exp(x) + ep[k-1]*log(x)/pow2(x);
361 }
362 /* linearly interpolate in charge */
363 if( fp_equal( zlow, zhigh ) )
364 {
365 rate = Ratelow;
366 }
367 else
368 {
369 slope = (Ratehigh - Ratelow)/(zhigh - zlow);
370 rate = slope*(Charge - zlow) + Ratelow;
371 }
372 rate = rate/ChargeSquared/Charge*1.0e-7;
373 /* must convert to cs and need to know the valid temp range */
374 templow = 0.025*27.211396*TE1RYD/EVRYD*ChargeSquared;
375 temphigh = 0.80*27.211396*TE1RYD/EVRYD*ChargeSquared;
376 TeUse = MIN2((double)phycon.te,temphigh);
377 temp = MAX2(TeUse,templow);
378 Hydcs123_v = rate*gHi*sqrt(temp)/COLL_CONST;
379
380 if( chType == 'p' )
381 {
382 /* COLL_CONST is incorrect for protons, correct here */
383 Hydcs123_v *= pow( PROTON_MASS/ELECTRON_MASS, 1.5 );
384 }
385 }
386 else if( ipHi == 4 || ipHi == 5 || ipHi == 6 )
387 {
388 /* n = 3
389 * for the rates involving n = 3, must do something different. */
390 if( nelem < 1 )
391 {
392 fprintf( ioQQQ, " insane charge given to Hydcs123\n" );
394 }
395 else if( nelem == 1 )
396 {
397 zlow = 2.;
398 Ratelow = He2cs123(ipLow,ipHi);
399 zhigh = 2.;
400 Ratehigh = Ratelow;
401 }
402 else if( nelem > 1 && nelem <= 5 )
403 {
404 zlow = 2.;
405 Ratelow = He2cs123(ipLow,ipHi);
406 zhigh = 6.;
407 Ratehigh = C6cs123(ipLow,ipHi);
408 }
409 else if( nelem > 5 && nelem <= 9 )
410 {
411 zlow = 6.;
412 Ratelow = C6cs123(ipLow,ipHi);
413 zhigh = 10.;
414 Ratehigh = Ne10cs123(ipLow,ipHi);
415 }
416 else if( nelem > 9 && nelem <= 19 )
417 {
418 zlow = 10.;
419 Ratelow = Ne10cs123(ipLow,ipHi);
420 zhigh = 20.;
421 Ratehigh = Ca20cs123(ipLow,ipHi);
422 }
423 else if( nelem > 19 && nelem <= 25 )
424 {
425 zlow = 20.;
426 Ratelow = Ca20cs123(ipLow,ipHi);
427 zhigh = 26.;
428 Ratehigh = Fe26cs123(ipLow,ipHi);
429 }
430 /*>>>chng 98 dec 17, to else to stop comment from good compilers*/
431 /*else if( nelem > 26 )*/
432 else
433 {
434 Charge = 26.;
435 zlow = 26.;
436 Ratelow = Fe26cs123(ipLow,ipHi);
437 zhigh = 26.;
438 Ratehigh = Ratelow;
439 }
440
441 /* linearly interpolate */
442 if( fp_equal( zlow, zhigh ) )
443 {
444 rate = Ratelow;
445 }
446 else
447 {
448 slope = (Ratehigh - Ratelow)/(zhigh - zlow);
449 rate = slope*(Charge - zlow) + Ratelow;
450 }
451
453 Hydcs123_v = rate;
454
455 /* the routines C6cs123, Ne10cs123, etc... do not resolve L for n>2 */
456 /* dividing by N should roughly recover the original l-resolved data */
457 Hydcs123_v /= 3.;
458 }
459 else
460 {
461 /* this branch 1-2s, 1-2p */
462 if( nelem == 1 )
463 {
464 /* this brance for helium, then return */
465 Hydcs123_v = He2cs123(ipLow,ipHi);
466 return( Hydcs123_v );
467 }
468
469 /* electron temperature in eV */
470 tev = phycon.te / EVDEGK;
471 /* energy in eV for hydrogenic species and these quantum numbers */
472 EE = ChargeSquared*EVRYD*(1./QuanNLo/QuanNLo - 1./QuanNUp/QuanNUp);
473 /* EE/kT for this transion */
474 q = EE/tev;
475 TeUse = MIN2(q,10.);
476 /* q is now EE/kT but between 1 and 10 */
477 q = MAX2(1.,TeUse);
478 expq = exp(q);
479
480 /* i must be 0 or 1 */
481 ASSERT( i==0 || i==1 );
482 C = C0[i] + C1[i]/Charge + C2[i]/ChargeSquared;
483 D = D0[i] + D1[i]/Charge + D2[i]/ChargeSquared;
484
485 /* following code changed so that ee1 always returns e1,
486 * orifinal version only returned e1 for x < 1 */
487 /* use disabled e1: */
488 /*if( q < 1. )*/
489 /*{*/
490 /*rate = (B[i-1] + D*q)*exp(-q) + (A[i-1] + C*q - D*q*q)**/
491 /*ee1(q);*/
492 /*}*/
493 /*else*/
494 /*{*/
495 /*rate = (B[i-1] + D*q) + (A[i-1] + C*q - D*q*q)*ee1(q);*/
496 /*}*/
497 /*rate *= 8.010e-8/2./ChargeSquared/tev*sqrt(tev);*/
498 /* convert to de-excitation */
499 /*if( q < 1. )*/
500 /*{*/
501 /*rate = rate*exp(q)*gLo/gHi;*/
502 /*}*/
503 /*else*/
504 /*{*/
505 /*rate = rate*gLo/gHi;*/
506 /*}*/
507
508 /*>>>chng 98 dec 17, ee1 always returns e1 */
509 rate = (B[i] + D*q)/expq + (A[i] + C*q - D*q*q)*
510 ee1(q);
511 rate *= 8.010e-8/2./ChargeSquared/tev*sqrt(tev);
512 /* convert to de-excitation */
513 rate *= expq*gLo/gHi;
514
515 /* convert to cs */
516 Hydcs123_v = rate*gHi*phycon.sqrte/COLL_CONST;
517 }
518 return( Hydcs123_v );
519}
520
521/*C6cs123 line collision rates for lower levels of hydrogenic carbon, n=1,2,3 */
522STATIC double C6cs123(long int i,
523 long int j)
524{
525 double C6cs123_v,
526 TeUse,
527 t,
528 x;
529 static const double a[3] = {-92.23774,-1631.3878,-6326.4947};
530 static const double b[3] = {-11.93818,-218.3341,-849.8927};
531 static const double c[3] = {0.07762914,1.50127,5.847452};
532 static const double d[3] = {78.401154,1404.8475,5457.9291};
533 static const double e[3] = {332.9531,5887.4263,22815.211};
534
535 DEBUG_ENTRY( "C6cs123()" );
536
537 /* These are fits to Table 5 of
538 * >>refer c6 cs Aggarwal, K.M., & Kingston, A.E. 1991, J Phys B, 24, 4583
539 * C VI collision rates for 1s-3l, 2s-3l, and 2p-3l,
540 * principal quantum numbers n and l.
541 *
542 * i is the lower level and runs from 1 to 3 (1s, 2s, 2p)
543 * j is the upper level and runs from 2 to 6 (2s, 2p, 3s, 3p, 3d)
544 * 1s-2s,2p is not done here.
545 * check temperature: fits only good between 3.8 < log Te < 6.2
546 */
547 /* arrays for fits of 3 transitions see the code below for key: */
548
549 TeUse = MAX2(phycon.te,6310.);
550 t = MIN2(TeUse,1.6e6);
551 x = log10(t);
552
553 if( i == 1 && j == 2 )
554 {
555 /* 1s - 2s (first entry) */
556 fprintf( ioQQQ, " Carbon VI 2s-1s not done in C6cs123\n" );
558 }
559
560 else if( i == 1 && j == 3 )
561 {
562 /* 1s - 2p (second entry) */
563 fprintf( ioQQQ, " Carbon VI 2p-1s not done in C6cs123\n" );
565 }
566
567 else if( i == 1 && ( j == 4 || j == 5 || j == 6 ) )
568 {
569 /* 1s - 3 (first entry) */
570 C6cs123_v = a[0] + b[0]*x + c[0]*pow2(x)*sqrt(x) + d[0]*log(x) +
571 e[0]*log(x)/pow2(x);
572 }
573 else if( i == 2 && ( j == 4 || j == 5 || j == 6 ) )
574 {
575 /* 2s - 3 (second entry) */
576 C6cs123_v = a[1] + b[1]*x + c[1]*pow2(x)*sqrt(x) + d[1]*log(x) +
577 e[1]*log(x)/pow2(x);
578 }
579 else if( i == 3 && ( j == 4 || j == 5 || j == 6 ) )
580 {
581 /* 2p - 3s (third entry) */
582 C6cs123_v = a[2] + b[2]*x + c[2]*pow2(x)*sqrt(x) + d[2]*log(x) +
583 e[2]*log(x)/pow2(x);
584 }
585 else
586 {
587 fprintf( ioQQQ, " insane levels for C VI n=1,2,3 !!!\n" );
589 }
590 return( C6cs123_v );
591}
592
593/*Ca20cs123 line collision rates for lower levels of hydrogenic calcium, n=1,2,3 */
594STATIC double Ca20cs123(long int i,
595 long int j)
596{
597 double Ca20cs123_v,
598 TeUse,
599 t,
600 x;
601 static const double a[3] = {-12.5007,-187.2303,-880.18896};
602 static const double b[3] = {-1.438749,-22.17467,-103.1259};
603 static const double c[3] = {0.008219688,0.1318711,0.6043752};
604 static const double d[3] = {10.116516,153.2650,717.4036};
605 static const double e[3] = {45.905343,685.7049,3227.2836};
606
607 DEBUG_ENTRY( "Ca20cs123()" );
608
609 /*
610 * These are fits to Table 5 of
611 * >>refer ca20 cs Aggarwal, K.M., & Kingston, A.E. 1992, J Phys B, 25, 751
612 * Ca XX collision rates for 1s-3l, 2s-3l, and 2p-3l,
613 * principal quantum numbers n and l.
614 *
615 * i is the lower level and runs from 1 to 3 (1s, 2s, 2p)
616 * j is the upper level and runs from 2 to 6 (2s, 2p, 3s, 3p, 3d)
617 * 1s-2s,2p is not done here.
618 * check temperature: fits only good between 5.0 < log Te < 7.2
619 */
620
621 /* arrays for fits of 3 transitions see the code below for key: */
622
623 TeUse = MAX2(phycon.te,1.0e5);
624 t = MIN2(TeUse,1.585e7);
625 x = log10(t);
626
627 if( i == 1 && j == 2 )
628 {
629 /* 1s - 2s (first entry) */
630 fprintf( ioQQQ, " Ca XX 2s-1s not done in Ca20cs123\n" );
632 }
633
634 else if( i == 1 && j == 3 )
635 {
636 /* 1s - 2p (second entry) */
637 fprintf( ioQQQ, " Ca XX 2p-1s not done in Ca20cs123\n" );
639 }
640
641 else if( i == 1 && ( j == 4 || j == 5 || j == 6 ))
642 {
643 /* 1s - 3 (first entry) */
644 Ca20cs123_v = a[0] + b[0]*x + c[0]*pow2(x)*sqrt(x) + d[0]*log(x) +
645 e[0]*log(x)/pow2(x);
646 }
647 else if( i == 2 && ( j == 4 || j == 5 || j == 6 ))
648 {
649 /* 2s - 3 (second entry) */
650 Ca20cs123_v = a[1] + b[1]*x + c[1]*pow2(x)*sqrt(x) + d[1]*log(x) +
651 e[1]*log(x)/pow2(x);
652 }
653 else if( i == 3 && ( j == 4 || j == 5 || j == 6 ))
654 {
655 /* 2p - 3s (third entry) */
656 Ca20cs123_v = a[2] + b[2]*x + c[2]*pow2(x)*sqrt(x) + d[2]*log(x) +
657 e[2]*log(x)/pow2(x);
658 }
659 else
660 {
661 fprintf( ioQQQ, " insane levels for Ca XX n=1,2,3 !!!\n" );
663 }
664 return( Ca20cs123_v );
665}
666
667/*Ne10cs123 line collision rates for lower levels of hydrogenic neon, n=1,2,3 */
668STATIC double Ne10cs123(long int i,
669 long int j)
670{
671 double Ne10cs123_v,
672 TeUse,
673 t,
674 x;
675 static const double a[3] = {3.346644,151.2435,71.7095};
676 static const double b[3] = {0.5176036,20.05133,13.1543};
677 static const double c[3] = {-0.00408072,-0.1311591,-0.1099238};
678 static const double d[3] = {-3.064742,-129.8303,-71.0617};
679 static const double e[3] = {-11.87587,-541.8599,-241.2520};
680
681 DEBUG_ENTRY( "Ne10cs123()" );
682
683 /* These are fits to Table 5 of
684 * >>refer ne10 cs Aggarwal, K.M., & Kingston, A.E. 1991, PhyS, 44, 517
685 * Ne X collision rates for 1s-3, 2s-3l, and 2p-3l,
686 * principal quantum numbers n and l.
687 *
688 * i is the lower level and runs from 1 to 3 (1s, 2s, 2p)
689 * j is the upper level and runs from 2 to 6 (2s, 2p, 3s, 3p, 3d)
690 * 1s-2s,2p is not done here.
691 * check temperature: fits only good between 3.8 < log Te < 6.2
692 * */
693 /* arrays for fits of 3 transitions see the code below for key: */
694
695 TeUse = MAX2(phycon.te,6310.);
696 t = MIN2(TeUse,1.6e6);
697 x = log10(t);
698
699 if( i == 1 && j == 2 )
700 {
701 /* 1s - 2s (first entry) */
702 fprintf( ioQQQ, " Neon X 2s-1s not done in Ne10cs123\n" );
704 }
705
706 else if( i == 1 && j == 3 )
707 {
708 /* 1s - 2p (second entry) */
709 fprintf( ioQQQ, " Neon X 2p-1s not done in Ne10cs123\n" );
711 }
712
713 else if( i == 1 && ( j == 4 || j == 5 || j == 6 ) )
714 {
715 /* 1s - 3 (first entry) */
716 Ne10cs123_v = a[0] + b[0]*x + c[0]*pow2(x)*sqrt(x) + d[0]*log(x) +
717 e[0]*log(x)/pow2(x);
718 }
719 else if( i == 2 && ( j == 4 || j == 5 || j == 6 ) )
720 {
721 /* 2s - 3 (second entry) */
722 Ne10cs123_v = a[1] + b[1]*x + c[1]*pow2(x)*sqrt(x) + d[1]*log(x) +
723 e[1]*log(x)/pow2(x);
724 }
725 else if( i == 3 && ( j == 4 || j == 5 || j == 6 ) )
726 {
727 /* 2p - 3s (third entry) */
728 Ne10cs123_v = a[2] + b[2]*x + c[2]*pow2(x)*sqrt(x) + d[2]*log(x) +
729 e[2]*log(x)/pow2(x);
730 }
731 else
732 {
733 fprintf( ioQQQ, " insane levels for Ne X n=1,2,3 !!!\n" );
735 }
736 return( Ne10cs123_v );
737}
738
739/*He2cs123 line collision strengths for lower levels of helium ion, n=1,2,3, by K Korista */
740STATIC double He2cs123(long int i,
741 long int j)
742{
743 double He2cs123_v,
744 t;
745 static const double a[11]={0.12176209,0.32916723,0.46546497,0.044501688,
746 0.040523277,0.5234889,1.4903214,1.4215094,1.0295881,4.769306,9.7226127};
747 static const double b[11]={0.039936166,2.9711166e-05,-0.020835863,3.0508137e-04,
748 -2.004485e-15,4.41475e-06,1.0622666e-05,2.0538877e-06,0.80638448,2.0967075e-06,
749 7.6089851e-05};
750 static const double c[11]={143284.77,0.73158545,-2.159172,0.43254802,2.1338557,
751 8.9899702e-06,-2.9001451e-12,1.762076e-05,52741.735,-2153.1219,-3.3996921e-11};
752
753 DEBUG_ENTRY( "He2cs123()" );
754
755 /* These are fits to Table 2
756 * >>refer he2 cs Aggarwal, K.M., Callaway, J., Kingston, A.E., Unnikrishnan, K.
757 * >>refercon 1992, ApJS, 80, 473
758 * He II collision rates for 1s-2s, 1s-2p, 1s-3s, 1s-3p, 1s-3d, 2s-3s, 2s-3p, 2s-3d,
759 * 2p-3s, 2p-3p, and 2p-3d.
760 * principal quantum numbers n and l.
761 *
762 * i is the lower level and runs from 1 to 3 (1s, 2s, 2p)
763 * j is the upper level and runs from 2 to 6 (2s, 2p, 3s, 3p, 3d)
764 * check temperature: fits only good between 5,000K and 500,000K
765 * */
766 /* array for fits of 11 transitions see the code below for key: */
767
768 t = phycon.te;
769 if( t < 5000. )
770 {
771 t = 5000.;
772 }
773 else if( t > 5.0e05 )
774 {
775 t = 5.0e05;
776 }
777
778 /**************fits begin here**************
779 * */
780 if( i == 1 && j == 2 )
781 {
782 /* 1s - 2s (first entry) */
783 He2cs123_v = a[0] + b[0]*exp(-t/c[0]);
784 }
785 else if( i == 1 && j == 3 )
786 {
787 /* 1s - 2p (second entry) */
788 He2cs123_v = a[1] + b[1]*pow(t,c[1]);
789 }
790 else if( i == 1 && j == 4 )
791 {
792 /* 1s - 3s (third entry) */
793 He2cs123_v = a[2] + b[2]*log(t) + c[2]/log(t);
794 }
795 else if( i == 1 && j == 5 )
796 {
797 /* 1s - 3p (fourth entry) */
798 He2cs123_v = a[3] + b[3]*pow(t,c[3]);
799 }
800 else if( i == 1 && j == 6 )
801 {
802 /* 1s - 3d (fifth entry) */
803 He2cs123_v = a[4] + b[4]*pow(t,c[4]);
804 }
805 else if( i == 2 && j == 4 )
806 {
807 /* 2s - 3s (sixth entry) */
808 He2cs123_v = (a[5] + c[5]*t)/(1 + b[5]*t);
809 }
810 else if( i == 2 && j == 5 )
811 {
812 /* 2s - 3p (seventh entry) */
813 He2cs123_v = a[6] + b[6]*t + c[6]*t*t;
814 }
815 else if( i == 2 && j == 6 )
816 {
817 /* 2s - 3d (eighth entry) */
818 He2cs123_v = (a[7] + c[7]*t)/(1 + b[7]*t);
819 }
820 else if( i == 3 && j == 4 )
821 {
822 /* 2p - 3s (ninth entry) */
823 He2cs123_v = a[8] + b[8]*exp(-t/c[8]);
824 }
825 else if( i == 3 && j == 5 )
826 {
827 /* 2p - 3p (tenth entry) */
828 He2cs123_v = a[9] + b[9]*t + c[9]/t;
829 }
830 else if( i == 3 && j == 6 )
831 {
832 /* 2p - 3d (eleventh entry) */
833 He2cs123_v = a[10] + b[10]*t + c[10]*t*t;
834 }
835 else
836 {
837 /**************fits end here************** */
838 fprintf( ioQQQ, " insane levels for He II n=1,2,3 !!!\n" );
840 }
841 return( He2cs123_v );
842}
843
844/*Fe26cs123 line collision rates for lower levels of hydrogenic iron, n=1,2,3 */
845STATIC double Fe26cs123(long int i,
846 long int j)
847{
848 double Fe26cs123_v,
849 TeUse,
850 t,
851 x;
852 static const double a[3] = {-4.238398,-238.2599,-1211.5237};
853 static const double b[3] = {-0.4448177,-27.06869,-136.7659};
854 static const double c[3] = {0.0022861,0.153273,0.7677703};
855 static const double d[3] = {3.303775,191.7165,972.3731};
856 static const double e[3] = {15.82689,878.1333,4468.696};
857
858 DEBUG_ENTRY( "Fe26cs123()" );
859
860 /* These are fits to Table 5 of
861 * >>refer fe26 cs Aggarwal, K.M., & Kingston, A.E. 1993, ApJS, 85, 187
862 * Fe XXVI collision rates for 1s-3, 2s-3, and 2p-3,
863 * principal quantum numbers n and l.
864 *
865 * i is the lower level and runs from 1 to 3 (1s, 2s, 2p)
866 * j is the upper level and runs from 2 to 6 (2s, 2p, 3s, 3p, 3d)
867 * 1s-2s,2p is not done here.
868 * check temperature: fits only good between 5.2 < log Te < 7.2
869 * */
870 /* arrays for fits of 3 transitions see the code below for key: */
871
872 TeUse = MAX2(phycon.te,1.585e5);
873 t = MIN2(TeUse,1.585e7);
874 x = log10(t);
875
876 if( i == 1 && j == 2 )
877 {
878 /* 1s - 2s (first entry) */
879 fprintf( ioQQQ, " Fe XXVI 2s-1s not done in Fe26cs123\n" );
881 }
882
883 else if( i == 1 && j == 3 )
884 {
885 /* 1s - 2p (second entry) */
886 fprintf( ioQQQ, " Fe XXVI 2p-1s not done in Fe26cs123\n" );
888 }
889
890 else if( i == 1 && ( j == 4 || j == 5 || j == 6 ) )
891 {
892 /* 1s - 3 (first entry) */
893 Fe26cs123_v = a[0] + b[0]*x + c[0]*pow2(x)*sqrt(x) + d[0]*log(x) +
894 e[0]*log(x)/pow2(x);
895 }
896 else if( i == 2 && ( j == 4 || j == 5 || j == 6 ) )
897 {
898 /* 2s - 3 (second entry) */
899 Fe26cs123_v = a[1] + b[1]*x + c[1]*pow2(x)*sqrt(x) + d[1]*log(x) +
900 e[1]*log(x)/pow2(x);
901 }
902 else if( i == 3 && ( j == 4 || j == 5 || j == 6 ) )
903 {
904 /* 2p - 3s (third entry) */
905 Fe26cs123_v = a[2] + b[2]*x + c[2]*pow2(x)*sqrt(x) + d[2]*log(x) +
906 e[2]*log(x)/pow2(x);
907 }
908 else
909 {
910 fprintf( ioQQQ, " insane levels for Ca XX n=1,2,3 !!!\n" );
912 }
913 return( Fe26cs123_v );
914}
915
916
917STATIC double CS_ThermAve_PR78(long ipISO, long nelem, long nHi, long nLo, double deltaE, double temp )
918{
919
920 double coll_str;
921
922 DEBUG_ENTRY( "CS_ThermAve_PR78()" );
923
924 global_ipISO = ipISO;
925 global_nelem = nelem;
926 global_nHi = nHi;
927 global_nLo = nLo;
928 global_deltaE = deltaE;
929
930 kTRyd = temp / TE1RYD;
931
932 if( !iso_ctrl.lgCS_therm_ave[ipISO] )
933 {
934 /* Must do some thermal averaging for densities greater
935 * than about 10000 and less than about 1e10,
936 * because kT gives significantly different results.
937 * Still, do sparser integration than is done below */
938 if( (dense.eden > 10000.) && (dense.eden < 1E10 ) )
939 {
940 coll_str = qg32( 0.0, 6.0, Therm_ave_coll_str_int_PR78);
941 }
942 else
943 {
944 /* Do NOT average over Maxwellian */
945 coll_str = CS_PercivalRichards78( kTRyd );
946 }
947 }
948 else
949 {
950 /* DO average over Maxwellian */
951 coll_str = qg32( 0.0, 1.0, Therm_ave_coll_str_int_PR78);
952 coll_str += qg32( 1.0, 10.0, Therm_ave_coll_str_int_PR78);
953 }
954
955 return coll_str;
956}
957
958/* The integrand for calculating the thermal average of collision strengths */
959STATIC double Therm_ave_coll_str_int_PR78( double EOverKT )
960{
961 double integrand;
962
963 DEBUG_ENTRY( "Therm_ave_coll_str_int_PR78()" );
964
965 integrand = exp( -1.*EOverKT ) * CS_PercivalRichards78( EOverKT * kTRyd );
966
967 return integrand;
968}
969
970STATIC double CS_PercivalRichards78( double Ebar )
971{
972 double cross_section, coll_str;
973 double stat_weight;
974 double A, D, L, F, G, H;
975 double np, n, s, Z_3, xPlus, xMinus, y;
976 long ipISO, nelem;
977
978 DEBUG_ENTRY( "CS_PercivalRichards78()" );
979
980 if( Ebar < global_deltaE )
981 {
982 DEBUG_ENTRY( "CS_PercivalRichards78()" );
983 return 0.;
984 }
985
986 ipISO = global_ipISO;
987 nelem = global_nelem;
988 np = (double)global_nHi;
989 n = (double)global_nLo;
990 s = np - n;
991 ASSERT( s > 0. );
992
993 A = (8./3./s) * pow(np/s/n, 3.) * (0.184 - 0.04 * pow( s, -0.66)) * pow( 1. - 0.2*s/n/np, 1.+2.*s);
994
995 Z_3 = (double)(nelem + 1. - ipISO);
996
997 D = exp( - Z_3 * Z_3 / n / np / Ebar / Ebar );
998
999 L = log( (1. + 0.53 * Ebar * Ebar * n * np / Z_3 / Z_3) / (1. + 0.4*Ebar) );
1000
1001 F = pow( 1. - 0.3 * s * D / n /np, 1. + 2.*s );
1002
1003 G = 0.5* POW3( Ebar * n * n / Z_3 / np );
1004
1005 xPlus = 2. * Z_3 / ( n * n * Ebar * ( sqrt( 2. - n*n/np/np ) + 1. ) );
1006 xMinus = 2. * Z_3 / ( n * n * Ebar * ( sqrt( 2. - n*n/np/np ) - 1. ) );
1007
1008 y = 1. / (1. - D * log ( 18. * s )/ 4. / s);
1009
1010 H = POW2( xMinus) * log( 1. + 2.*xMinus/3. ) / ( 2.*y + 1.5*xMinus );
1011 H -= POW2( xPlus ) * log( 1. + 2.* xPlus/3. ) / ( 2.*y + 1.5*xPlus );
1012
1013 /* this is the LHS of equation 1 of PR78 */
1014 cross_section = (A*D*L + F*G*H);
1015 /* this is the result after solving equation 1 for the cross section */
1016 cross_section *= PI * POW2( n * n * BOHR_RADIUS_CM / Z_3 ) / Ebar;
1017
1018 if( ipISO == ipH_LIKE )
1019 stat_weight = 2. * n * n;
1020 else if( ipISO == ipHE_LIKE )
1021 stat_weight = 4. * n * n;
1022 else
1023 TotalInsanity();
1024
1025 /* convert to collision strength */
1026 coll_str = cross_section * stat_weight * Ebar / ( PI * POW2( BOHR_RADIUS_CM ) );
1027 return coll_str;
1028}
1029
1030#if 0
1031STATIC void TestPercivalRichards( void )
1032{
1033 double CStemp;
1034
1035 /* this reproduces Table 1 of PR78 */
1036 for( long i=0; i<5; i++ )
1037 {
1038 double Ebar[5] = {0.1, 0.4, 0.8, 1.0, 10.};
1039
1040 CStemp = CS_PercivalRichards78( 0, 2, 12, 10, Ebar[i] );
1041 }
1042
1043 /* this reproduces Table 2 of PR78 */
1044 for( long i=0; i<5; i++ )
1045 {
1046 double Ebar[5] = {0.1, 0.4, 0.8, 1.0, 10.};
1047
1048 CStemp = CS_ThermAve_PR78( ipISO, 0, N_(ipHi), N_(ipLo), phycon.te );
1049 }
1050
1051 return;
1052}
1053#endif
1054
1055realnum HydroCSInterp(long int nelem,
1056 long int ipHi,
1057 long int ipLo,
1058 long int ipCollider )
1059{
1060 double CStemp;
1061 long ipISO = ipH_LIKE;
1062
1063 DEBUG_ENTRY( "HydroCSInterp()" );
1064
1065 if( !iso_ctrl.lgColl_excite[ipISO] )
1066 return 0.;
1067
1068 /* This set of collision strengths should only be used
1069 * if the Storey and Hummer flag is set */
1070 if( opac.lgCaseB_HummerStorey )
1071 {
1072 if( N_(ipLo) == N_(ipHi) )
1073 {
1074 if( N_(ipHi) <= iso_sp[ipH_LIKE][nelem].n_HighestResolved_max &&
1075 abs( L_(ipLo) - L_(ipHi) ) != 1 )
1076 {
1077 /* if delta L is not +/- 1, set collision strength to zero. */
1078 CStemp = 0.;
1079 }
1080 else
1081 {
1082 CStemp = CS_l_mixing_PS64(
1083 nelem,
1084 iso_sp[ipH_LIKE][nelem].st[ipLo].lifetime(),
1085 nelem+1.-ipH_LIKE,
1086 iso_sp[ipH_LIKE][nelem].st[ipLo].n(),
1087 iso_sp[ipH_LIKE][nelem].st[ipLo].l(),
1088 iso_sp[ipH_LIKE][nelem].st[ipHi].g(),
1089 ipCollider);
1090 }
1091 }
1092
1093 else
1094 {
1095 if( N_(ipHi) <= iso_sp[ipH_LIKE][nelem].n_HighestResolved_max &&
1096 abs( L_(ipLo) - L_(ipHi) ) != 1 )
1097 {
1098 /* if delta L is not +/- 1, set collision strength to zero. */
1099 CStemp = 0.;
1100 }
1101 else if( ipCollider == ipELECTRON )
1102 {
1103 CStemp = CS_ThermAve_PR78( ipH_LIKE, nelem, N_(ipHi), N_(ipLo),
1104 iso_sp[ipH_LIKE][nelem].trans(ipHi,ipLo).EnergyErg() / EN1RYD , phycon.te );
1105 }
1106 else
1107 CStemp = 0.;
1108 }
1109 }
1110 else
1111 {
1112 /* HCSAR_interp interpolates on a table to return R-matrix collision strengths
1113 * for hydrogen only */
1114 if( nelem==ipHYDROGEN && ipCollider==ipELECTRON && N_(ipHi) <= 5 && ( N_(ipHi) != N_(ipLo) ) )
1115 {
1116 CStemp = HCSAR_interp(ipLo,ipHi);
1117 }
1118 else if( nelem==ipHYDROGEN && ipCollider==ipELECTRON && ( N_(ipHi) != N_(ipLo) ) )
1119 {
1120 CStemp = hydro_vs_deexcit( ipH_LIKE, nelem, ipHi, ipLo, iso_sp[ipH_LIKE][nelem].trans(ipHi,ipLo).Emis().Aul() );
1121 }
1122 else if( nelem>ipHYDROGEN && ipCollider==ipELECTRON && N_(ipHi) <= 3 && N_(ipLo) < 3 )
1123 {
1124 CStemp = Hydcs123(ipLo + 1,ipHi + 1,nelem,'e');
1125 }
1126 else if( nelem>ipHYDROGEN && ipCollider==ipPROTON && ipHi==ipH2p && ipLo==ipH2s )
1127 {
1128 CStemp = Hydcs123(ipLo + 1,ipHi + 1,nelem,'p');
1129 }
1130 else if( N_(ipLo) == N_(ipHi) )
1131 {
1132 if( iso_ctrl.lgCS_Vrinceanu[ipH_LIKE] )
1133 {
1134 CStemp = CS_l_mixing_VF01(ipH_LIKE, nelem,
1135 iso_sp[ipH_LIKE][nelem].st[ipLo].n(),
1136 iso_sp[ipH_LIKE][nelem].st[ipLo].l(),
1137 iso_sp[ipH_LIKE][nelem].st[ipHi].l(),
1138 iso_sp[ipH_LIKE][nelem].st[ipLo].S(),
1139 phycon.te,
1140 ipCollider );
1141 }
1142 else if (abs( L_(ipLo) - L_(ipHi) ) == 1)
1143 CStemp = CS_l_mixing_PS64(
1144 nelem,
1145 iso_sp[ipH_LIKE][nelem].st[ipLo].lifetime(),
1146 nelem+1.-ipH_LIKE,
1147 iso_sp[ipH_LIKE][nelem].st[ipLo].n(),
1148 iso_sp[ipH_LIKE][nelem].st[ipLo].l(),
1149 iso_sp[ipH_LIKE][nelem].st[ipHi].g(),
1150 ipCollider);
1151 else
1152 CStemp = 0.;
1153 }
1154 else
1155 {
1156 ASSERT( N_(ipHi) != N_(ipLo) );
1157 /* highly excited levels */
1158 CStemp = CS_VS80( ipH_LIKE, nelem, ipHi, ipLo,
1159 iso_sp[ipH_LIKE][nelem].trans(ipHi,ipLo).Emis().Aul(),
1160 phycon.te,
1161 ipCollider );
1162 }
1163 }
1164
1165 return (realnum)CStemp;
1166}
#define NHCSTE
Definition atmdat.h:126
FILE * ioQQQ
Definition cddefines.cpp:7
#define ASSERT(exp)
Definition cddefines.h:578
#define STATIC
Definition cddefines.h:97
#define POW3
Definition cddefines.h:936
T pow2(T a)
Definition cddefines.h:931
#define MIN2
Definition cddefines.h:761
double qg32(double, double, double(*)(double))
Definition service.cpp:1053
double ee1(double x)
Definition service.cpp:312
const int LIMELM
Definition cddefines.h:258
#define POW2
Definition cddefines.h:929
#define EXIT_FAILURE
Definition cddefines.h:140
#define cdEXIT(FAIL)
Definition cddefines.h:434
float realnum
Definition cddefines.h:103
#define MAX2
Definition cddefines.h:782
bool fp_equal(sys_float x, sys_float y, int n=3)
Definition cddefines.h:812
NORETURN void TotalInsanity(void)
Definition service.cpp:886
const int ipHYDROGEN
Definition cddefines.h:305
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
static t_ADfA & Inst()
Definition cddefines.h:175
realnum h_coll_str(long ipLo, long ipHi, long ipTe)
@ ipPROTON
Definition collision.h:10
@ ipELECTRON
Definition collision.h:9
t_dense dense
Definition dense.cpp:24
double CS_l_mixing_VF01(long int ipISO, long int nelem, long int n, long int l, long int lp, long int s, double temp, long int Collider)
double CS_l_mixing_PS64(long nelem, double tau, double target_charge, long n, long l, double gHi, long Collider)
STATIC double cross_section(double EgammaRyd, double EthRyd, long nelem, long n, long l, long s)
double hydro_vs_deexcit(long ipISO, long nelem, long ipHi, long ipLo, double Aul)
double CS_VS80(long int ipISO, long int nelem, long int ipHi, long int ipLo, double Aul, double temp, long int Collider)
realnum HydroCSInterp(long int nelem, long int ipHi, long int ipLo, long int ipCollider)
static const realnum HCSTE[NHCSTE]
static double global_deltaE
STATIC double He2cs123(long int i, long int j)
static long global_nLo
static long global_ipISO
STATIC double Therm_ave_coll_str_int_PR78(double EOverKT)
static double kTRyd
STATIC double Ca20cs123(long int i, long int j)
STATIC double CS_ThermAve_PR78(long ipISO, long nelem, long nHi, long nLo, double deltaE, double temp)
static long global_nHi
STATIC double CS_PercivalRichards78(double Ebar)
static long global_nelem
STATIC double Hydcs123(long int ilow, long int ihigh, long int iz, long int chType)
STATIC realnum HCSAR_interp(int ipLo, int ipHi)
STATIC double Fe26cs123(long int i, long int j)
STATIC double C6cs123(long int i, long int j)
STATIC double Ne10cs123(long int i, long int j)
t_iso_sp iso_sp[NISO][LIMELM]
Definition iso.cpp:8
t_isoCTRL iso_ctrl
Definition iso.cpp:6
#define N_(A_)
Definition iso.h:20
const int ipHE_LIKE
Definition iso.h:63
const int ipH2p
Definition iso.h:29
const int ipH2s
Definition iso.h:28
const int ipH_LIKE
Definition iso.h:62
#define L_(A_)
Definition iso.h:21
t_opac opac
Definition opacity.cpp:5
#define S(I_, J_)
t_phycon phycon
Definition phycon.cpp:6
UNUSED const double PI
Definition physconst.h:29
UNUSED const double BOHR_RADIUS_CM
Definition physconst.h:222
UNUSED const double ELECTRON_MASS
Definition physconst.h:91
UNUSED const double EN1RYD
Definition physconst.h:179
UNUSED const double EVRYD
Definition physconst.h:189
UNUSED const double PROTON_MASS
Definition physconst.h:94
UNUSED const double EVDEGK
Definition physconst.h:186
UNUSED const double COLL_CONST
Definition physconst.h:229
UNUSED const double EE
Definition physconst.h:23
UNUSED const double TE1RYD
Definition physconst.h:183
static double * g
Definition species2.cpp:28
static const double C1