cloudy
trunk
Main Page
Related Pages
Data Structures
Files
File List
Globals
All
Data Structures
Files
Functions
Variables
Typedefs
Enumerations
Enumerator
Friends
Macros
Pages
source
iso_continuum_lower.cpp
Go to the documentation of this file.
1
/* This file is part of Cloudy and is copyright (C)1978-2008 by Gary J. Ferland and
2
* others. For conditions of distribution and use see copyright notice in license.txt */
3
/*iso_continuum_lower - limit max prin. quan. no. due to continuum lowering processes */
4
#include "
cddefines.h
"
5
#include "
phycon.h
"
6
#include "
dense.h
"
7
#include "
iso.h
"
8
#include "
hydrogenic.h
"
9
#include "
trace.h
"
10
11
void
iso_continuum_lower
(
long
ipISO,
long
nelem )
12
{
13
double
a;
14
long
int
np, nd, ns, nc;
15
long
eff_charge;
16
17
/* size of rate matrices will be defined according to the n calculated here */
18
19
ASSERT
( nelem <
LIMELM
);
20
/* this may change at a future date. */
21
ASSERT
( ipISO <= 1 );
22
23
eff_charge = nelem + 1 - ipISO;
24
25
/* Particle packing - the density here should be density of all nuclei in the plasma */
26
/* This one is just nuclear charge, which is independent of iso, and always nelem+1. */
27
a = sqrt( 1.8887E8 * (nelem+1.) / pow((
double
)
dense
.
xNucleiTotal
, 0.333) );
28
ASSERT
( a > 0. );
29
ASSERT
( a < 1.e15 );
30
if
( a > (
double
)
iso
.
n_HighestResolved_max
[ipISO][nelem]+(
double
)
iso
.
nCollapsed_max
[ipISO][nelem] )
31
{
32
np =
iso
.
n_HighestResolved_max
[ipISO][nelem]+
iso
.
nCollapsed_max
[ipISO][nelem] + 1;
33
}
34
else
35
np = (long)a;
36
37
/* Debye shielding - the density here is electron density */
38
/* This one depends on effective charge. */
39
a = 2.6E7 * eff_charge * eff_charge * pow(
phycon
.
te
/
dense
.
eden
, 0.25);
40
ASSERT
( a > 0. );
41
ASSERT
( a < 1.e15 );
42
if
( a > (
double
)
iso
.
n_HighestResolved_max
[ipISO][nelem]+(
double
)
iso
.
nCollapsed_max
[ipISO][nelem] )
43
{
44
nd =
iso
.
n_HighestResolved_max
[ipISO][nelem]+
iso
.
nCollapsed_max
[ipISO][nelem] + 1;
45
}
46
else
47
nd = (long)a;
48
49
/* Stark broadening - this should be the density of singly charged ions,
50
* both positive and negative. The sum of protons, electrons, and HeII should be
51
* good enough. */
52
/* This one depends on effective charge. */
53
a = 3171. * pow( (
double
)eff_charge, 0.8 ) * pow(
dense
.
eden
+ (
double
)
dense
.
xIonDense
[
ipHYDROGEN
][1]
54
+ (
double
)
dense
.
xIonDense
[
ipHELIUM
][1], -0.1333);
55
ASSERT
( a > 0. );
56
ASSERT
( a < 1.e15 );
57
if
( a > (
double
)
iso
.
n_HighestResolved_max
[ipISO][nelem]+(
double
)
iso
.
nCollapsed_max
[ipISO][nelem] )
58
{
59
ns =
iso
.
n_HighestResolved_max
[ipISO][nelem]+
iso
.
nCollapsed_max
[ipISO][nelem] + 1;
60
}
61
else
62
ns = (long)a;
63
64
ASSERT
( np > 3 );
65
ASSERT
( nd > 3 );
66
ASSERT
( ns > 3 );
67
68
nc =
MIN3
(np, nd, ns);
69
70
/* I assert greater than three because the code depends upon having at least up to n=3, and
71
* because it would take EXTREMELY dense conditions to lower the continuum that much, and
72
* the code would probably get a wrong answer then anyway. */
73
ASSERT
( nc > 3 );
74
75
if
( nc <
iso
.
n_HighestResolved_max
[ipISO][nelem])
76
{
77
iso
.
lgLevelsLowered
[ipISO][nelem] =
true
;
78
iso
.
lgLevelsEverLowered
[ipISO][nelem] =
true
;
79
hydro
.
lgReevalRecom
=
true
;
80
iso
.
n_HighestResolved_local
[ipISO][nelem] = nc;
81
iso
.
nCollapsed_local
[ipISO][nelem] = 0;
82
iso
.
numLevels_local
[ipISO][nelem] =
iso_get_total_num_levels
( ipISO, nc, 0 );
83
}
84
/* Here is the case where the critical n lies among the one or more collapsed levels */
85
/* we just get rid of any that are too high. */
86
else
if
( nc <=
iso
.
n_HighestResolved_max
[ipISO][nelem] +
iso
.
nCollapsed_max
[ipISO][nelem] )
87
{
88
iso
.
lgLevelsLowered
[ipISO][nelem] =
true
;
89
iso
.
lgLevelsEverLowered
[ipISO][nelem] =
true
;
90
hydro
.
lgReevalRecom
=
true
;
91
iso
.
n_HighestResolved_local
[ipISO][nelem] =
iso
.
n_HighestResolved_max
[ipISO][nelem];
92
iso
.
nCollapsed_local
[ipISO][nelem] = nc -
iso
.
n_HighestResolved_local
[ipISO][nelem];
93
iso
.
numLevels_local
[ipISO][nelem] =
94
iso_get_total_num_levels
( ipISO,
iso
.
n_HighestResolved_max
[ipISO][nelem],
iso
.
nCollapsed_local
[ipISO][nelem] );
95
}
96
/* This is usually where control will flow, because in most conditions the continuum will not be lowered.
97
* Nothing changes in this case. */
98
else
99
{
100
iso
.
numLevels_local
[ipISO][nelem] =
iso
.
numLevels_max
[ipISO][nelem];
101
iso
.
nCollapsed_local
[ipISO][nelem] =
iso
.
nCollapsed_max
[ipISO][nelem];
102
iso
.
n_HighestResolved_local
[ipISO][nelem] =
iso
.
n_HighestResolved_max
[ipISO][nelem];
103
104
/* if levels were lowered on last pass but are not now, must reeval */
105
if
(
iso
.
lgLevelsLowered
[ipISO][nelem] )
106
{
107
hydro
.
lgReevalRecom
=
true
;
108
}
109
else
110
{
111
hydro
.
lgReevalRecom
=
false
;
112
}
113
114
iso
.
lgLevelsLowered
[ipISO][nelem] =
false
;
115
}
116
117
/* None of these can be greater than that which was originally malloc'd. */
118
ASSERT
(
iso
.
numLevels_local
[ipISO][nelem] <=
iso
.
numLevels_max
[ipISO][nelem] );
119
ASSERT
(
iso
.
nCollapsed_local
[ipISO][nelem] <=
iso
.
nCollapsed_max
[ipISO][nelem] );
120
ASSERT
(
iso
.
n_HighestResolved_local
[ipISO][nelem] <=
iso
.
n_HighestResolved_max
[ipISO][nelem] );
121
122
/* don't print levels in continuum. */
123
iso
.
numPrintLevels
[ipISO][nelem] =
MIN2
(
iso
.
numPrintLevels
[ipISO][nelem],
iso
.
numLevels_local
[ipISO][nelem] );
124
125
/* Lyman lines can not be greater than original malloc or critical pqn. */
126
iso
.
nLyman
[ipISO] =
MIN2
( nc,
iso
.
nLyman
[ipISO]);
127
128
if
(
trace
.
lgTrace
)
129
{
130
fprintf(
ioQQQ
,
" iso_continuum_lower: ipISO %li nelem %li nc %li numLevels %li nCollapsed %li n_HighestResolved %li \n"
,
131
ipISO,
132
nelem,
133
nc,
134
iso
.
numLevels_local
[ipISO][nelem],
135
iso
.
nCollapsed_local
[ipISO][nelem],
136
iso
.
n_HighestResolved_local
[ipISO][nelem]
137
);
138
}
139
140
return
;
141
}
Generated for cloudy by
1.8.1.2