cloudy trunk
Loading...
Searching...
No Matches
rt_line_one_tauinc.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/*RT_line_one_tauinc increment optical depths for all heavy element lines, zone by zone,
4 * mainly called by RT_tau_inc, but also by FeII */
5#include "cddefines.h"
6#include "doppvel.h"
7#include "geometry.h"
8#include "rfield.h"
9#include "radius.h"
10#include "wind.h"
11#include "rt.h"
12#include "physconst.h"
13#include "cosmology.h"
14#include "transition.h"
15
17 /* following four are flags to generate info if some species has extreme maser */
18 long int maser_flag_species,
19 long int maser_flag_ion,
20 long int maser_flag_hi,
21 long int maser_flag_lo,
22 realnum DopplerWidth )
23
24{
25 DEBUG_ENTRY( "RT_line_one_tauinc()" );
26
27 /* routine increments optical depths for static or expanding atmosphere */
28
29 /* this is line center frequency, including bulk motion of gas */
30 long int ipLineCenter = t.Emis().ipFine() + rfield.ipFineConVelShift;
31 double OpacityEffective, EffectiveThickness;
32 realnum dTau_total;
33
34 /* find line center opacity - use fine opacity if array indices are OK */
35 if( t.Emis().ipFine()>=0 && ipLineCenter>0 && ipLineCenter<rfield.nfine && rfield.lgOpacityFine )
36 {
37 /* use fine opacities fine grid fine mesh to get optical depth
38 * to continuum source */
39 /* total line center optical depth, all overlapping lines included */
40 OpacityEffective = rfield.fine_opac_zone[ipLineCenter];
41 }
42 else
43 {
44 OpacityEffective = t.Emis().PopOpc() * t.Emis().opacity() / DopplerWidth;
45 }
46
47#if 0
48 if( rfield.anu[ t.ipCont-1 ] < rfield.plsfrq )
49 {
50 /* transition is below plasma frequency - make optical depth huge */
51 dTau_total = 1.e10;
52
53 t.Emis().TauIn() = dTau_total;
54 t.Emis().TauCon() = dTau_total;
55 t.Emis().TauTot() = dTau_total;
56 }
57 else
58#endif
59 if( cosmology.lgDo )
60 {
61 /* dv/dr (s-1), equal to dv/dt / v */
62 /* in this case, dv/dr is just the Hubble factor */
63 wind.dvdr = GetHubbleFactor(cosmology.redshift_current);
64
65 fixit(); //This doppler width is sqrt(2kt/m), but Seager et al use sqrt(3kt/m). Resolve
66 EffectiveThickness = DopplerWidth / wind.dvdr;
67 dTau_total = (realnum)(OpacityEffective * EffectiveThickness);
68
69 t.Emis().TauIn() = dTau_total;
70 t.Emis().TauCon() = dTau_total;
71 t.Emis().TauTot() = dTau_total;
72 }
73
74 /* use cumulated fine optical depth for both d-critical and static,
75 * for d-critical speeds are only roughly sonic
76 * optical depth is computed including velocity shift */
77 else if( ! wind.lgBallistic() )
78 {
79 /* static and negative velocity solutions */
80 EffectiveThickness = radius.drad_x_fillfac;
81 dTau_total = (realnum)(OpacityEffective * EffectiveThickness);
82
83 t.Emis().TauIn() += dTau_total;
84 t.Emis().TauCon() += dTau_total;
85
86 }
87
88 else
89 {
90 /* ballistic outflowing wind
91 * effective length scale for Sobolev or LVG approximation, eqn 3 of
92 * >>refer RT wind Castor, J.I., Abbott, D.C., & Klein, R.I., 1975, ApJ, 195, 157
93 */
94
95 /* dv/dr (s-1), equal to dv/dt / v */
96 wind.dvdr = fabs(wind.AccelTotalOutward - wind.AccelGravity) / wind.windv;
97 /* depth (cm) over which wind accelerates by one velocity width
98 * include filling factor */
99 EffectiveThickness = DopplerWidth / SDIV(wind.dvdr) * geometry.FillFac;
100
101 /* min2 is to not let the physical scale exceed the current depth */
102 EffectiveThickness = MIN2( radius.depth, EffectiveThickness );
103 dTau_total = (realnum)(OpacityEffective * EffectiveThickness);
104
105 t.Emis().TauIn() = dTau_total;
106 t.Emis().TauCon() = dTau_total;
107 t.Emis().TauTot() = dTau_total;
108 }
109
110 /* keep track of any masers */
111 if( dTau_total < rt.dTauMase )
112 {
113 rt.dTauMase = dTau_total;
114 rt.mas_species = maser_flag_species;
115 rt.mas_ion = maser_flag_ion;
116 rt.mas_hi = maser_flag_hi;
117 rt.mas_lo = maser_flag_lo;
118 if( rt.dTauMase < -1. )
119 rt.lgMaserCapHit = true;
120 }
121
122 return;
123}
#define MIN2
Definition cddefines.h:761
float realnum
Definition cddefines.h:103
sys_float SDIV(sys_float x)
Definition cddefines.h:952
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
void fixit(void)
Definition service.cpp:991
long int & ipFine() const
Definition emission.h:413
double & PopOpc() const
Definition emission.h:603
realnum & TauCon() const
Definition emission.h:453
realnum & TauIn() const
Definition emission.h:423
realnum & opacity() const
Definition emission.h:593
realnum & TauTot() const
Definition emission.h:433
long & ipCont() const
Definition transition.h:450
EmissionList::reference Emis() const
Definition transition.h:408
t_cosmology cosmology
Definition cosmology.cpp:11
realnum GetHubbleFactor(realnum z)
Definition cosmology.cpp:13
t_geometry geometry
Definition geometry.cpp:5
t_radius radius
Definition radius.cpp:5
t_rfield rfield
Definition rfield.cpp:8
t_rt rt
Definition rt.cpp:5
void RT_line_one_tauinc(const TransitionProxy &t, long int maser_flag_species, long int maser_flag_ion, long int maser_flag_hi, long int maser_flag_lo, realnum DopplerWidth)
Wind wind
Definition wind.cpp:5