This class represents 1D flow domains that satisfy the one-dimensional similarity solution for chemically-reacting, axisymmetric flows. More...
#include <Flow1D.h>
This class represents 1D flow domains that satisfy the one-dimensional similarity solution for chemically-reacting, axisymmetric flows.
Public Member Functions | |
Flow1D (ThermoPhase *ph=0, size_t nsp=1, size_t points=1) | |
Create a new flow domain. | |
Flow1D (shared_ptr< ThermoPhase > th, size_t nsp=1, size_t points=1) | |
Delegating constructor. | |
Flow1D (shared_ptr< Solution > sol, const string &id="", size_t points=1) | |
Create a new flow domain. | |
string | domainType () const override |
Domain type flag. | |
string | componentName (size_t n) const override |
Name of component n . May be overloaded. | |
size_t | componentIndex (const string &name) const override |
index of component with name name . | |
virtual bool | componentActive (size_t n) const |
Returns true if the specified component is an active part of the solver state. | |
void | show (const double *x) override |
Print the solution. | |
shared_ptr< SolutionArray > | asArray (const double *soln) const override |
Save the state of this domain as a SolutionArray. | |
void | fromArray (SolutionArray &arr, double *soln) override |
Restore the solution for this domain from a SolutionArray. | |
void | setFreeFlow () |
Set flow configuration for freely-propagating flames, using an internal point with a fixed temperature as the condition to determine the inlet mass flux. | |
void | setAxisymmetricFlow () |
Set flow configuration for axisymmetric counterflow flames, using specified inlet mass fluxes. | |
void | setUnstrainedFlow () |
Set flow configuration for burner-stabilized flames, using specified inlet mass fluxes. | |
void | solveEnergyEqn (size_t j=npos) |
Specify that the energy equation should be solved at point j . | |
virtual size_t | getSolvingStage () const |
Get the solving stage (used by IonFlow specialization) | |
virtual void | setSolvingStage (const size_t stage) |
Solving stage mode for handling ionized species (used by IonFlow specialization) | |
virtual void | solveElectricField (size_t j=npos) |
Set to solve electric field in a point (used by IonFlow specialization) | |
virtual void | fixElectricField (size_t j=npos) |
Set to fix voltage in a point (used by IonFlow specialization) | |
virtual bool | doElectricField (size_t j) const |
Retrieve flag indicating whether electric field is solved or not (used by IonFlow specialization) | |
void | enableRadiation (bool doRadiation) |
Turn radiation on / off. | |
bool | radiationEnabled () const |
Returns true if the radiation term in the energy equation is enabled. | |
double | radiativeHeatLoss (size_t j) const |
Return radiative heat loss at grid point j. | |
void | setBoundaryEmissivities (double e_left, double e_right) |
Set the emissivities for the boundary values. | |
double | leftEmissivity () const |
Return emissivity at left boundary. | |
double | rightEmissivity () const |
Return emissivity at right boundary. | |
void | fixTemperature (size_t j=npos) |
Specify that the the temperature should be held fixed at point j . | |
bool | doEnergy (size_t j) |
true if the energy equation is solved at point j or false if a fixed temperature condition is imposed. | |
void | resize (size_t components, size_t points) override |
Change the grid size. Called after grid refinement. | |
void | setGas (const double *x, size_t j) |
Set the gas object state to be consistent with the solution at point j. | |
void | setGasAtMidpoint (const double *x, size_t j) |
Set the gas state to be consistent with the solution at the midpoint between j and j + 1. | |
double | density (size_t j) const |
Get the density [kg/m³] at point j | |
bool | isFree () const |
Retrieve flag indicating whether flow is freely propagating. | |
bool | isStrained () const |
Retrieve flag indicating whether flow uses radial momentum. | |
void | setViscosityFlag (bool dovisc) |
Specify if the viscosity term should be included in the momentum equation. | |
void | eval (size_t jGlobal, double *xGlobal, double *rsdGlobal, integer *diagGlobal, double rdt) override |
Evaluate the residual functions for axisymmetric stagnation flow. | |
size_t | leftExcessSpecies () const |
Index of the species on the left boundary with the largest mass fraction. | |
size_t | rightExcessSpecies () const |
Index of the species on the right boundary with the largest mass fraction. | |
Problem Specification | |
void | setupGrid (size_t n, const double *z) override |
called to set up initial grid, and after grid refinement | |
void | resetBadValues (double *xg) override |
When called, this function should reset "bad" values in the state vector such as negative species concentrations. | |
ThermoPhase & | phase () |
Access the phase object used to compute thermodynamic properties for points in this domain. | |
Kinetics & | kinetics () |
Access the Kinetics object used to compute reaction rates for points in this domain. | |
void | setKinetics (shared_ptr< Kinetics > kin) override |
Set the Kinetics object used for reaction rate calculations. | |
void | setTransport (shared_ptr< Transport > trans) override |
Set the transport manager used for transport property calculations. | |
void | setTransportModel (const string &trans) |
Set the transport model. | |
string | transportModel () const |
Retrieve transport model. | |
void | enableSoret (bool withSoret) |
Enable thermal diffusion, also known as Soret diffusion. | |
bool | withSoret () const |
Indicates if thermal diffusion (Soret effect) term is being calculated. | |
void | setFluxGradientBasis (ThermoBasis fluxGradientBasis) |
Compute species diffusive fluxes with respect to their mass fraction gradients (fluxGradientBasis = ThermoBasis::mass) or mole fraction gradients (fluxGradientBasis = ThermoBasis::molar, default) when using the mixture-averaged transport model. | |
ThermoBasis | fluxGradientBasis () const |
Compute species diffusive fluxes with respect to their mass fraction gradients (fluxGradientBasis = ThermoBasis::mass) or mole fraction gradients (fluxGradientBasis = ThermoBasis::molar, default) when using the mixture-averaged transport model. | |
void | setPressure (double p) |
Set the pressure. | |
double | pressure () const |
The current pressure [Pa]. | |
void | _getInitialSoln (double *x) override |
Write the initial solution estimate into array x. | |
void | _finalize (const double *x) override |
In some cases, a domain may need to set parameters that depend on the initial solution estimate. | |
void | setFixedTempProfile (vector< double > &zfixed, vector< double > &tfixed) |
Sometimes it is desired to carry out the simulation using a specified temperature profile, rather than computing it by solving the energy equation. | |
void | setTemperature (size_t j, double t) |
Set the temperature fixed point at grid point j, and disable the energy equation so that the solution will be held to this value. | |
double | T_fixed (size_t j) const |
The fixed temperature value at point j. | |
Two-Point control method | |
In this method two control points are designated in the 1D domain, and the value of the temperature at these points is fixed. The values of the control points are imposed and thus serve as a boundary condition that affects the solution of the governing equations in the 1D domain. The imposition of fixed points in the domain means that the original set of governing equations' boundary conditions would over-determine the problem. Thus, the boundary conditions are changed to reflect the fact that the control points are serving as internal boundary conditions. The imposition of the two internal boundary conditions requires that two other boundary conditions be changed. The first is the boundary condition for the continuity equation at the left boundary, which is changed to be a value that is derived from the solution at the left boundary. The second is the continuity boundary condition at the right boundary, which is also determined from the flow solution by using the oxidizer axial velocity equation variable to compute the mass flux at the right boundary. This method is based on the work of Nishioka et al. [32] . | |
double | leftControlPointTemperature () const |
Returns the temperature at the left control point. | |
double | leftControlPointCoordinate () const |
Returns the z-coordinate of the left control point. | |
void | setLeftControlPointTemperature (double temperature) |
Sets the temperature of the left control point. | |
void | setLeftControlPointCoordinate (double z_left) |
Sets the coordinate of the left control point. | |
double | rightControlPointTemperature () const |
Returns the temperature at the right control point. | |
double | rightControlPointCoordinate () const |
Returns the z-coordinate of the right control point. | |
void | setRightControlPointTemperature (double temperature) |
Sets the temperature of the right control point. | |
void | setRightControlPointCoordinate (double z_right) |
Sets the coordinate of the right control point. | |
void | enableTwoPointControl (bool twoPointControl) |
Sets the status of the two-point control. | |
bool | twoPointControlEnabled () const |
Returns the status of the two-point control. | |
![]() | |
Domain1D (size_t nv=1, size_t points=1, double time=0.0) | |
Constructor. | |
Domain1D (const Domain1D &)=delete | |
Domain1D & | operator= (const Domain1D &)=delete |
virtual string | domainType () const |
Domain type flag. | |
string | type () const |
String indicating the domain implemented. | |
size_t | domainIndex () |
The left-to-right location of this domain. | |
virtual bool | isConnector () |
True if the domain is a connector domain. | |
void | setSolution (shared_ptr< Solution > sol) |
Set the solution manager. | |
virtual void | setKinetics (shared_ptr< Kinetics > kin) |
Set the kinetics manager. | |
virtual void | setTransport (shared_ptr< Transport > trans) |
Set transport model to existing instance. | |
const OneDim & | container () const |
The container holding this domain. | |
void | setContainer (OneDim *c, size_t index) |
Specify the container object for this domain, and the position of this domain in the list. | |
void | setBandwidth (int bw=-1) |
Set the Jacobian bandwidth. See the discussion of method bandwidth(). | |
size_t | bandwidth () |
Set the Jacobian bandwidth for this domain. | |
virtual void | init () |
Initialize. | |
virtual void | setInitialState (double *xlocal=0) |
virtual void | setState (size_t point, const double *state, double *x) |
virtual void | resetBadValues (double *xg) |
When called, this function should reset "bad" values in the state vector such as negative species concentrations. | |
virtual void | resize (size_t nv, size_t np) |
Resize the domain to have nv components and np grid points. | |
Refiner & | refiner () |
Return a reference to the grid refiner. | |
size_t | nComponents () const |
Number of components at each grid point. | |
void | checkComponentIndex (size_t n) const |
Check that the specified component index is in range. | |
void | checkComponentArraySize (size_t nn) const |
Check that an array size is at least nComponents(). | |
size_t | nPoints () const |
Number of grid points in this domain. | |
void | checkPointIndex (size_t n) const |
Check that the specified point index is in range. | |
void | checkPointArraySize (size_t nn) const |
Check that an array size is at least nPoints(). | |
virtual string | componentName (size_t n) const |
Name of component n . May be overloaded. | |
void | setComponentName (size_t n, const string &name) |
Set the name of the component n to name . | |
virtual size_t | componentIndex (const string &name) const |
index of component with name name . | |
void | setBounds (size_t n, double lower, double upper) |
Set the upper and lower bounds for a solution component, n. | |
void | setTransientTolerances (double rtol, double atol, size_t n=npos) |
Set tolerances for time-stepping mode. | |
void | setSteadyTolerances (double rtol, double atol, size_t n=npos) |
Set tolerances for steady-state mode. | |
double | rtol (size_t n) |
Relative tolerance of the nth component. | |
double | atol (size_t n) |
Absolute tolerance of the nth component. | |
double | steady_rtol (size_t n) |
Steady relative tolerance of the nth component. | |
double | steady_atol (size_t n) |
Steady absolute tolerance of the nth component. | |
double | transient_rtol (size_t n) |
Transient relative tolerance of the nth component. | |
double | transient_atol (size_t n) |
Transient absolute tolerance of the nth component. | |
double | upperBound (size_t n) const |
Upper bound on the nth component. | |
double | lowerBound (size_t n) const |
Lower bound on the nth component. | |
void | initTimeInteg (double dt, const double *x0) |
Performs the setup required before starting a time-stepping solution. | |
void | setSteadyMode () |
Set the internally-stored reciprocal of the time step to 0.0, which is used to indicate that the problem is in steady-state mode. | |
bool | steady () |
True if in steady-state mode. | |
bool | transient () |
True if not in steady-state mode. | |
void | needJacUpdate () |
Set this if something has changed in the governing equations (for example, the value of a constant has been changed, so that the last-computed Jacobian is no longer valid. | |
virtual void | eval (size_t j, double *x, double *r, integer *mask, double rdt=0.0) |
Evaluate the residual function at point j. | |
size_t | index (size_t n, size_t j) const |
Returns the index of the solution vector, which corresponds to component n at grid point j. | |
double | value (const double *x, size_t n, size_t j) const |
Returns the value of solution component n at grid point j of the solution vector x. | |
virtual void | setJac (MultiJac *jac) |
virtual shared_ptr< SolutionArray > | asArray (const double *soln) const |
Save the state of this domain as a SolutionArray. | |
shared_ptr< SolutionArray > | toArray (bool normalize=false) const |
Save the state of this domain to a SolutionArray. | |
virtual void | fromArray (SolutionArray &arr, double *soln) |
Restore the solution for this domain from a SolutionArray. | |
void | fromArray (const shared_ptr< SolutionArray > &arr) |
Restore the solution for this domain from a SolutionArray. | |
shared_ptr< Solution > | solution () const |
Return thermo/kinetics/transport manager used in the domain. | |
size_t | size () const |
Return the size of the solution vector (the product of m_nv and m_points). | |
void | locate () |
Find the index of the first grid point in this domain, and the start of its variables in the global solution vector. | |
virtual size_t | loc (size_t j=0) const |
Location of the start of the local solution vector in the global solution vector. | |
size_t | firstPoint () const |
The index of the first (that is, left-most) grid point belonging to this domain. | |
size_t | lastPoint () const |
The index of the last (that is, right-most) grid point belonging to this domain. | |
void | linkLeft (Domain1D *left) |
Set the left neighbor to domain 'left. | |
void | linkRight (Domain1D *right) |
Set the right neighbor to domain 'right.'. | |
void | append (Domain1D *right) |
Append domain 'right' to this one, and update all links. | |
Domain1D * | left () const |
Return a pointer to the left neighbor. | |
Domain1D * | right () const |
Return a pointer to the right neighbor. | |
double | prevSoln (size_t n, size_t j) const |
Value of component n at point j in the previous solution. | |
void | setID (const string &s) |
Specify an identifying tag for this domain. | |
string | id () const |
Returns the identifying tag for this domain. | |
virtual void | show (std::ostream &s, const double *x) |
Print the solution. | |
virtual void | show (const double *x) |
Print the solution. | |
double | z (size_t jlocal) const |
Get the coordinate [m] of the point with local index jlocal | |
double | zmin () const |
Get the coordinate [m] of the first (leftmost) grid point in this domain. | |
double | zmax () const |
Get the coordinate [m] of the last (rightmost) grid point in this domain. | |
void | setProfile (const string &name, double *values, double *soln) |
Set initial values for a component at each grid point. | |
vector< double > & | grid () |
Access the array of grid coordinates [m]. | |
const vector< double > & | grid () const |
Access the array of grid coordinates [m]. | |
double | grid (size_t point) const |
virtual void | setupGrid (size_t n, const double *z) |
called to set up initial grid, and after grid refinement | |
virtual void | _getInitialSoln (double *x) |
Writes some or all initial solution values into the global solution array, beginning at the location pointed to by x. | |
virtual double | initialValue (size_t n, size_t j) |
Initial value of solution component n at grid point j. | |
virtual void | _finalize (const double *x) |
In some cases, a domain may need to set parameters that depend on the initial solution estimate. | |
void | forceFullUpdate (bool update) |
In some cases, for computational efficiency some properties (such as transport coefficients) may not be updated during Jacobian evaluations. | |
void | setData (shared_ptr< vector< double > > &data) |
Set shared data pointer. | |
Public Attributes | |
double | m_zfixed = Undef |
Location of the point where temperature is fixed. | |
double | m_tfixed = -1.0 |
Temperature at the point used to fix the flame location. | |
Protected Member Functions | |
AnyMap | getMeta () const override |
Retrieve meta data. | |
void | setMeta (const AnyMap &state) override |
Retrieve meta data. | |
virtual void | evalContinuity (size_t j, double *x, double *r, int *diag, double rdt) |
Alternate version of evalContinuity with legacy signature. | |
virtual void | evalUo (double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) |
Evaluate the oxidizer axial velocity equation residual. | |
double | shear (const double *x, size_t j) const |
Compute the shear term from the momentum equation using a central three-point differencing scheme. | |
double | conduction (const double *x, size_t j) const |
Compute the conduction term from the energy equation using a central three-point differencing scheme. | |
size_t | mindex (size_t k, size_t j, size_t m) |
Array access mapping for a 3D array stored in a 1D vector. | |
virtual void | grad_hk (const double *x, size_t j) |
Compute the spatial derivative of species specific molar enthalpies using upwind differencing. | |
Updates of cached properties | |
These methods are called by eval() to update cached properties and data that are used for the evaluation of the governing equations. | |
void | updateThermo (const double *x, size_t j0, size_t j1) |
Update the thermodynamic properties from point j0 to point j1 (inclusive), based on solution x. | |
virtual void | updateTransport (double *x, size_t j0, size_t j1) |
Update the transport properties at grid points in the range from j0 to j1 , based on solution x . | |
virtual void | updateDiffFluxes (const double *x, size_t j0, size_t j1) |
Update the diffusive mass fluxes. | |
virtual void | updateProperties (size_t jg, double *x, size_t jmin, size_t jmax) |
Update the properties (thermo, transport, and diffusion flux). | |
void | computeRadiation (double *x, size_t jmin, size_t jmax) |
Computes the radiative heat loss vector over points jmin to jmax and stores the data in the qdotRadiation variable. | |
Governing Equations | |
Methods called by eval() to calculate residuals for individual governing equations. | |
virtual void | evalContinuity (double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) |
Evaluate the continuity equation residual. | |
virtual void | evalMomentum (double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) |
Evaluate the momentum equation residual. | |
virtual void | evalLambda (double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) |
Evaluate the lambda equation residual. | |
virtual void | evalEnergy (double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) |
Evaluate the energy equation residual. | |
virtual void | evalSpecies (double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) |
Evaluate the species equations' residuals. | |
virtual void | evalElectricField (double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) |
Evaluate the electric field equation residual to be zero everywhere. | |
Solution components | |
double | T (const double *x, size_t j) const |
Get the temperature at point j from the local state vector x . | |
double & | T (double *x, size_t j) |
Get the temperature at point j from the local state vector x . | |
double | T_prev (size_t j) const |
Get the temperature at point j from the previous time step. | |
double | rho_u (const double *x, size_t j) const |
Get the axial mass flux [kg/m²/s] at point j from the local state vector x . | |
double | u (const double *x, size_t j) const |
Get the axial velocity [m/s] at point j from the local state vector x . | |
double | V (const double *x, size_t j) const |
Get the spread rate (tangential velocity gradient) [1/s] at point j from the local state vector x . | |
double | V_prev (size_t j) const |
Get the spread rate [1/s] at point j from the previous time step. | |
double | lambda (const double *x, size_t j) const |
Get the radial pressure gradient [N/m⁴] at point j from the local state vector x | |
double | Uo (const double *x, size_t j) const |
Get the oxidizer inlet velocity [m/s] linked to point j from the local state vector x . | |
double | Y (const double *x, size_t k, size_t j) const |
Get the mass fraction of species k at point j from the local state vector x . | |
double & | Y (double *x, size_t k, size_t j) |
Get the mass fraction of species k at point j from the local state vector x . | |
double | Y_prev (size_t k, size_t j) const |
Get the mass fraction of species k at point j from the previous time step. | |
double | X (const double *x, size_t k, size_t j) const |
Get the mole fraction of species k at point j from the local state vector x . | |
double | flux (size_t k, size_t j) const |
Get the diffusive mass flux [kg/m²/s] of species k at point j | |
Convective spatial derivatives | |
These methods use upwind differencing to calculate spatial derivatives for velocity, species mass fractions, and temperature. Upwind differencing is a numerical discretization method that considers the direction of the flow to improve stability. | |
double | dVdz (const double *x, size_t j) const |
Calculates the spatial derivative of velocity V with respect to z at point j using upwind differencing. | |
double | dYdz (const double *x, size_t k, size_t j) const |
Calculates the spatial derivative of the species mass fraction Y_k with respect to z for species k at point j using upwind differencing. | |
double | dTdz (const double *x, size_t j) const |
Calculates the spatial derivative of temperature T with respect to z at point j using upwind differencing. | |
virtual AnyMap | getMeta () const |
Retrieve meta data. | |
virtual void | setMeta (const AnyMap &meta) |
Retrieve meta data. | |
Protected Attributes | |
double | m_press = -1.0 |
pressure [Pa] | |
vector< double > | m_dz |
Grid spacing. Element j holds the value of z(j+1) - z(j) . | |
vector< double > | m_rho |
Density at each grid point. | |
vector< double > | m_wtm |
Mean molecular weight at each grid point. | |
vector< double > | m_wt |
Molecular weight of each species. | |
vector< double > | m_cp |
Specific heat capacity at each grid point. | |
vector< double > | m_visc |
Dynamic viscosity at each grid point [Pa∙s]. | |
vector< double > | m_tcon |
Thermal conductivity at each grid point [W/m/K]. | |
vector< double > | m_diff |
Coefficient used in diffusion calculations for each species at each grid point. | |
vector< double > | m_multidiff |
Vector of size m_nsp × m_nsp × m_points for saving multicomponent diffusion coefficients. | |
Array2D | m_dthermal |
Array of size m_nsp by m_points for saving thermal diffusion coefficients. | |
Array2D | m_flux |
Array of size m_nsp by m_points for saving diffusive mass fluxes. | |
Array2D | m_hk |
Array of size m_nsp by m_points for saving molar enthalpies. | |
Array2D | m_dhk_dz |
Array of size m_nsp by m_points-1 for saving enthalpy fluxes. | |
Array2D | m_wdot |
Array of size m_nsp by m_points for saving species production rates. | |
size_t | m_nsp |
Number of species in the mechanism. | |
ThermoPhase * | m_thermo = nullptr |
Phase object used for calculating thermodynamic properties. | |
Kinetics * | m_kin = nullptr |
Kinetics object used for calculating species production rates. | |
Transport * | m_trans = nullptr |
Transport object used for calculating transport properties. | |
double | m_epsilon_left = 0.0 |
Emissivity of the surface to the left of the domain. | |
double | m_epsilon_right = 0.0 |
Emissivity of the surface to the right of the domain. | |
vector< size_t > | m_kRadiating |
Indices within the ThermoPhase of the radiating species. | |
vector< double > | m_qdotRadiation |
radiative heat loss vector | |
vector< double > | m_fixedtemp |
Fixed values of the temperature at each grid point that are used when solving with the energy equation disabled. | |
vector< double > | m_zfix |
Relative coordinates used to specify a fixed temperature profile. | |
vector< double > | m_tfix |
Fixed temperature values at the relative coordinates specified in m_zfix. | |
size_t | m_kExcessLeft = 0 |
Index of species with a large mass fraction at the left boundary, for which the mass fraction may be calculated as 1 minus the sum of the other mass fractions. | |
size_t | m_kExcessRight = 0 |
Index of species with a large mass fraction at the right boundary, for which the mass fraction may be calculated as 1 minus the sum of the other mass fractions. | |
double | m_zLeft = Undef |
Location of the left control point when two-point control is enabled. | |
double | m_tLeft = Undef |
Temperature of the left control point when two-point control is enabled. | |
double | m_zRight = Undef |
Location of the right control point when two-point control is enabled. | |
double | m_tRight = Undef |
Temperature of the right control point when two-point control is enabled. | |
flags | |
vector< bool > | m_do_energy |
For each point in the domain, true if energy equation is solved or false if temperature is held constant. | |
bool | m_do_soret = false |
true if the Soret diffusion term should be calculated. | |
ThermoBasis | m_fluxGradientBasis = ThermoBasis::molar |
Determines whether diffusive fluxes are computed using gradients of mass fraction or mole fraction. | |
bool | m_do_multicomponent = false |
true if transport fluxes are computed using the multicomponent diffusion coefficients, or false if mixture-averaged diffusion coefficients are used. | |
bool | m_do_radiation = false |
Determines whether radiative heat loss is calculated. | |
bool | m_dovisc |
Determines whether the viscosity term in the momentum equation is calculated. | |
bool | m_isFree |
Flag that is true for freely propagating flames anchored by a temperature fixed point. | |
bool | m_usesLambda |
Flag that is true for counterflow configurations that use the pressure eigenvalue \Lambda in the radial momentum equation. | |
bool | m_twoPointControl = false |
Flag for activating two-point flame control. | |
![]() | |
shared_ptr< vector< double > > | m_state |
data pointer shared from OneDim | |
double | m_rdt = 0.0 |
Reciprocal of the time step. | |
size_t | m_nv = 0 |
Number of solution components. | |
size_t | m_points |
Number of grid points. | |
vector< double > | m_slast |
Solution vector at the last time step. | |
vector< double > | m_max |
Upper bounds on solution components. | |
vector< double > | m_min |
Lower bounds on solution components. | |
vector< double > | m_rtol_ss |
Relative tolerances for steady mode. | |
vector< double > | m_rtol_ts |
Relative tolerances for transient mode. | |
vector< double > | m_atol_ss |
Absolute tolerances for steady mode. | |
vector< double > | m_atol_ts |
Absolute tolerances for transient mode. | |
vector< double > | m_z |
1D spatial grid coordinates | |
OneDim * | m_container = nullptr |
Parent OneDim simulation containing this and adjacent domains. | |
size_t | m_index |
Left-to-right location of this domain. | |
size_t | m_iloc = 0 |
Starting location within the solution vector for unknowns that correspond to this domain. | |
size_t | m_jstart = 0 |
Index of the first point in this domain in the global point list. | |
Domain1D * | m_left = nullptr |
Pointer to the domain to the left. | |
Domain1D * | m_right = nullptr |
Pointer to the domain to the right. | |
string | m_id |
Identity tag for the domain. | |
unique_ptr< Refiner > | m_refiner |
Refiner object used for placing grid points. | |
vector< string > | m_name |
Names of solution components. | |
int | m_bw = -1 |
See bandwidth() | |
bool | m_force_full_update = false |
see forceFullUpdate() | |
shared_ptr< Solution > | m_solution |
Composite thermo/kinetics/transport handler. | |
Private Attributes | |
vector< double > | m_ybar |
Holds the average of the species mass fractions between grid points j and j+1. | |
Flow1D | ( | ThermoPhase * | ph = 0 , |
size_t | nsp = 1 , |
||
size_t | points = 1 |
||
) |
Create a new flow domain.
ph | Object representing the gas phase. This object will be used to evaluate all thermodynamic, kinetic, and transport properties. |
nsp | Number of species. |
points | Initial number of grid points |
Definition at line 19 of file Flow1D.cpp.
Flow1D | ( | shared_ptr< ThermoPhase > | th, |
size_t | nsp = 1 , |
||
size_t | points = 1 |
||
) |
Delegating constructor.
Definition at line 90 of file Flow1D.cpp.
Create a new flow domain.
sol | Solution object used to evaluate all thermodynamic, kinetic, and transport properties |
id | name of flow domain |
points | initial number of grid points |
Definition at line 98 of file Flow1D.cpp.
~Flow1D | ( | ) |
Definition at line 116 of file Flow1D.cpp.
|
overridevirtual |
Domain type flag.
string
. Reimplemented from Domain1D.
Reimplemented in IonFlow.
Definition at line 123 of file Flow1D.cpp.
|
overridevirtual |
called to set up initial grid, and after grid refinement
Reimplemented from Domain1D.
Definition at line 185 of file Flow1D.cpp.
|
overridevirtual |
When called, this function should reset "bad" values in the state vector such as negative species concentrations.
This function may be called after a failed solution attempt.
Reimplemented from Domain1D.
Definition at line 200 of file Flow1D.cpp.
|
inline |
|
inline |
|
overridevirtual |
Set the Kinetics object used for reaction rate calculations.
Reimplemented from Domain1D.
Definition at line 133 of file Flow1D.cpp.
|
overridevirtual |
Set the transport manager used for transport property calculations.
Reimplemented from Domain1D.
Definition at line 139 of file Flow1D.cpp.
void setTransportModel | ( | const string & | trans | ) |
string transportModel | ( | ) | const |
|
inline |
|
inline |
void setFluxGradientBasis | ( | ThermoBasis | fluxGradientBasis | ) |
Compute species diffusive fluxes with respect to their mass fraction gradients (fluxGradientBasis = ThermoBasis::mass) or mole fraction gradients (fluxGradientBasis = ThermoBasis::molar, default) when using the mixture-averaged transport model.
fluxGradientBasis | set flux computation to mass or mole basis |
Definition at line 219 of file Flow1D.cpp.
|
inline |
Compute species diffusive fluxes with respect to their mass fraction gradients (fluxGradientBasis = ThermoBasis::mass) or mole fraction gradients (fluxGradientBasis = ThermoBasis::molar, default) when using the mixture-averaged transport model.
|
inline |
|
inline |
|
overridevirtual |
Write the initial solution estimate into array x.
Reimplemented from Domain1D.
Definition at line 229 of file Flow1D.cpp.
|
overridevirtual |
In some cases, a domain may need to set parameters that depend on the initial solution estimate.
In such cases, the parameters may be set in method _finalize. This method is called just before the Newton solver is called, and the x array is guaranteed to be the local solution vector for this domain that will be used as the initial guess. If no such parameters need to be set, then method _finalize does not need to be overloaded.
Reimplemented from Domain1D.
Reimplemented in IonFlow.
Definition at line 258 of file Flow1D.cpp.
|
inline |
|
inline |
|
inline |
|
overridevirtual |
Name of component n
. May be overloaded.
Reimplemented from Domain1D.
Definition at line 799 of file Flow1D.cpp.
|
overridevirtual |
index of component with name name
.
Reimplemented from Domain1D.
Definition at line 823 of file Flow1D.cpp.
|
virtual |
Returns true if the specified component is an active part of the solver state.
Reimplemented in IonFlow.
Definition at line 848 of file Flow1D.cpp.
|
overridevirtual |
Print the solution.
x | Pointer to the local portion of the system state vector |
Reimplemented from Domain1D.
Definition at line 782 of file Flow1D.cpp.
|
overridevirtual |
Save the state of this domain as a SolutionArray.
soln | local solution vector for this domain |
Reimplemented from Domain1D.
Definition at line 915 of file Flow1D.cpp.
|
overridevirtual |
Restore the solution for this domain from a SolutionArray.
[in] | arr | SolutionArray defining the state of this domain |
[out] | soln | Value of the solution vector, local to this domain |
Reimplemented from Domain1D.
Definition at line 949 of file Flow1D.cpp.
|
inline |
|
inline |
|
inline |
void solveEnergyEqn | ( | size_t | j = npos | ) |
Specify that the energy equation should be solved at point j
.
The converse of this method is fixTemperature().
j | Point at which to enable the energy equation. npos means all points. |
Definition at line 1046 of file Flow1D.cpp.
|
virtual |
Get the solving stage (used by IonFlow specialization)
Reimplemented in IonFlow.
Definition at line 1070 of file Flow1D.cpp.
|
virtual |
Solving stage mode for handling ionized species (used by IonFlow specialization)
stage=1
: the fluxes of charged species are set to zerostage=2
: the electric field equation is solved, and the drift flux for ionized species is evaluated Reimplemented in IonFlow.
Definition at line 1076 of file Flow1D.cpp.
|
virtual |
Set to solve electric field in a point (used by IonFlow specialization)
Reimplemented in IonFlow.
Definition at line 1082 of file Flow1D.cpp.
|
virtual |
Set to fix voltage in a point (used by IonFlow specialization)
Reimplemented in IonFlow.
Definition at line 1088 of file Flow1D.cpp.
|
virtual |
Retrieve flag indicating whether electric field is solved or not (used by IonFlow specialization)
Reimplemented in IonFlow.
Definition at line 1094 of file Flow1D.cpp.
|
inline |
|
inline |
|
inline |
void setBoundaryEmissivities | ( | double | e_left, |
double | e_right | ||
) |
Set the emissivities for the boundary values.
Reads the emissivities for the left and right boundary values in the radiative term and writes them into the variables, which are used for the calculation.
Definition at line 1100 of file Flow1D.cpp.
|
inline |
|
inline |
void fixTemperature | ( | size_t | j = npos | ) |
Specify that the the temperature should be held fixed at point j
.
The converse of this method is enableEnergyEqn().
j | Point at which to specify a fixed temperature. npos means all points. |
Definition at line 1114 of file Flow1D.cpp.
double leftControlPointTemperature | ( | ) | const |
Returns the temperature at the left control point.
Definition at line 1147 of file Flow1D.cpp.
double leftControlPointCoordinate | ( | ) | const |
Returns the z-coordinate of the left control point.
Definition at line 1162 of file Flow1D.cpp.
void setLeftControlPointTemperature | ( | double | temperature | ) |
Sets the temperature of the left control point.
Definition at line 1177 of file Flow1D.cpp.
void setLeftControlPointCoordinate | ( | double | z_left | ) |
Sets the coordinate of the left control point.
Definition at line 1192 of file Flow1D.cpp.
double rightControlPointTemperature | ( | ) | const |
Returns the temperature at the right control point.
Definition at line 1202 of file Flow1D.cpp.
double rightControlPointCoordinate | ( | ) | const |
Returns the z-coordinate of the right control point.
Definition at line 1217 of file Flow1D.cpp.
void setRightControlPointTemperature | ( | double | temperature | ) |
Sets the temperature of the right control point.
Definition at line 1232 of file Flow1D.cpp.
void setRightControlPointCoordinate | ( | double | z_right | ) |
Sets the coordinate of the right control point.
Definition at line 1247 of file Flow1D.cpp.
void enableTwoPointControl | ( | bool | twoPointControl | ) |
Sets the status of the two-point control.
Definition at line 1257 of file Flow1D.cpp.
|
inline |
|
inline |
|
overridevirtual |
Change the grid size. Called after grid refinement.
Reimplemented from Domain1D.
Reimplemented in IonFlow.
Definition at line 159 of file Flow1D.cpp.
void setGas | ( | const double * | x, |
size_t | j | ||
) |
Set the gas object state to be consistent with the solution at point j.
Definition at line 238 of file Flow1D.cpp.
void setGasAtMidpoint | ( | const double * | x, |
size_t | j | ||
) |
Set the gas state to be consistent with the solution at the midpoint between j and j + 1.
Definition at line 246 of file Flow1D.cpp.
|
inline |
|
inline |
|
inline |
|
inline |
|
overridevirtual |
Evaluate the residual functions for axisymmetric stagnation flow.
If jGlobal == npos, the residual function is evaluated at all grid points. Otherwise, the residual function is only evaluated at grid points j-1, j, and j+1. This option is used to efficiently evaluate the Jacobian numerically.
These residuals at all the boundary grid points are evaluated using a default boundary condition that may be modified by a boundary object that is attached to the domain. The boundary object connected will modify these equations by subtracting the boundary object's values for V, T, mdot, etc. As a result, these residual equations will force the solution variables to the values of the connected boundary object.
jGlobal | Global grid point at which to update the residual | |
[in] | xGlobal | Global state vector |
[out] | rsdGlobal | Global residual vector |
[out] | diagGlobal | Global boolean mask indicating whether each solution component has a time derivative (1) or not (0). |
[in] | rdt | Reciprocal of the timestep (rdt=0 implies steady-state.) |
Reimplemented from Domain1D.
Reimplemented in StFlow.
Definition at line 305 of file Flow1D.cpp.
|
inline |
|
inline |
|
overrideprotectedvirtual |
|
overrideprotectedvirtual |
|
inlineprotected |
Update the thermodynamic properties from point j0 to point j1 (inclusive), based on solution x.
The gas state is set to be consistent with the solution at the points from j0 to j1.
Properties that are computed and cached are:
|
protectedvirtual |
Update the transport properties at grid points in the range from j0
to j1
, based on solution x
.
Evaluates the solution at the midpoint between j
and j + 1
to compute the transport properties. For example, the viscosity at element j
is the viscosity evaluated at the midpoint between j
and j + 1
.
Reimplemented in IonFlow.
Definition at line 368 of file Flow1D.cpp.
|
protectedvirtual |
Update the diffusive mass fluxes.
Reimplemented in IonFlow.
Definition at line 416 of file Flow1D.cpp.
|
protectedvirtual |
Update the properties (thermo, transport, and diffusion flux).
This function is called in eval after the points which need to be updated are defined.
Definition at line 344 of file Flow1D.cpp.
|
protected |
Computes the radiative heat loss vector over points jmin to jmax and stores the data in the qdotRadiation variable.
The simple radiation model used was established by Liu and Rogg [23]. This model considers the radiation of CO2 and H2O.
This model uses the optically thin limit and the gray-gas approximation to simply calculate a volume specified heat flux out of the Planck absorption coefficients, the boundary emissivities and the temperature. Polynomial lines calculate the species Planck coefficients for H2O and CO2. The data for the lines are taken from the RADCAL program [9]. The coefficients for the polynomials are taken from TNF Workshop material.
Definition at line 462 of file Flow1D.cpp.
|
protectedvirtual |
Evaluate the continuity equation residual.
\frac{d(\rho u)}{dz} + 2\rho V = 0
Axisymmetric flame: The continuity equation propagates information from right-to-left. The \rho u at point 0 is dependent on \rho u at point 1, but not on \dot{m} from the inlet.
Freely-propagating flame: The continuity equation propagates information away from a fixed temperature point that is set in the domain.
Unstrained flame: A specified mass flux; the main example being burner-stabilized flames.
The default boundary condition for the continuity equation is ( u = 0 ) at the right boundary. Because the equation is a first order equation, only one boundary condition is needed.
[in] | x | Local domain state vector, includes variables like temperature, density, etc. |
[out] | rsd | Local domain residual vector that stores the continuity equation residuals. |
[out] | diag | Local domain diagonal matrix that controls whether an entry has a time-derivative (used by the solver). |
[in] | rdt | Reciprocal of the timestep. |
[in] | jmin | The index for the starting point in the local domain grid. |
[in] | jmax | The index for the ending point in the local domain grid. |
Definition at line 509 of file Flow1D.cpp.
|
protectedvirtual |
Evaluate the momentum equation residual.
\rho u \frac{dV}{dz} + \rho V^2 = \frac{d}{dz}\left( \mu \frac{dV}{dz} \right) - \Lambda
The radial momentum equation is used for axisymmetric flows, and incorporates terms for time and spatial variations of radial velocity ( V ). The default boundary condition is zero radial velocity ( V ) at the left and right boundary.
For argument explanation, see evalContinuity().
Definition at line 566 of file Flow1D.cpp.
|
protectedvirtual |
Evaluate the lambda equation residual.
\frac{d\Lambda}{dz} = 0
The lambda equation serves as an eigenvalue that allows the momentum equation and continuity equations to be simultaneously satisfied in axisymmetric flows. The lambda equation propagates information from left-to-right. The default boundary condition is \Lambda = 0 at the left boundary. The equation is first order and so only one boundary condition is needed.
For argument explanation, see evalContinuity().
Definition at line 602 of file Flow1D.cpp.
|
protectedvirtual |
Evaluate the energy equation residual.
\rho c_p u \frac{dT}{dz} = \frac{d}{dz}\left( \lambda \frac{dT}{dz} \right) - \sum_k h_kW_k\dot{\omega}_k - \sum_k j_k \frac{dh_k}{dz}
The energy equation includes contributions from chemical reactions and diffusion. Default is zero temperature ( T ) at the left and right boundaries. These boundary values are updated by the specific boundary object connected to the domain.
For argument explanation, see evalContinuity().
Definition at line 645 of file Flow1D.cpp.
|
protectedvirtual |
Evaluate the species equations' residuals.
\rho u \frac{dY_k}{dz} + \frac{dj_k}{dz} = W_k\dot{\omega}_k
The species equations include terms for temporal and spatial variations of species mass fractions ( Y_k ). The default boundary condition is zero flux for species at the left and right boundary.
For argument explanation, see evalContinuity().
Reimplemented in IonFlow.
Definition at line 727 of file Flow1D.cpp.
|
protectedvirtual |
Evaluate the electric field equation residual to be zero everywhere.
The electric field equation is implemented in the IonFlow class. The default boundary condition is zero electric field ( E ) at the boundary, and E is zero within the domain.
For argument explanation, see evalContinuity().
Reimplemented in IonFlow.
Definition at line 767 of file Flow1D.cpp.
|
protectedvirtual |
Alternate version of evalContinuity with legacy signature.
Implemented by StFlow; included here to prevent compiler warnings about shadowed virtual functions.
Reimplemented in StFlow.
Definition at line 776 of file Flow1D.cpp.
|
protectedvirtual |
Evaluate the oxidizer axial velocity equation residual.
The function calculates the oxidizer axial velocity equation as
\frac{dU_{o}}{dz} = 0
This equation serves as a dummy equation that is used only in the context of two-point flame control, and serves as the way for two interior control points to be specified while maintaining block tridiagonal structure. The default boundary condition is U_o = 0 at the right and zero flux at the left boundary.
For argument explanation, see evalContinuity().
Definition at line 687 of file Flow1D.cpp.
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
|
inlineprotected |
Calculates the spatial derivative of velocity V with respect to z at point j using upwind differencing.
For more details on the upwinding scheme, see the science reference documentation.
\frac{\partial V}{\partial z} \bigg|_{j} \approx \frac{V_{\ell} - V_{\ell-1}}{z_{\ell} - z_{\ell-1}}
Where the value of \ell is determined by the sign of the axial velocity. If the axial velocity is positive, the value of \ell is j. If the axial velocity is negative, the value of \ell is j + 1. A positive velocity means that the flow is moving left-to-right.
[in] | x | The local domain state vector. |
[in] | j | The grid point index at which the derivative is computed. |
|
inlineprotected |
Calculates the spatial derivative of the species mass fraction Y_k with respect to z for species k at point j using upwind differencing.
For details on the upwinding scheme, see dVdz().
[in] | x | The local domain state vector. |
[in] | k | The species index. |
[in] | j | The grid point index at which the derivative is computed. |
|
inlineprotected |
Calculates the spatial derivative of temperature T with respect to z at point j using upwind differencing.
For details on the upwinding scheme, see dVdz().
[in] | x | The local domain state vector. |
[in] | j | The grid point index at which the derivative is computed. |
|
inlineprotected |
Compute the shear term from the momentum equation using a central three-point differencing scheme.
The term to be discretized is:
\frac{d}{dz}\left(\mu \frac{dV}{dz}\right)
For more details on the discretization scheme used for the second derivative, see the documentation.
\frac{d}{dz}\left(\mu \frac{dV}{dz}\right) \approx \frac{\mu_{j+1/2} \frac{V_{j+1} - V_j}{z_{j+1} - z_j} - \mu_{j-1/2} \frac{V_j - V_{j-1}}{z_j - z_{j-1}}}{\frac{z_{j+1} - z_{j-1}}{2}}
[in] | x | The local domain state vector. |
[in] | j | The grid point index at which the derivative is computed. |
|
inlineprotected |
Compute the conduction term from the energy equation using a central three-point differencing scheme.
For the details about the discretization, see shear().
[in] | x | The local domain state vector. |
[in] | j | The grid point index at which the derivative is computed. |
|
inlineprotected |
Array access mapping for a 3D array stored in a 1D vector.
Used for accessing data in the m_multidiff member variable.
[in] | k | First species index. |
[in] | j | The grid point index. |
[in] | m | The second species index. |
|
protectedvirtual |
Compute the spatial derivative of species specific molar enthalpies using upwind differencing.
Updates all species molar enthalpies for all species at point j. Updates the m_dhk_dz 2D array.
For details on the upwinding scheme, see dVdz().
[in] | x | The local domain state vector. |
[in] | j | The index at which the derivative is computed. |
Definition at line 1138 of file Flow1D.cpp.
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
Coefficient used in diffusion calculations for each species at each grid point.
The value stored is different depending on the transport model (multicomponent versus mixture averaged) and flux gradient basis (mass or molar). Vector size is m_nsp × m_points, where m_diff[k + j*m_nsp]
contains the value for species k
at point j
.
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
Indices within the ThermoPhase of the radiating species.
First index is for CO2, second is for H2O.
|
protected |
For each point in the domain, true
if energy equation is solved or false
if temperature is held constant.
|
protected |
|
protected |
Determines whether diffusive fluxes are computed using gradients of mass fraction or mole fraction.
|
protected |
|
protected |
Determines whether radiative heat loss is calculated.
|
protected |
Determines whether the viscosity term in the momentum equation is calculated.
|
protected |
Flag that is true
for freely propagating flames anchored by a temperature fixed point.
|
protected |
Flag that is true
for counterflow configurations that use the pressure eigenvalue \Lambda in the radial momentum equation.
|
protected |
|
protected |
|
protected |
Fixed values of the temperature at each grid point that are used when solving with the energy equation disabled.
Values are interpolated from profiles specified with the setFixedTempProfile method as part of _finalize().
|
protected |
Relative coordinates used to specify a fixed temperature profile.
0 corresponds to the left edge of the domain and 1 corresponds to the right edge of the domain. Length is the same as the m_tfix array.
|
protected |
Fixed temperature values at the relative coordinates specified in m_zfix.
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
|
protected |
double m_zfixed = Undef |
double m_tfixed = -1.0 |
|
private |