cloudy trunk
Loading...
Searching...
No Matches
parse_constant.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/*ParseConstant parse parameters from the 'constant ...' command */
4#include "cddefines.h"
5#include "physconst.h"
6#include "optimize.h"
7#include "thermal.h"
8#include "dense.h"
9#include "pressure.h"
10#include "phycon.h"
11#include "input.h"
12#include "parse.h"
13#include "parser.h"
14
16{
17 DEBUG_ENTRY( "ParseConstant()" );
18
19 if( p.nMatch("GRAI") && p.nMatch("TEMP") )
20 {
21 /* constant grain temperature command */
22 thermal.ConstGrainTemp = (realnum)p.FFmtRead();
23
24 /* if linear option is not on the line, convert to exponent if <= 10 */
25 if( !p.nMatch("LINE") )
26 {
27 if( thermal.ConstGrainTemp <= 10. )
28 thermal.ConstGrainTemp = (realnum)pow((realnum)10.f,thermal.ConstGrainTemp);
29 }
30
31 if( p.lgEOL() )
32 {
33 p.NoNumb("grain temperature");
34 }
35 }
36
37 else if( p.nMatch("TEMP") )
38 {
39 /* a constant temperature model */
40 thermal.lgTemperatureConstant = true;
41 thermal.lgTemperatureConstantCommandParsed = true;
42
43 /* this is an option to specify the temperature in different units
44 * keV, eV now supported */
45 realnum convert_to_Kelvin = 1;
46 if( p.nMatch(" EV ") )
47 {
48 convert_to_Kelvin = (realnum)EVDEGK;
49 }
50 else if( p.nMatch(" KEV") )
51 {
52 convert_to_Kelvin = (realnum)(EVDEGK * 1000.);
53 }
54
55 /* this is the "force" temperature. same var used in force temp
56 * command, but lgTSetOn is not set so then allowed to vary
57 * so constant temperature requires both lgTSetOn true and ConstTemp > 0 */
58 thermal.ConstTemp = (realnum)p.FFmtRead();
59 if( p.lgEOL() )
60 p.NoNumb("temperature");
61
62 /* if linear option is not on the line and T<=10, assume number is log */
63 if( p.nMatch(" LOG") || (thermal.ConstTemp <= 10. && !p.nMatch("LINE")) )
64 {
65 if( thermal.ConstTemp > log10(BIGFLOAT) )
66 {
67 fprintf(ioQQQ," PROBLEM temperature entered as a log but is too large "\
68 "for this processor. I am interpreting it as the linear temperature.\n");
69 }
70 else
71 thermal.ConstTemp = (realnum)pow((realnum)10.f,thermal.ConstTemp);
72 }
73 /* do units conversion here */
74 thermal.ConstTemp *= convert_to_Kelvin;
75
76 /* check temperature bounds */
77 if( thermal.ConstTemp < phycon.TEMP_LIMIT_LOW )
78 {
79 thermal.ConstTemp = (realnum)(1.0001*phycon.TEMP_LIMIT_LOW);
80 fprintf( ioQQQ, " PROBLEM Te too low, reset to %g K.\n",
81 thermal.ConstTemp );
82 }
83 if( thermal.ConstTemp > phycon.TEMP_LIMIT_HIGH )
84 {
85 thermal.ConstTemp = (realnum)(0.9999*phycon.TEMP_LIMIT_HIGH);
86 fprintf( ioQQQ, " PROBLEM Te too high, reset to %g K.\n",
87 thermal.ConstTemp );
88 }
89
90 /* set the real electron temperature to the forced value */
91 TempChange( thermal.ConstTemp, false );
92
93 /* vary option */
94 if( optimize.lgVarOn )
95 {
96 /* no luminosity options on vary */
97 optimize.nvarxt[optimize.nparm] = 1;
98 // the keyword LOG is not used above, but is checked elsewhere
99 strcpy( optimize.chVarFmt[optimize.nparm], "CONSTANT TEMP %f LOG" );
100
101 /* pointer to where to write */
102 optimize.nvfpnt[optimize.nparm] = input.nRead;
103
104 /* log of temp will be pointer */
105 optimize.vparm[0][optimize.nparm] = (realnum)log10(thermal.ConstTemp);
106 optimize.vincr[optimize.nparm] = 0.1f;
107 optimize.varang[optimize.nparm][0] = (realnum)log10(1.00001*phycon.TEMP_LIMIT_LOW);
108 optimize.varang[optimize.nparm][1] = (realnum)log10(0.99999*phycon.TEMP_LIMIT_HIGH);
109 ++optimize.nparm;
110 }
111 }
112
113 else if( p.nMatch("DENS") )
114 {
115 /* constant density */
116 strcpy( dense.chDenseLaw, "CDEN" );
117 /* turn off radiation pressure */
118 pressure.lgPres_radiation_ON = false;
119 pressure.lgPres_magnetic_ON = false;
120 pressure.lgPres_ram_ON = false;
121 }
122
123 else if( p.nMatch("PRES") )
124 {
125 /* constant pressure */
126 strcpy( dense.chDenseLaw, "CPRE" );
127
128 /* >>chng 06 jun 20, add reset option, to reset density to keep
129 * initial pressure itself constant from iteration to iteration,
130 * rather than initial density */
131 if( p.nMatch("RESE") )
132 {
133 /* this says not to keep initial density constant,
134 * reset it to keep pressure const */
135 dense.lgDenseInitConstant = false;
136 }
137 else
138 {
139 /* this is default, says keep initial density constant,
140 * so pressure from iter to iter not really const */
141 dense.lgDenseInitConstant = true;
142 }
143
144 if( p.nMatch("TIME") )
145 {
146 // pressure varies as a function of time
147 /* this says not to keep initial density constant,
148 * reset it to keep pressure const */
149 dense.lgDenseInitConstant = false;
150 dense.lgPressureVaryTime = true;
151
152 // required number is timescale for time variation
153 dense.PressureVaryTimeTimescale = (realnum)p.FFmtRead();
154 if( dense.PressureVaryTimeTimescale<SMALLFLOAT )
155 {
156 // need two numbers
157 fprintf(ioQQQ," PROBLEM the constant pressure time command requires"
158 " a positive timescale.\n");
160 }
161 // optional number is index for time variation
162 dense.PressureVaryTimeIndex = (realnum)p.FFmtRead();
163 // pressure will be initial pressure - pressure.PresTotlInit - multiplied by
164 // (time / time scale ) ^ index
165 if( p.lgEOL() )
166 {
167 // need two numbers
168 fprintf(ioQQQ," PROBLEM the constant pressure time command requires"
169 " two numbers, the timescale for the variation and an index.\n");
171 }
172 }
173
174 if( p.nMatch(" GAS") )
175 {
176 /* constant gas pressure (no radiation)
177 * turn off radiation pressure */
178 pressure.lgPres_radiation_ON = false;
179
180 /* turn off incident continuum */
181 pressure.lgContRadPresOn = false;
182
183 /* turn off magnetic and ram pressure */
184 pressure.lgPres_magnetic_ON = false;
185 pressure.lgPres_ram_ON = false;
186
187 /* optional number is power law index */
188 pressure.PresPowerlaw = (realnum)p.FFmtRead();
189 }
190
191 else
192 {
193 /* constant total pressure, gas+rad+incident continuum
194 * turn on radiation pressure */
195 pressure.lgPres_radiation_ON = true;
196 pressure.lgPres_magnetic_ON = true;
197 pressure.lgPres_ram_ON = true;
198
199 /* option to turn off continuum pressure */
200 if( p.nMatch("NO CO") )
201 {
202 pressure.lgContRadPresOn = false;
203 }
204 else
205 {
206 pressure.lgContRadPresOn = true;
207 }
208
209 /* option to not abort when too much radiation pressure, no abort */
210 if( p.nMatch("NO AB") )
211 {
212 pressure.lgRadPresAbortOK = false;
213 }
214 else
215 {
216 pressure.lgRadPresAbortOK = true;
217 }
218 /* there is no optional power law option for constant total pressure */
219 pressure.PresPowerlaw = 0.;
220
221 /* option to set pressure */
222 if( p.nMatch(" SET") )
223 {
224 /* number on line is log of nT - option to specify initial pressure */
225 pressure.lgPressureInitialSpecified = true;
226 /* this is log of nT product - if not present then set zero */
227 pressure.PressureInitialSpecified = p.FFmtRead();
228 if( p.lgEOL() )
229 p.NoNumb("initial pressure" );
230 else
231 /* pressure in nkT units */
232 pressure.PressureInitialSpecified = pow(10.,pressure.PressureInitialSpecified)*
233 BOLTZMANN;
234 }
235 else
236 pressure.lgPressureInitialSpecified = false;
237 }
238
239 /* vary option */
240 if( optimize.lgVarOn && pressure.lgPressureInitialSpecified )
241 {
242 /* no options on vary */
243 optimize.nvarxt[optimize.nparm] = 1;
244 // the keyword LOG is not used above, but is checked elsewhere
245 strcpy( optimize.chVarFmt[optimize.nparm], "CONSTANT PRESSURE SET %f LOG" );
246 if( !dense.lgDenseInitConstant )
247 strcat( optimize.chVarFmt[optimize.nparm], " RESET" );
248 if( !pressure.lgContRadPresOn )
249 strcat( optimize.chVarFmt[optimize.nparm], " NO CONTINUUM" );
250 if( !pressure.lgRadPresAbortOK )
251 strcat( optimize.chVarFmt[optimize.nparm], " NO ABORT" );
252
253 /* pointer to where to write */
254 optimize.nvfpnt[optimize.nparm] = input.nRead;
255
256 /* log of temp will be pointer */
257 optimize.vparm[0][optimize.nparm] = (realnum)log10(pressure.PressureInitialSpecified/BOLTZMANN);
258 optimize.vincr[optimize.nparm] = 0.1f;
259 ++optimize.nparm;
260 }
261 }
262
263 else
264 {
265 /* no keys were recognized */
266 fprintf( ioQQQ, " The keyword should be TEMPerature, DENSity, GAS or PRESsure, sorry.\n" );
268 }
269 return;
270}
FILE * ioQQQ
Definition cddefines.cpp:7
#define EXIT_FAILURE
Definition cddefines.h:140
#define cdEXIT(FAIL)
Definition cddefines.h:434
float realnum
Definition cddefines.h:103
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
double FFmtRead(void)
Definition parser.cpp:353
bool nMatch(const char *chKey) const
Definition parser.h:135
bool lgEOL(void) const
Definition parser.h:98
NORETURN void NoNumb(const char *chDesc) const
Definition parser.cpp:233
UNUSED const realnum BIGFLOAT
Definition cpu.h:189
const realnum SMALLFLOAT
Definition cpu.h:191
t_dense dense
Definition dense.cpp:24
t_input input
Definition input.cpp:12
t_optimize optimize
Definition optimize.cpp:5
void ParseConstant(Parser &p)
t_phycon phycon
Definition phycon.cpp:6
UNUSED const double BOLTZMANN
Definition physconst.h:97
UNUSED const double EVDEGK
Definition physconst.h:186
t_pressure pressure
Definition pressure.cpp:5
void TempChange(double TempNew, bool lgForceUpdate)
t_thermal thermal
Definition thermal.cpp:5