cloudy trunk
Loading...
Searching...
No Matches
parse_coronal.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/*ParseCoronal parse parameters off coronal equilibrium command */
4#include "cddefines.h"
5#include "rfield.h"
6#include "thermal.h"
7#include "input.h"
8#include "optimize.h"
9#include "phycon.h"
10#include "radius.h"
11#include "dynamics.h"
12#include "parser.h"
13#include "atmdat.h"
14
15/*ParseCoronal parse parameters off coronal equilibrium command */
17{
18 double a;
19
20 DEBUG_ENTRY( "ParseCoronal()" );
21
22 /* use coronal command to establish initial conditions in a cooling
23 * time-varying cloud */
24 if( p.nMatch( "INIT" ) && p.nMatch( "TIME" ) )
25 {
26 dynamics.lg_coronal_time_init = true;
27 if( p.nMatch( "TRAC" ) )
28 dynamics.lgTracePrint = true;
29 }
30
31 /* coronal equilibrium; set constant temperature to number on line */
32 thermal.lgTemperatureConstant = true;
33 thermal.lgTemperatureConstantCommandParsed = true;
34
35 // kinetic temperatures for a given ion are higher for coronal equilibrium
36 // simulations - large Fe chianti models are needed to get the full cooling
37 // tests show that this mainly affects the cooling around 1e7 K.
38 // the other elements included in Chianti hardly affect the cooling so
39 // are not changed
40 if( !atmdat.lgChiantiLevelsSet )
41 {
42 atmdat.nChiantiMaxLevelsFe = atmdat.nChiantiCollLevelsFe;
43 atmdat.nChiantiMaxLevels = atmdat.nChiantiCollLevels;
44 }
45
46 a = p.FFmtRead();
47 if( p.lgEOL() )
48 {
49 fprintf( ioQQQ, " There should be a temperature on this line.\n" );
51 }
52
53 /* numbers less than or equal to 10 are the log of the temperature */
54 if( (a <= 10. && !p.nMatch("LINE")) || p.nMatch(" LOG") )
55 {
56 thermal.ConstTemp = (realnum)pow(10.,a);
57 }
58 else
59 {
60 thermal.ConstTemp = (realnum)a;
61 }
62
63 /* check temperature bounds */
64 if( thermal.ConstTemp < phycon.TEMP_LIMIT_LOW )
65 {
66 thermal.ConstTemp = (realnum)(1.0001*phycon.TEMP_LIMIT_LOW);
67 fprintf( ioQQQ, " PROBLEM Te too low, reset to %g K.\n",
68 thermal.ConstTemp );
69 }
70 if( thermal.ConstTemp > phycon.TEMP_LIMIT_HIGH )
71 {
72 thermal.ConstTemp = (realnum)(0.9999*phycon.TEMP_LIMIT_HIGH);
73 fprintf( ioQQQ, " PROBLEM Te too high, reset to %g K.\n",
74 thermal.ConstTemp );
75 }
76
77 /* now simulate a LASER line */
78 strcpy( rfield.chSpType[rfield.nShape], "LASER" );
79 /* scan off the laser's energy ion Rydbergs */
80 rfield.slope[rfield.nShape] = rfield.egamry*0.99;
81 /* default width is 0.05 */
82 rfield.cutoff[rfield.nShape][0] = 0.05;
83
84 /* simulate an ionization parameter line */
85 strcpy( rfield.chRSpec[p.m_nqh], "SQCM" );
86 strcpy( rfield.chSpNorm[p.m_nqh], "IONI" );
87
88 /* >>chng 96 jun 17, to stop mole network from crashing */
89 /* >>chng 05 aug 15, this sets ionization parameter, in test case ism_hot_brems the
90 * value of 1e-10 was enough to dominate the ionization of he-like N - it's ionization
91 * then jumped due to large optical depth in the continuum - change U from -10 to -15 */
92 /* >>chng 05 aug 16, this very strongly affected the coll_t4 sim - apparently there
93 * was a significant photoionization contribution from the -10 continuum,
94 * this was close to a 'no photoionization' case, but lower further to insure
95 * no photo contribution
96 * chang from -15 to -20 */
97 rfield.totpow[p.m_nqh] = -20.f;
98
99 ++rfield.nShape;
100 if( rfield.nShape >= LIMSPC )
101 {
102 /* too many continua were entered */
103 fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
105 }
106 ++p.m_nqh;
107 if( p.m_nqh >= LIMSPC )
108 {
109 /* too many continua were entered */
110 fprintf( ioQQQ, " Too many continua entered; increase LIMSPC\n" );
112 }
113
114 /* set R to large value if U specified but R is not */
115 /* set radius to very large value if not already set */
116 /* >>chng 01 jul 24, from Radius == 0 to this, as per PvH comments */
117 if( !radius.lgRadiusKnown )
118 {
119 radius.Radius = pow(10.,radius.rdfalt);
120 }
121
122 /* vary option */
123 if( optimize.lgVarOn )
124 {
125 /* no luminosity options on vary */
126 optimize.nvarxt[optimize.nparm] = 1;
127 strcpy( optimize.chVarFmt[optimize.nparm], "COROnal equilibrium %f LOG" );
128 if( dynamics.lg_coronal_time_init )
129 strcat( optimize.chVarFmt[optimize.nparm], " TIME INIT" );
130
131 /* pointer to where to write */
132 optimize.nvfpnt[optimize.nparm] = input.nRead;
133
134 /* log of temp will be pointer */
135 optimize.vparm[0][optimize.nparm] = (realnum)log10(thermal.ConstTemp);
136 optimize.vincr[optimize.nparm] = 0.1f;
137 optimize.varang[optimize.nparm][0] = (realnum)log10(1.00001*phycon.TEMP_LIMIT_LOW);
138 optimize.varang[optimize.nparm][1] = (realnum)log10(0.99999*phycon.TEMP_LIMIT_HIGH);
139 ++optimize.nparm;
140 }
141 return;
142}
t_atmdat atmdat
Definition atmdat.cpp:6
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
long int m_nqh
Definition parser.h:41
t_dynamics dynamics
Definition dynamics.cpp:44
t_input input
Definition input.cpp:12
t_optimize optimize
Definition optimize.cpp:5
void ParseCoronal(Parser &p)
t_phycon phycon
Definition phycon.cpp:6
t_radius radius
Definition radius.cpp:5
t_rfield rfield
Definition rfield.cpp:8
const int LIMSPC
Definition rfield.h:18
t_thermal thermal
Definition thermal.cpp:5