cloudy trunk
Loading...
Searching...
No Matches
opacity_add1subshell.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/*OpacityAdd1Subshell add opacity due to single shell to main opacity array*/
4/*OpacityAdd1SubshellInduc add opacity of individual species, including stimulated emission */
5#include "cddefines.h"
6#include "rfield.h"
7#include "hydrogenic.h"
8#include "opacity.h"
9
11 /*ipOpac is opacity index within opac opacity offset for this species */
12 long int ipOpac,
13 /* lower freq limit to opacity range on energy mesh */
14 long int ipLowLim,
15 /* upper limit to opacity range on energy mesh */
16 long int ipUpLim,
17 /* abundance, we bail if zero */
18 realnum abundance,
19 /* either static 's' or volitile 'v' */
20 char chStat )
21{
22 long int i ,
23 ipOffset,
24 limit;
25
26 DEBUG_ENTRY( "OpacityAdd1Subshell()" );
27
28 /* code spends roughly 20% of its time in this loop*/
29
30 ASSERT( chStat == 's' || chStat == 'v' );
31
32 ipOffset = ipOpac - ipLowLim;
33 ASSERT( ipLowLim > 0 );
34 /* >>chng 02 aug 13, negative offset is ok, it is only this plus ipLowLim that matters */
35 /*ASSERT( ipOffset >= 0 );*/
36
37 limit = MIN2(ipUpLim,rfield.nflux);
38
39 /* do nothing if abundance is zero, or if static opacities do not
40 * need to be redone */
41 if( abundance <= 0. || (chStat=='s' && !opac.lgRedoStatic) )
42 {
43 return;
44 }
45
46 /* volative (outer shell, constantly reevaluated) or static opacity? */
47 if( chStat=='v' )
48 {
49 for( i=ipLowLim-1; i < limit; i++ )
50 {
51 opac.opacity_abs[i] += opac.OpacStack[i+ipOffset]*abundance;
52 }
53 }
54 else
55 {
56 for( i=ipLowLim-1; i < limit; i++ )
57 {
58 opac.OpacStatic[i] += opac.OpacStack[i+ipOffset]*abundance;
59 }
60 }
61 return;
62}
63
64/*OpacityAdd1SubshellInduc add opacity of individual species, including stimulated emission */
66 /* pointer to opacity within opacity stack */
67 long int ipOpac,
68 /* pointer to low energy in continuum array for this opacity band */
69 long int ipLowEnergy,
70 /* pointer to high energy in continuum array for this opacity */
71 long int ipHiEnergy,
72 /* this abundance of this species, may be zero */
73 double abundance,
74 /* the departure coef, may be infinite or zero */
75 double DepartCoef ,
76 /* either 'v' for volitile or 's' for static opacities */
77 char chStat )
78{
79 long int i,
80 iup,
81 k;
82
83 DEBUG_ENTRY( "OpacityAdd1SubshellInduc()" );
84
85 /* add opacity of individual species, including stimulated emission
86 * abundance is the density of the lower level (cm^-3)
87 * DepartCoef is its departure coefficient, can be zero */
88
89 /* this is opacity offset, must be positive */
90 ASSERT( ipOpac > 0 );
91
92 /* check that chStat is either 'v' or 's' */
93 ASSERT( chStat == 'v' || chStat == 's' );
94
95 /* do nothing if abundance is zero, or if static opacities do not
96 * need to be redone */
97 if( abundance <= 0. || (chStat=='s' && !opac.lgRedoStatic) )
98 {
99 return;
100 }
101
102 k = ipOpac - ipLowEnergy;
103
104 /* DepartCoef is dep coef, rfield.lgInducProcess is turned off with 'no indcued' command */
105 if( (DepartCoef > 1e-35 && rfield.lgInducProcess) && hydro.lgHInducImp )
106 {
107 iup = MIN2(ipHiEnergy,rfield.nflux);
108 /* >>>chng 99 apr 29, following was present, caused pdr to make opac at impossible energy*/
109 /*iup = MAX2(ipLowEnergy,iup);*/
110 double DepartCoefInv = 1./DepartCoef;
111 if( chStat == 'v' )
112 {
113 /* volitile opacities, always reevaluate */
114 for( i=ipLowEnergy-1; i < iup; i++ )
115 {
116 opac.opacity_abs[i] += opac.OpacStack[i+k]*abundance*
117 MAX2(0. , 1.- rfield.ContBoltz[i]*DepartCoefInv);
118 }
119 }
120 else
121 {
122 /* static opacities, save in special array */
123 for( i=ipLowEnergy-1; i < iup; i++ )
124 {
125 opac.OpacStatic[i] += opac.OpacStack[i+k]*abundance*
126 MAX2(0. , 1.- rfield.ContBoltz[i]*DepartCoefInv);
127 }
128 }
129 }
130
131 else
132 {
133 /* DepartCoef is the departure coef, which can go to zero,
134 * neglect stimulated emission in this case */
135 iup = MIN2(ipHiEnergy,rfield.nflux);
136 /* >>>chng 99 apr 29, following was present, caused pdr to make opac at impossible energy*/
137 /*iup = MAX2(ipLowEnergy,iup);*/
138 if( chStat == 'v' )
139 {
140 for( i=ipLowEnergy-1; i < iup; i++ )
141 {
142 opac.opacity_abs[i] += opac.OpacStack[i+k]*abundance;
143 }
144 }
145 else
146 {
147 for( i=ipLowEnergy-1; i < iup; i++ )
148 {
149 opac.OpacStatic[i] += opac.OpacStack[i+k]*abundance;
150 }
151 }
152 }
153
154 return;
155}
#define ASSERT(exp)
Definition cddefines.h:578
#define MIN2
Definition cddefines.h:761
float realnum
Definition cddefines.h:103
#define MAX2
Definition cddefines.h:782
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
t_hydro hydro
Definition hydrogenic.cpp:5
t_opac opac
Definition opacity.cpp:5
void OpacityAdd1SubshellInduc(long int ipOpac, long int ipLowEnergy, long int ipHiEnergy, double abundance, double DepartCoef, char chStat)
void OpacityAdd1Subshell(long int ipOpac, long int ipLowLim, long int ipUpLim, realnum abundance, char chStat)
t_rfield rfield
Definition rfield.cpp:8