cloudy trunk
Loading...
Searching...
No Matches
pressure_change.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/*PressureChange called by ConvPresTempEdenIoniz
4 * evaluate the current pressure, density change needed to get it to converge,
5 * applies this correction factor to all gas constituents,
6 * sets conv.lgConvPres true if good pressure, false if pressure change capped */
7/*lgConvPres finds amount by which density should change to move towards pressure equilibrium
8 * returns true if pressure is converged */
9#include "cddefines.h"
10
11#include "pressure_change.h"
12
13#include "colden.h"
14#include "conv.h"
15#include "cosmology.h"
16#include "dark_matter.h"
17#include "dense.h"
18#include "dynamics.h"
19#include "geometry.h"
20#include "mole.h"
21#include "phycon.h"
22#include "pressure.h"
23#include "radius.h"
24#include "struc.h"
25#include "thermal.h"
26#include "trace.h"
27#include "wind.h"
28
29/* zoneDensity finds density where prescribed */
31{
32 double new_density;
33
34 DEBUG_ENTRY( "zoneDensity()" );
35
36 pressure.PresTotlError = 0.;
37
38 /* evaluate a series of possible options, and set new_density */
39 /* inside out globule */
40 if( strcmp(dense.chDenseLaw,"GLOB") == 0 )
41 {
42 /* GLBDST is distance from globule, or glbrad-DEPTH */
43 if( radius.glbdst < 0. )
44 {
45 fprintf( ioQQQ, " Globule distance is negative, internal overflow has occured, sorry.\n" );
46 fprintf( ioQQQ, " This is routine lgConvPres, GLBDST is%10.2e\n",
47 radius.glbdst );
49 }
50 new_density =
51 (radius.glbden*pow(radius.glbrad/(radius.glbdst),radius.glbpow))*
52 (scalingDensity()/dense.gas_phase[ipHYDROGEN]);
53 }
54 else if( cosmology.lgDo )
55 {
56 /* cosmological - density varies because of expansion of universe */
57 double dnew = GetDensity( cosmology.redshift_current );
58 new_density = dnew*(scalingDensity()/dense.gas_phase[ipHYDROGEN]);
59 }
60 else if( (strcmp(dense.chDenseLaw,"WIND") == 0) ) {
61
62 if ( wind.lgStatic() )
63 {
64 fprintf( ioQQQ, " PROBLEM WIND called with zero velocity - this is impossible.\n Sorry.\n" );
66 }
67
68 /* this is positive wind velocity the outflowing wind beyond sonic point */
69 else if( wind.lgBallistic() )
70 {
71
72 /* this is logic for supersonic wind solution,
73 * well above sonic point. */
74 if( nzone > 1 )
75 {
76 /* Wind model */
77
78 if( trace.lgTrace && trace.lgWind )
79 {
80 fprintf(ioQQQ," lgConvPres sets AccelGravity %.3e lgDisk?%c\n",
81 wind.AccelGravity ,
82 TorF(wind.lgDisk) );
83 }
84
85 /* following is form of constant acceleration equation, v^2 = u^2 + 2 a s
86 * struc.windv[nzone-2] is velocity of previous zone
87 * this increments that velocity to form square of new wind
88 * velocity for outer edge of this zone */
89 double term = POW2(struc.windv[nzone-2]) + 2.*(wind.AccelTotalOutward - wind.AccelGravity)* radius.drad;
90
91 /* increment velocity if it is substantially positive */
92 fixit(); // RHS of comparison should be ~ sound speed squared
93 if( term <= 1e3 )
94 {
95 /* wind velocity is well below sonic point, give up,
96 * do not change velocity */
97 wind.lgVelPos = false;
98 }
99 else
100 {
101 /* wind.windv is velocity at OUTER edge of this zone */
102 term = sqrt(term);
103 if (wind.windv > 0)
104 {
105 wind.windv = (realnum) term;
106 }
107 else
108 {
109 wind.windv = -(realnum) term;
110 }
111 wind.lgVelPos = true;
112 }
113
114 if( trace.lgTrace && trace.lgWind )
115 {
116 fprintf(ioQQQ," lgConvPres new wind V zn%li %.3e AccelTotalOutward %.3e AccelGravity %.3e\n",
117 nzone,wind.windv, wind.AccelTotalOutward, wind.AccelGravity );
118 }
119 }
120
121 /* conservation of mass sets density here */
122 new_density = wind.emdot/(wind.windv*radius.r1r0sq) *
123 (scalingDensity()/dense.gas_phase[ipHYDROGEN]);
124 }
125 else
126 {
127 fprintf(ioQQQ,"chDenseLaw==\"WIND\" must now be ballistic or static\n");
129 }
130 }
131
132 else if( strcmp(dense.chDenseLaw,"SINE") == 0 )
133 {
134 /* rapid density fluctuation */
135 if( dense.lgDenFlucRadius )
136 {
137 new_density =
138 (dense.cfirst*cos(radius.depth*dense.flong+dense.flcPhase) +
139 dense.csecnd)*(scalingDensity()/dense.gas_phase[ipHYDROGEN]);
140 }
141 else
142 {
143 new_density =
144 (dense.cfirst*cos(colden.colden[ipCOL_HTOT]*dense.flong+dense.flcPhase) +
145 dense.csecnd)*(scalingDensity()/dense.gas_phase[ipHYDROGEN]);
146 }
147 }
148
149 else if( strcmp(dense.chDenseLaw,"POWR") == 0 )
150 {
151 /* power law function of radius */
152 double dnew = dense.den0*pow(radius.Radius/radius.rinner,(double)dense.DensityPower);
153 new_density = dnew*(scalingDensity()/dense.gas_phase[ipHYDROGEN]);
154 }
155
156 else if( strcmp(dense.chDenseLaw,"POWD") == 0 )
157 {
158 /* power law function of depth */
159 double dnew = dense.den0*pow(1. + radius.depth/dense.rscale,(double)dense.DensityPower);
160 new_density = dnew*(scalingDensity()/dense.gas_phase[ipHYDROGEN]);
161 }
162
163 else if( strcmp(dense.chDenseLaw,"POWC") == 0 )
164 {
165 /* power law function of column density */
166 double dnew = dense.den0*pow(1.f + colden.colden[ipCOL_HTOT]/
167 dense.rscale,dense.DensityPower);
168 new_density = dnew*(scalingDensity()/dense.gas_phase[ipHYDROGEN]);
169 }
170
171
172 else if( strncmp( dense.chDenseLaw ,"DLW" , 3) == 0 )
173 {
174 if( strcmp(dense.chDenseLaw,"DLW1") == 0 )
175 {
176 /* call ACF sub */
177 new_density =
178 dense_fabden(radius.Radius,radius.depth)*
179 (scalingDensity()/dense.gas_phase[ipHYDROGEN]);
180 }
181 else if( strcmp(dense.chDenseLaw,"DLW2") == 0 )
182 {
183 /* call table interpolation subroutine
184 * >>chng 96 nov 29, added dense_tabden */
185 new_density = dense_tabden(radius.Radius,radius.depth) *
186 (scalingDensity()/dense.gas_phase[ipHYDROGEN]);
187 }
188 else if( strcmp(dense.chDenseLaw,"DLW3") == 0 )
189 {
190 /* call parametrized wind subroutine */
191 new_density = dense_parametric_wind(radius.Radius) *
192 (scalingDensity()/dense.gas_phase[ipHYDROGEN]);
193 }
194 else
195 {
196 fprintf( ioQQQ, " Insanity, lgConvPres gets chCPres=%4.4s\n",
197 dense.chDenseLaw );
199 }
200
201 if( new_density/scalingDensity() > 3. || new_density/scalingDensity()< 1./3 )
202 {
203 static bool lgWARN2BIG=false;
204 if( !lgWARN2BIG )
205 {
206 lgWARN2BIG = true;
207 fprintf(ioQQQ,"\n\n >========== Warning! The tabulated or functional change in density as a function of depth was VERY large. This is zone %li.\n",nzone);
208 fprintf(ioQQQ," >========== Warning! This will cause convergence problems. \n");
209 fprintf(ioQQQ," >========== Warning! The current radius is %.3e. \n",radius.Radius);
210 fprintf(ioQQQ," >========== Warning! The current depth is %.3e. \n",radius.depth);
211 fprintf(ioQQQ," >========== Warning! Consider using a more modest change in density vs radius. \n\n\n");
212 }
213 }
214 }
215
216 else if( strcmp(dense.chDenseLaw,"CDEN") == 0 )
217 {
218 /* this is the default, constant density */
219 new_density = scalingDensity();
220 }
221
222 else
223 {
224 fprintf( ioQQQ, " Unknown density law=%s= This is a major internal error.\n",
225 dense.chDenseLaw );
226 ShowMe();
228 }
229
230 return new_density;
231}
232
233enum {
242};
243
245STATIC double stepDensity(const PresMode & presmode, solverState &st);
246
248{
249 /* >>chng 04 feb 11, add option to remember current density and pressure */
250 conv.hist_pres_density.push_back( dense.gas_phase[ipHYDROGEN] );
251 conv.hist_pres_current.push_back( pressure.PresTotlCurr );
252 conv.hist_pres_error.push_back( pressure.PresTotlError );
253}
254
255STATIC bool lgTestPressureConvergence( double new_density)
256{
257 double density_change_factor = new_density/scalingDensity();
258
259 /* now see whether current pressure is within error bounds */
260 if( density_change_factor > 1. + conv.PressureErrorAllowed ||
261 density_change_factor < 1. - conv.PressureErrorAllowed )
262 {
263 return false;
264 }
265 else
266 {
267 return true;
268 }
269}
270
271STATIC double limitedDensityScaling( double new_density, double dP_chng_factor )
272{
273 double density_change_factor = new_density/scalingDensity();
274
275 /* dP_chng_factor is initially 1, becomes smaller if sign of pressure change, changes */
276 double pdelta = conv.MaxFractionalDensityStepPerIteration * dP_chng_factor;
277
278 /* make sure that change is not too extreme */
279 density_change_factor = MIN2(density_change_factor,1.+pdelta);
280 density_change_factor = MAX2(density_change_factor,1.-pdelta);
281 return density_change_factor;
282}
283
284/*PressureChange evaluate the current pressure, and change needed to
285 * get it to PresTotlInit,
286 * return value is true if density was changed, false if no changes were necessary */
288 /* this is change factor, 1 at first, becomes smaller as oscillations occur */
289 double dP_chng_factor, const PresMode & presmode, solverState &st )
290{
291 DEBUG_ENTRY( "PressureChange()" );
292
293 /* first evaluate total pressure for this location, and current conditions
294 * CurrentPressure is just sum of gas plus local line radiation pressure */
295 /* this sets values of pressure.PresTotlCurr, also wind velocity for dynamics */
298
299 double new_density = stepDensity(presmode, st);
300
301 conv.lgConvPres = lgTestPressureConvergence(new_density);
302
303 double density_change_factor = limitedDensityScaling(new_density, dP_chng_factor);
304
305 {
306 /*@-redef@*/
307 enum{DEBUG_LOC=false};
308 static long int nsave=-1;
309 /*@+redef@*/
310 if( DEBUG_LOC /*&& nzone > 150 && iteration > 1*/ )
311 {
312 if( nsave-nzone ) fprintf(ioQQQ,"\n");
313 nsave = nzone;
314 fprintf(ioQQQ,"nnzzone\t%li\t%.2f%%\t%.3f\t%.2e\t%.2e\t%.2e\n",
315 nzone,
316 density_change_factor,
317 /* when this is negative we need to raise the density */
318 pressure.PresTotlError*100.,
319 pressure.PresTotlCurr,
320 pressure.PresGasCurr,
321 pressure.pres_radiation_lines_curr );
322 }
323 }
324
325 if( trace.lgTrace )
326 {
327 fprintf( ioQQQ,
328 " PressureChange called, changing HDEN from %10.3e to %10.3e Set fill fac to %10.3e\n",
329 dense.gas_phase[ipHYDROGEN], density_change_factor*dense.gas_phase[ipHYDROGEN],
330 geometry.FillFac );
331 }
332
333 bool lgChange = ( density_change_factor != 1. );
334
335 if( lgChange )
336 {
337 ScaleAllDensities((realnum) density_change_factor);
338
339 /* must call TempChange since ionization has changed, there are some
340 * terms that affect collision rates (H0 term in electron collision) */
341 TempChange(phycon.te , false);
342 }
343
344 {
345 /*@-redef@*/
346 enum {DEBUG_LOC=false};
347 /*@+redef@*/
348 if( DEBUG_LOC && (nzone>215)/**/ )
349 {
350 fprintf( ioQQQ,
351 "%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%.2e\t%c\n",
352 radius.depth,
353 pressure.PresTotlCurr,
354 pressure.PresTotlInit + pressure.PresInteg,
355 pressure.PresTotlInit,
356 pressure.PresGasCurr,
357 pressure.PresRamCurr,
358 pressure.pres_radiation_lines_curr,
359 /* subtract continuum rad pres which has already been added on */
360 pressure.PresInteg - pressure.pinzon,
361 wind.windv/1e5,
362 sqrt(5.*pressure.PresGasCurr/3./dense.xMassDensity)/1e5,
363 TorF(conv.lgConvPres) );
364 }
365 }
366
367 return lgChange;
368}
369
370/* ============================================================================== */
371/* DynaPresChngFactor, called from PressureChange to evaluate factor needed
372 * to find new density needed for
373 * current conditions and wind solution, returns ratio of new to old density */
374
375/* object is to get the local ram pressure
376 * RamNeeded = pressure.PresTotlInit + pressure.PresGasCurr + pressure.PresInteg;
377 * to equal the sum of the inital pressur at the illuminated face, the local gas pressure,
378 * and the integrate radiative acceleration from the incident continuum
379 *
380 * the local gas pressure is linear in density if the temperature is constant,
381 *
382 * the local ram pressure is inversely linear in density because of the relationship
383 * between velocity and density introduced by the mass flux conservation
384 */
385
386/* stepDensity finds a density which should move towards pressure equilibrium
387 * sets pressure.PresTotlError */
388STATIC double stepDensity(const PresMode &presmode, solverState &st)
389{
390 DEBUG_ENTRY( "stepDensity()" );
391
392 double PresTotlCorrect = st.press;
393
394 double densityCurrent = scalingDensity();
395
396 double er = pressure.PresTotlCurr-PresTotlCorrect;
397 pressure.PresTotlError = er/PresTotlCorrect;
398
399 if( trace.lgTrace && presmode.global != CPRE )
400 {
401 fprintf( ioQQQ,
402 " DynaPresChngFactor set PresTotlCorrect=%.3e PresTotlInit=%.3e PresInteg=%.3e DivergePresInteg=%.3e\n",
403 PresTotlCorrect , pressure.PresTotlInit , pressure.PresInteg*pressure.lgContRadPresOn,
404 dynamics.DivergePresInteg );
405 }
406
407 if( dynamics.lgTracePrint && presmode.global != CPRE)
408 fprintf(ioQQQ,"Pressure: init %g +rad %g +diverge %g = %g cf %g\n",
409 pressure.PresTotlInit, pressure.PresInteg*pressure.lgContRadPresOn,
410 dynamics.DivergePresInteg, PresTotlCorrect, pressure.PresTotlCurr);
411 /*fprintf(ioQQQ,"DEBUG\t%.2f\thden\t%.4e\tPtot\t%.4e\tPgas\t%.4e\n",
412 fnzone, scalingDensity(),pressure.PresTotlCurr,pressure.PresGasCurr );*/
413
414 double rhohat;
415 if( presmode.global == CPRE || presmode.global == ORIGINAL ||
416 st.lastzone != nzone || fabs(er-st.erp) < SMALLFLOAT
417 || fabs(densityCurrent-st.dp) < SMALLFLOAT)
418 {
419 double dpdrho;
420 if ( presmode.zone == SUBSONIC )
421 {
422 dpdrho = pressure.PresTotlCurr/densityCurrent;
423 }
424 else
425 {
426 dpdrho = -pressure.PresTotlCurr/densityCurrent;
427 }
428 rhohat = densityCurrent - er/dpdrho;
429 }
430 else
431 {
432 /* second iteration on this zone, do linear fit to previous Pres - rho curve
433 * Linear approximation to guess root with two independent samples */
434 double dpdrho = (st.erp-er)/(st.dp-densityCurrent);
435 rhohat = densityCurrent - er/dpdrho;
436
437 /* Subsonic case: pressure ^ with density ^ => increase density further */
438 /* Super " case: pressure ^ with density v => decrease density further */
439
440 /* Force the solution towards the required branch when it
441 * appears to have the "wrong" value of dP/drho */
442 if(presmode.zone == SUBSONIC && dpdrho <= 0)
443 {
444 rhohat = 1.03*densityCurrent;
445 }
446 else if(presmode.zone == SUPERSONIC && dpdrho >= 0)
447 {
448 rhohat = 0.97*densityCurrent;
449 }
450 }
451
452 st.dp = densityCurrent;
453 st.erp = er;
454 st.lastzone = nzone;
455
456 if( presmode.global != CPRE && dynamics.lgTracePrint )
457 fprintf(ioQQQ,"windv %li r %g v %g f %g\n",
458 nzone,radius.depth,wind.windv,DynaFlux(radius.depth));
459
460 /* convergence trace at this level */
461 if( presmode.global != CPRE && trace.nTrConvg>=1 )
462 {
463 fprintf( ioQQQ,
464 " DynaPresChngFactor mode %s new density %.3f vel %.3e\n",
465 dynamics.chPresMode , rhohat , wind.windv );
466 }
467
468 return rhohat;
469}
470
472{
473 DEBUG_ENTRY("PresMode::set()");
474 /* dynamics.lgSetPresMode is flag to indicate sane value of dynamics.chPresMode.
475 * If set true with SET DYNAMICS PRESSURE MODE
476 * then do not override that choice */
477 if( !dynamics.lgSetPresMode )
478 {
479 /* above set true if pressure mode was set with
480 * SET DYNAMICS PRESSURE MODE - if we got here
481 * it was not set, and must make a guess */
482 if(pressure.PresGasCurr < pressure.PresRamCurr)
483 strcpy( dynamics.chPresMode , "supersonic" );
484 else
485 strcpy( dynamics.chPresMode , "subsonic" );
486 /* clear the flag - pressure mode has been set */
487 dynamics.lgSetPresMode = true;
488 }
489
490 /* if globally looking for transonic solution, then locally sometimes
491 * supersonic, sometimes subsonic - this branch sets global flag,
492 * which can also be set with SET DYNAMICS PRESSURE MODE.
493 * Under default conditions, ChPresMode was set in previous branch
494 * to sub or super sonic depending on current velocity on first time*/
495
496 if (strcmp(dense.chDenseLaw,"CPRE") != 0)
497 {
498 ASSERT (strcmp(dense.chDenseLaw,"DYNA") == 0);
499 }
500
501 if (strcmp(dense.chDenseLaw,"CPRE") == 0)
502 {
503 global = CPRE;
504 pressure.lgSonicPointAbortOK = false;
505 }
506 else if( strcmp( dynamics.chPresMode , "original" ) == 0 )
507 {
509 pressure.lgSonicPointAbortOK = true;
510 }
511 else if( strcmp( dynamics.chPresMode , "subsonic" ) == 0 )
512 {
514 pressure.lgSonicPointAbortOK = true;
515 }
516 /* supersonic */
517 else if( strcmp( dynamics.chPresMode , "supersonic" ) == 0 )
518 {
520 pressure.lgSonicPointAbortOK = true;
521 }
522 /* strong d */
523 else if( strcmp( dynamics.chPresMode , "strongd" ) == 0 )
524 {
525 global = STRONGD;
526 pressure.lgSonicPointAbortOK = false;
527 }
528 else if( strcmp( dynamics.chPresMode , "shock" ) == 0 )
529 {
530 global = SHOCK;
531 pressure.lgSonicPointAbortOK = false;
532 }
533 else if( strcmp( dynamics.chPresMode , "antishock-by-mach" ) == 0 )
534 {
536 pressure.lgSonicPointAbortOK = false;
537 }
538 else if( strcmp( dynamics.chPresMode , "antishock" ) == 0 )
539 {
541 pressure.lgSonicPointAbortOK = false;
542 }
543
544 /* this sets pressure mode for the current zone based on global mode
545 * and local conditions */
546 if(global == CPRE)
547 {
548 zone = SUBSONIC;
549 }
550 else if(global == ORIGINAL)
551 {
552 /* in this mode use comparison between ram and gas pressure to determine
553 * whether sub or super sonic */
554 if(pressure.PresGasCurr > pressure.PresRamCurr)
555 {
556 zone = SUBSONIC;
557 }
558 else
559 {
561 }
562 }
563 else if(global == STRONGD)
564 {
565 //if(nzone <= 1)
567 }
568 else if(global == SUBSONIC)
569 {
570 zone = SUBSONIC;
571 }
572 else if(global == SUPERSONIC)
573 {
575 }
576 else if(global == SHOCK)
577 {
578 if(radius.depth < dynamics.ShockDepth)
579 {
580 zone = SUBSONIC;
581 }
582 else
583 {
585 }
586 }
587 else if(global == ANTISHOCK)
588 {
589 if(radius.depth < dynamics.ShockDepth)
590 {
592 }
593 else
594 {
595 zone = SUBSONIC;
596 }
597 }
598 else if(global == ANTISHOCK2)
599 {
600 /* WJH 19 May 2004: This version of the antishock mode will
601 insert the antishock at the point where the isothermal Mach
602 number falls below a certain value, dynamics.ShockMach */
603 if( fabs(wind.windv) > dynamics.ShockMach*sqrt(pressure.PresGasCurr/dense.xMassDensity))
604 {
606 }
607 else
608 {
609 zone = SUBSONIC;
610 }
611 }
612 else
613 {
614 printf("Need to set global pressure mode\n");
616 }
617}
618
619double pressureZone(const PresMode &presmode)
620{
621 double press;
622 if( presmode.global == CPRE )
623 {
624 /* constant pressure */
625 if( pressure.lgContRadPresOn )
626 {
627 /* >>chng 01 oct 31, replace pneed with CorretPressure */
628 /* this has pressure due to incident continuum */
629 press = pressure.PresTotlInit + pressure.PresInteg;
630 }
631 else
632 {
633 /* this does not have pressure due to incident continuum*/
634 press = pressure.PresTotlInit*
635 /* following term normally unity, power law set with option par on cmmnd*/
636 pow(radius.Radius/radius.rinner,(double)pressure.PresPowerlaw);
637 }
638
639 if( dark.lgNFW_Set || pressure.gravity_symmetry >=0 )
640 {
641 fixit(); // doesn't this exclude current zone gravity pressure.RhoGravity?
642 // and doesn't the above exclude current rad press, pressure.pinzon?
643 // That would be a second order correction, which would need
644 // a consistent second order treatment elsewhere to make it
645 // worthwhile.
646 press += pressure.IntegRhoGravity;
647 }
648 }
649 else
650 {
651 /* update the current desired pressure */
652 press = pressure.PresTotlInit + pressure.PresInteg*pressure.lgContRadPresOn
653 + dynamics.DivergePresInteg;
654 }
655 return press;
656}
657
long int nzone
Definition cddefines.cpp:14
FILE * ioQQQ
Definition cddefines.cpp:7
#define ASSERT(exp)
Definition cddefines.h:578
#define STATIC
Definition cddefines.h:97
#define MIN2
Definition cddefines.h:761
#define POW2
Definition cddefines.h:929
#define EXIT_FAILURE
Definition cddefines.h:140
char TorF(bool l)
Definition cddefines.h:710
#define cdEXIT(FAIL)
Definition cddefines.h:434
float realnum
Definition cddefines.h:103
#define MAX2
Definition cddefines.h:782
NORETURN void TotalInsanity(void)
Definition service.cpp:886
const int ipHYDROGEN
Definition cddefines.h:305
#define DEBUG_ENTRY(funcname)
Definition cddefines.h:684
void ShowMe(void)
Definition service.cpp:181
void fixit(void)
Definition service.cpp:991
t_colden colden
Definition colden.cpp:5
#define ipCOL_HTOT
Definition colden.h:12
t_conv conv
Definition conv.cpp:5
t_cosmology cosmology
Definition cosmology.cpp:11
realnum GetDensity(realnum z)
Definition cosmology.cpp:31
const realnum SMALLFLOAT
Definition cpu.h:191
t_dark_matter dark
void ScaleAllDensities(realnum factor)
Definition dense.cpp:37
t_dense dense
Definition dense.cpp:24
realnum scalingDensity(void)
Definition dense.cpp:378
double dense_tabden(double r0, double depth)
double dense_fabden(double radius, double depth)
double dense_parametric_wind(double rad)
t_dynamics dynamics
Definition dynamics.cpp:44
realnum DynaFlux(double depth)
t_geometry geometry
Definition geometry.cpp:5
t_phycon phycon
Definition phycon.cpp:6
t_pressure pressure
Definition pressure.cpp:5
void PresTotCurrent(void)
STATIC bool lgTestPressureConvergence(double new_density)
@ STRONGD
@ SUBSONIC
@ SUPERSONIC
@ ANTISHOCK
@ ANTISHOCK2
@ ORIGINAL
STATIC void logPressureState()
double zoneDensity()
double pressureZone(const PresMode &presmode)
STATIC double stepDensity(const PresMode &presmode, solverState &st)
STATIC double limitedDensityScaling(double new_density, double dP_chng_factor)
bool PressureChange(double dP_chng_factor, const PresMode &presmode, solverState &st)
t_radius radius
Definition radius.cpp:5
t_struc struc
Definition struc.cpp:6
void TempChange(double TempNew, bool lgForceUpdate)
t_trace trace
Definition trace.cpp:5
Wind wind
Definition wind.cpp:5