cloudy trunk
Loading...
Searching...
No Matches
cool_calc.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/*CoolCalc compute calcium cooling */
4#include "cddefines.h"
5#include "taulines.h"
6#include "doppvel.h"
7#include "phycon.h"
8#include "ca.h"
9#include "dense.h"
10#include "thermal.h"
11#include "opacity.h"
12#include "rfield.h"
13#include "ligbar.h"
14#include "lines_service.h"
15#include "atoms.h"
16#include "cooling.h"
17#include "iso.h"
18
19void CoolCalc(void)
20{
21 double a21,
22 a31,
23 a41,
24 a42,
25 a51,
26 a52,
27 a53,
28 c21,
29 Ca2pop[5] ,
30 cs,
31 cs14,
32 cs15,
33 d41,
34 d42,
35 d51,
36 d52,
37 d53,
38 hlgam,
39 op41,
40 op51,
41 opckh,
42 opcxyz,
43 PhotoRate2,
44 PhotoRate3,
45 PhotoRate4,
46 PhotoRate5,
47 r21,
48 r31,
49 r41,
50 r42,
51 r51,
52 r52,
53 r53;
54 static double gCa2[5]={2.,4.,6.,2.,4.};
55 static double exCa2[4]={13650.2,60.7,11480.6,222.9};
56 static realnum opCax = 0.f;
57 static realnum opCay = 0.f;
58 static realnum opCaz = 0.f;
59
60 DEBUG_ENTRY( "CoolCalc()" );
61
62 /* Ca I 4228 */
65
66 /* photoionization of evcited levels by Ly-alpha */
67 hlgam = rfield.otslin[ iso_sp[ipH_LIKE][ipHYDROGEN].trans(ipH2p,ipH1s).ipCont() -1];
68 PhotoRate5 = 1.7e-18*hlgam;
69 PhotoRate4 = 8.4e-19*hlgam;
70 PhotoRate3 = 7.0e-18*hlgam;
71 PhotoRate2 = 4.8e-18*hlgam;
72
73 /* spontaneous decays
74 * frist two trans prob from
75 * >>refer Ca2 AS Zeippen, C.J. 1990, A&A, 229, 248 */
76 a21 = 1.02*TauLines[ipT7324].Emis().Pesc();
77 a31 = 1.05*TauLines[ipT7291].Emis().Pesc();
78 a41 = 1.4e8*TauLines[ipT3969].Emis().Pesc();
79 a51 = 1.4e8*TauLines[ipT3934].Emis().Pesc();
80 a42 = 7.9e6*TauLines[ipT8662].Emis().Pesc();
81 a52 = 8.2e5*TauLines[ipT8498].Emis().Pesc();
82 a53 = 7.48e6*TauLines[ipT8542].Emis().Pesc();
83
84 /* destruction of IR triplet by continuous opacities */
85 opcxyz = opac.opacity_abs[ TauLines[ipT7324].ipCont() -1];
86
87 /* opcxyz = opac(icaxyz) */
88 if( opcxyz > 0. )
89 {
90 d52 = 5.6*opcxyz/(opcxyz + opCax)*(1. - TauLines[ipT8498].Emis().Pesc());
91 d53 = 5.6*opcxyz/(opcxyz + opCay)*(1. - TauLines[ipT8542].Emis().Pesc());
92 d42 = 5.6*opcxyz/(opcxyz + opCaz)*(1. - TauLines[ipT8662].Emis().Pesc());
93 }
94 else
95 {
96 d52 = 0.;
97 d53 = 0.;
98 d42 = 0.;
99 }
100
101 /* near UV dest of KH by background continuum */
102 opckh = opac.opacity_abs[ TauLines[ipT3969].ipCont() -1];
103
104 /* opckh = opac(icakh) */
105 if( opckh > 0. )
106 {
107 op51 = dense.xIonDense[ipCALCIUM][1]*3.89e-7/GetDopplerWidth(dense.AtomicWeight[ipCALCIUM]);
108 d51 = 5.6*opckh/(opckh + op51);
109 op41 = dense.xIonDense[ipCALCIUM][1]*1.96e-7/GetDopplerWidth(dense.AtomicWeight[ipCALCIUM]);
110 d41 = 5.6*opckh/(opckh + op41);
111 }
112 else
113 {
114 op51 = 0.;
115 d51 = 0.;
116 op41 = 0.;
117 d41 = 0.;
118 }
119 /* net rates */
120 r21 = PhotoRate2 + a21;
121 r31 = PhotoRate3 + a31;
122 r41 = a41 + PhotoRate4 + d41;
123 r51 = a51 + PhotoRate5 + d51;
124 r42 = a42 + d42;
125 r52 = a52 + d52;
126 r53 = a53 + d53;
127 cs14 = 0.923*phycon.te10*phycon.te10;
128 cs15 = cs14*2.;
129 TauLines[ipT3969].Coll().col_str() = (realnum)cs14;
130 TauLines[ipT3934].Coll().col_str() = (realnum)cs15;
131
132 /* following used to correct rec contribution
133 * fcakh = a51 / ( a51 + eden*1.5e-5 / sqrte )
134 * cs 1-2 from
135 * >>refer Ca2 CS Saraph, H.E. 1970, J. Phys. B, 3, 952
136 * other
137 * >>refer Ca2 CS Chidichimo, M.C. 1981, J. Phys. B, 14, 4149 */
138 double Cooling , CoolingDeriv;
139 atom_pop5(gCa2,exCa2,5.8,8.6,cs14,cs15,20.6,22.9,9.8,3.4,44.4,1.0,
140 r21,r31,r41,r51,0.,r42,r52,0.,r53,0.,Ca2pop,
141 dense.xIonDense[ipCALCIUM][1],&Cooling , &CoolingDeriv, 0.,0.,0.,0.);
142
143 /* CDSQTE = 8.629E-6*EDEN/SQRTE */
144 c21 = 5.8/4.*dense.cdsqte;
145
146 /* remember largest ratio of Ly-al removal to total */
147 if( dense.xIonDense[ipCALCIUM][1] > 0. )
148 ca.Ca2RmLya = MAX2(ca.Ca2RmLya,(realnum)(PhotoRate2/(PhotoRate2+a21+c21)));
149
150 ca.Cak = (realnum)(Ca2pop[4]*a51*5.06e-12);
151 ca.Cah = (realnum)(Ca2pop[3]*a41*5.01e-12);
152 ca.Cax = (realnum)(Ca2pop[4]*a52*2.34e-12);
153 ca.Cay = (realnum)(Ca2pop[4]*a53*2.33e-12);
154 ca.Caz = (realnum)(Ca2pop[3]*a42*2.30e-12);
155 ca.Caf1 = (realnum)(Ca2pop[2]*a31*2.73e-12);
156 ca.Caf2 = (realnum)(Ca2pop[1]*a21*2.72e-12);
157 ca.popca2ex = (realnum)(Ca2pop[1] + Ca2pop[2] + Ca2pop[3] + Ca2pop[4]);
158
159 /* this is the total cooling due to the model atom */
160 ca.Cair = ca.Cax + ca.Cay + ca.Caz;
161 ca.c7306 = ca.Caf1 + ca.Caf2;
162 ca.Cakh = ca.Cak + ca.Cah;
163
164 // total cooling from 5-level atom
165 CoolAdd("Ca 2",7306,Cooling);
166 thermal.dCooldT += CoolingDeriv;
167
168 /*fprintf(ioQQQ,"DEBUG ca2\t%.2f\t%.5e\t%.4e\t%.4e\n",
169 fnzone, phycon.te,ca.Cakh,dense.xIonDense[ipCALCIUM][1]);*/
170
171 /* level populations that will be used for excited state photoionization */
172 ca.dstCala = (realnum)(Ca2pop[4]*PhotoRate5 + Ca2pop[3]*PhotoRate4);
173 ca.dCakh = (realnum)(ca.dstCala*5.03e-12);
174 ca.dCaf12 = (realnum)((Ca2pop[2]*PhotoRate3 + Ca2pop[1]*PhotoRate2)*2.31e-12);
175 opCax = (realnum)(Ca2pop[1]*1.13e-8/GetDopplerWidth(dense.AtomicWeight[ipCALCIUM]));
176 opCay = (realnum)(Ca2pop[2]*6.87e-8/GetDopplerWidth(dense.AtomicWeight[ipCALCIUM]));
177 opCaz = (realnum)(Ca2pop[1]*5.74e-8/GetDopplerWidth(dense.AtomicWeight[ipCALCIUM]));
178
179 /* total rate Lalpha destroys CaII,
180 * this is only used in ioncali to increase ionization rate by
181 * adding it to the ct vector */
182 if( dense.xIonDense[ipCALCIUM][1] > 0. )
183 {
184 ca.dstCala = (realnum)(
185 (ca.dstCala + ca.dCaf12/2.31e-12)/dense.xIonDense[ipCALCIUM][1]);
186 {
187 /*@-redef@*/
188 enum {DEBUG_LOC=false};
189 /*@+redef@*/
190 if( DEBUG_LOC )
191 {
192 fprintf(ioQQQ," hlgam is %e\n", hlgam);
193 }
194 }
195 }
196 else
197 {
198 ca.dstCala = 0.;
199 }
200 ca.Ca3d = (realnum)(Ca2pop[1] + Ca2pop[2]);
201 ca.Ca4p = (realnum)(Ca2pop[3] + Ca2pop[4]);
202
203 /* incl stimulated emission for Calcium II 5-level atom */
204 TauLines[ipT3934].Emis().PopOpc() = (Ca2pop[0] - Ca2pop[4]/2.);
205 (*TauLines[ipT3934].Hi()).Pop() = Ca2pop[4];
206 (*TauLines[ipT3934].Lo()).Pop() = Ca2pop[0];
207 TauLines[ipT3969].Emis().PopOpc() = (Ca2pop[0] - Ca2pop[3]);
208 (*TauLines[ipT3969].Hi()).Pop() = Ca2pop[3];
209 (*TauLines[ipT3969].Lo()).Pop() = Ca2pop[0];
210
211 TauLines[ipT8498].Emis().PopOpc() = (Ca2pop[1] - Ca2pop[4]);
212 (*TauLines[ipT8498].Hi()).Pop() = Ca2pop[4];
213 (*TauLines[ipT8498].Lo()).Pop() = Ca2pop[1];
214 TauLines[ipT8542].Emis().PopOpc() = (Ca2pop[2] - Ca2pop[4]*1.5);
215 (*TauLines[ipT8542].Hi()).Pop() = Ca2pop[4];
216 (*TauLines[ipT8542].Lo()).Pop() = Ca2pop[2];
217 TauLines[ipT8662].Emis().PopOpc() = (Ca2pop[1] - Ca2pop[3]*2.);
218 (*TauLines[ipT8662].Hi()).Pop() = Ca2pop[3];
219 (*TauLines[ipT8662].Lo()).Pop() = Ca2pop[1];
220 TauLines[ipT7291].Emis().PopOpc() = dense.xIonDense[ipCALCIUM][1];
221 (*TauLines[ipT7291].Hi()).Pop() = 0.;
222 (*TauLines[ipT7291].Lo()).Pop() = dense.xIonDense[ipCALCIUM][1];
223 TauLines[ipT7324].Emis().PopOpc() = dense.xIonDense[ipCALCIUM][1];
224 (*TauLines[ipT7324].Hi()).Pop() = 0.;
225 (*TauLines[ipT7324].Lo()).Pop() = dense.xIonDense[ipCALCIUM][1];
226
227 /* Ca IV 3.2 micron; data from
228 * >>refer Ca4 AS Mendoza, C. 1982, in IAU Symp. 103, Planetary
229 * >>refercon Nebulae, ed. D.R. Flower, (Dordrecht, Holland: D. Reidel Publishing Co.), 143
230 * Y(ik) from
231 * >>refer Ca4 CS Pelan, J., & Berrington, K. A. 1995, A&AS, 110, 209 */
232 if( phycon.te <= 1e5 )
233 {
234 cs = MAX2(1.0,8.854e-3*phycon.sqrte);
235 }
236 else if( phycon.te < 2.512e5 )
237 {
238 cs = 2.8;
239 }
240 else
241 {
242 cs = 641.1/(phycon.te30*phycon.te10*phycon.te02*phycon.te02/
243 phycon.te003);
244 }
245 PutCS(cs,TauLines[ipTCa3]);
247
248 return;
249}
long ipT8498
long ipT8662
long ipTCa3
long ipT8542
long ipT7324
long ipCaI4228
long ipT7291
long ipT3969
long ipT3934
void atom_level2(const TransitionProxy &t)
void atom_pop5(const double g[], const double EnerWN[], double cs12, double cs13, double cs14, double cs15, double cs23, double cs24, double cs25, double cs34, double cs35, double cs45, double a21, double a31, double a41, double a51, double a32, double a42, double a52, double a43, double a53, double a54, double p[], realnum abund, double *Cooling, double *CoolingDeriv, double pump01, double pump02, double pump03, double pump04)
Definition atom_pop5.cpp:13
t_ca ca
Definition ca.cpp:5
FILE * ioQQQ
Definition cddefines.cpp:7
float realnum
Definition cddefines.h:103
#define MAX2
Definition cddefines.h:782
const int ipCALCIUM
Definition cddefines.h:324
const int ipHYDROGEN
Definition cddefines.h:305
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
void CoolCalc(void)
Definition cool_calc.cpp:19
void CoolAdd(const char *chLabel, realnum lambda, double cool)
Definition cool_etc.cpp:13
t_dense dense
Definition dense.cpp:24
realnum GetDopplerWidth(realnum massAMU)
t_iso_sp iso_sp[NISO][LIMELM]
Definition iso.cpp:8
const int ipH1s
Definition iso.h:27
const int ipH2p
Definition iso.h:29
const int ipH_LIKE
Definition iso.h:62
t_opac opac
Definition opacity.cpp:5
t_phycon phycon
Definition phycon.cpp:6
t_rfield rfield
Definition rfield.cpp:8
TransitionList TauLines("TauLines", &AnonStates)
t_thermal thermal
Definition thermal.cpp:5
void MakeCS(const TransitionProxy &t)
void PutCS(double cs, const TransitionProxy &t)