cloudy trunk
Loading...
Searching...
No Matches
parse_dlaw.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/*ParseDLaw parse parameters on the dlaw command */
4#include "cddefines.h"
5#include "dense.h"
6#include "optimize.h"
7#include "abund.h"
8#include "input.h"
9#include "radius.h"
10#include "parser.h"
11
13{
14 bool lgEnd;
15 long int j;
16
17 DEBUG_ENTRY( "ParseDLaw()" );
18
19 if( dense.gas_phase[ipHYDROGEN] > 0. )
20 {
21 fprintf( ioQQQ, " PROBLEM DISASTER More than one density command was entered.\n" );
23 }
24
25 /* call fcn dense_fabden(RADIUS) which uses the ten parameters
26 * N.B.; existing version of dense_fabden must be deleted
27 * >>chng 96 nov 29, added table option */
28 if( p.nMatch("TABL") )
29 {
30 /* when called, read in densities from input stream */
31 strcpy( dense.chDenseLaw, "DLW2" );
32 if( p.nMatch("DEPT") )
33 {
34 dense.lgDLWDepth = true;
35 }
36 else
37 {
38 dense.lgDLWDepth = false;
39 }
40
41 p.getline();
42 dense.frad[0] = (realnum)p.FFmtRead();
43 dense.fhden[0] = (realnum)p.FFmtRead();
44 if( p.lgEOL() )
45 {
46 fprintf( ioQQQ, " No pairs entered - can\'t interpolate.\n Sorry.\n" );
48 }
49
50 dense.nvals = 2;
51 lgEnd = false;
52
53 /* read pairs of numbers until we find line starting with END */
54 /* >>chng 04 jan 27, loop to LIMTABDLAW from LIMTABD, as per
55 * var definitions, caught by Will Henney */
56 while( !lgEnd && dense.nvals < LIMTABDLAW )
57 {
58 p.getline();
59 lgEnd = p.m_lgEOF;
60 if( !lgEnd )
61 {
62 if( p.strcmp("END") == 0 )
63 lgEnd = true;
64 }
65
66 if( !lgEnd )
67 {
68 dense.frad[dense.nvals-1] = (realnum)p.FFmtRead();
69 dense.fhden[dense.nvals-1] = (realnum)p.FFmtRead();
70 dense.nvals += 1;
71 }
72 }
73 --dense.nvals;
74
75 for( long i=1; i < dense.nvals; i++ )
76 {
77 /* the radius values are assumed to be strictly increasing */
78 if( dense.frad[i] <= dense.frad[i-1] )
79 {
80 fprintf( ioQQQ, " Radii must be in increasing order. Sorry.\n" );
82 }
83 }
84 }
85 else if( p.nMatch("WIND") )
86 {
87 strcpy( dense.chDenseLaw, "DLW3" );
88 /* This sets up a steady-state "wind" profile parametrized as in Springmann (1994):
89 *
90 * v(r) = v_star + (v_inf - v_0) * sqrt( Beta1 x + (1-Beta1) x^Beta2 )
91 *
92 * A mass loss rate into 4pi sterradians Mdot then allows the density via continuity:
93 *
94 * n(r) = Mdot / ( 4Pi m_H * mu * r^2 * v(r) )
95 */
96
97 /* The parameters must be specified in this order:
98 *
99 * Mdot, v_inf, Beta2, Beta1, v_0, v_star.
100 *
101 * Only the first three are required. The final three may be omitted right to left and
102 * take default values Beta1 = v_0 = v_star = 0. */
103
104 for( j=0; j < 6; j++ )
105 {
106 dense.DensityLaw[j] = p.FFmtRead();
107 if( j <= 2 && p.lgEOL() )
108 p.NoNumb("density law element");
109 }
110 }
111 else
112 {
113 /* this is usual case, call dense_fabden to get density */
114 for( j=0; j < 10; j++ )
115 {
116 dense.DensityLaw[j] = p.FFmtRead();
117 if( j == 0 && p.lgEOL() )
118 p.NoNumb("density law element");
119 }
120
121 /* set flag so we know which law to use later */
122 strcpy( dense.chDenseLaw, "DLW1" );
123
124 /* vary option */
125 if( optimize.lgVarOn )
126 {
127 ostringstream oss;
128 oss << "DLAW %f" << setprecision(7);
129 for( j=1; j < 10; j++ )
130 oss << " " << dense.DensityLaw[j];
131 strcpy( optimize.chVarFmt[optimize.nparm], oss.str().c_str() );
132 optimize.lgOptimizeAsLinear[optimize.nparm] = true;
133
134 /* index for where to write */
135 optimize.nvfpnt[optimize.nparm] = input.nRead;
136 optimize.vparm[0][optimize.nparm] = (realnum)dense.DensityLaw[0];
137 optimize.vincr[optimize.nparm] = 0.5f;
138 optimize.nvarxt[optimize.nparm] = 1;
139 ++optimize.nparm;
140 }
141 }
142
143 /* set fake density to signal that density command was entered */
144 /* real density will be set once all input commands have been read */
145 /* this is necessary since density may depend on subsequent commands */
146 dense.SetGasPhaseDensity( ipHYDROGEN, 1.f );
147
148 return;
149}
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
const int ipHYDROGEN
Definition cddefines.h:305
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
bool getline(void)
Definition parser.cpp:164
double FFmtRead(void)
Definition parser.cpp:353
bool nMatch(const char *chKey) const
Definition parser.h:135
int strcmp(const char *s2)
Definition parser.h:177
bool lgEOL(void) const
Definition parser.h:98
NORETURN void NoNumb(const char *chDesc) const
Definition parser.cpp:233
bool m_lgEOF
Definition parser.h:42
t_dense dense
Definition dense.cpp:24
const int LIMTABDLAW
Definition dense.h:26
t_input input
Definition input.cpp:12
t_optimize optimize
Definition optimize.cpp:5
void ParseDLaw(Parser &p)