Cantera
2.4.0
|
This class models the ion transportation in a flame. More...
#include <IonFlow.h>
Public Member Functions | |
IonFlow (IdealGasPhase *ph=0, size_t nsp=1, size_t points=1) | |
virtual void | setSolvingStage (const size_t phase) |
set the solving stage More... | |
virtual void | resize (size_t components, size_t points) |
Change the grid size. Called after grid refinement. More... | |
virtual void | _finalize (const double *x) |
In some cases, a domain may need to set parameters that depend on the initial solution estimate. More... | |
void | solveElectricField (size_t j=npos) |
set to solve electric field on a point More... | |
void | fixElectricField (size_t j=npos) |
set to fix voltage on a point More... | |
bool | doElectricField (size_t j) |
void | setElectronTransport (vector_fp &tfix, vector_fp &diff_e, vector_fp &mobi_e) |
Sometimes it is desired to carry out the simulation using a specified electron transport profile, rather than assuming it as a constant (0.4). More... | |
Public Member Functions inherited from StFlow | |
StFlow (IdealGasPhase *ph=0, size_t nsp=1, size_t points=1) | |
Create a new flow domain. More... | |
virtual std::string | componentName (size_t n) const |
Name of the nth component. May be overloaded. More... | |
virtual size_t | componentIndex (const std::string &name) const |
virtual void | showSolution (const doublereal *x) |
Print the solution. More... | |
virtual XML_Node & | save (XML_Node &o, const doublereal *const sol) |
Save the current solution for this domain into an XML_Node. More... | |
virtual void | restore (const XML_Node &dom, doublereal *soln, int loglevel) |
Restore the solution for this domain from an XML_Node. More... | |
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. More... | |
void | setAxisymmetricFlow () |
Set flow configuration for axisymmetric counterflow or burner-stabilized flames, using specified inlet mass fluxes. More... | |
virtual std::string | flowType () |
Return the type of flow domain being represented, either "Free Flame" or "Axisymmetric Stagnation". More... | |
void | solveEnergyEqn (size_t j=npos) |
void | enableRadiation (bool doRadiation) |
Turn radiation on / off. More... | |
bool | radiationEnabled () const |
Returns true if the radiation term in the energy equation is enabled. More... | |
void | setBoundaryEmissivities (doublereal e_left, doublereal e_right) |
Set the emissivities for the boundary values. More... | |
void | fixTemperature (size_t j=npos) |
bool | doEnergy (size_t j) |
virtual void | setFixedPoint (int j0, doublereal t0) |
void | setGas (const doublereal *x, size_t j) |
Set the gas object state to be consistent with the solution at point j. More... | |
void | setGasAtMidpoint (const doublereal *x, size_t j) |
Set the gas state to be consistent with the solution at the midpoint between j and j + 1. More... | |
doublereal | density (size_t j) const |
virtual bool | fixed_mdot () |
void | setViscosityFlag (bool dovisc) |
virtual void | eval (size_t j, doublereal *x, doublereal *r, integer *mask, doublereal rdt) |
virtual void | evalRightBoundary (double *x, double *res, int *diag, double rdt) |
Evaluate all residual components at the right boundary. More... | |
virtual void | evalContinuity (size_t j, double *x, double *r, int *diag, double rdt) |
Evaluate the residual corresponding to the continuity equation at all interior grid points. More... | |
size_t | leftExcessSpecies () const |
Index of the species on the left boundary with the largest mass fraction. More... | |
size_t | rightExcessSpecies () const |
Index of the species on the right boundary with the largest mass fraction. More... | |
virtual void | setupGrid (size_t n, const doublereal *z) |
called to set up initial grid, and after grid refinement More... | |
virtual void | resetBadValues (double *xg) |
thermo_t & | phase () |
Kinetics & | kinetics () |
void | setThermo (IdealGasPhase &th) |
Set the thermo manager. More... | |
void | setKinetics (Kinetics &kin) |
Set the kinetics manager. The kinetics manager must. More... | |
void | setTransport (Transport &trans) |
set the transport manager More... | |
void | enableSoret (bool withSoret) |
Enable thermal diffusion, also known as Soret diffusion. More... | |
bool | withSoret () const |
void | setPressure (doublereal p) |
Set the pressure. More... | |
doublereal | pressure () const |
The current pressure [Pa]. More... | |
virtual void | _getInitialSoln (double *x) |
Write the initial solution estimate into array x. More... | |
void | setFixedTempProfile (vector_fp &zfixed, vector_fp &tfixed) |
Sometimes it is desired to carry out the simulation using a specified temperature profile, rather than computing it by solving the energy equation. More... | |
void | setTemperature (size_t j, doublereal t) |
doublereal | T_fixed (size_t j) const |
The fixed temperature value at point j. More... | |
Public Member Functions inherited from Domain1D | |
Domain1D (size_t nv=1, size_t points=1, double time=0.0) | |
Constructor. More... | |
Domain1D (const Domain1D &)=delete | |
Domain1D & | operator= (const Domain1D &)=delete |
int | domainType () |
Domain type flag. More... | |
size_t | domainIndex () |
The left-to-right location of this domain. More... | |
bool | isConnector () |
True if the domain is a connector domain. More... | |
const OneDim & | container () const |
The container holding this domain. More... | |
void | setContainer (OneDim *c, size_t index) |
Specify the container object for this domain, and the position of this domain in the list. More... | |
void | setBandwidth (int bw=-1) |
Set the Jacobian bandwidth. See the discussion of method bandwidth(). More... | |
size_t | bandwidth () |
Set the Jacobian bandwidth for this domain. More... | |
virtual void | init () |
virtual void | setInitialState (doublereal *xlocal=0) |
virtual void | setState (size_t point, const doublereal *state, doublereal *x) |
Refiner & | refiner () |
Return a reference to the grid refiner. More... | |
size_t | nComponents () const |
Number of components at each grid point. More... | |
void | checkComponentIndex (size_t n) const |
Check that the specified component index is in range. More... | |
void | checkComponentArraySize (size_t nn) const |
Check that an array size is at least nComponents(). More... | |
size_t | nPoints () const |
Number of grid points in this domain. More... | |
void | checkPointIndex (size_t n) const |
Check that the specified point index is in range. More... | |
void | checkPointArraySize (size_t nn) const |
Check that an array size is at least nPoints(). More... | |
void | setComponentName (size_t n, const std::string &name) |
size_t | componentIndex (const std::string &name) const |
index of component with name name. More... | |
void | setBounds (size_t n, doublereal lower, doublereal upper) |
void | setTransientTolerances (doublereal rtol, doublereal atol, size_t n=npos) |
Set tolerances for time-stepping mode. More... | |
void | setSteadyTolerances (doublereal rtol, doublereal atol, size_t n=npos) |
Set tolerances for steady-state mode. More... | |
doublereal | rtol (size_t n) |
Relative tolerance of the nth component. More... | |
doublereal | atol (size_t n) |
Absolute tolerance of the nth component. More... | |
double | steady_rtol (size_t n) |
Steady relative tolerance of the nth component. More... | |
double | steady_atol (size_t n) |
Steady absolute tolerance of the nth component. More... | |
double | transient_rtol (size_t n) |
Transient relative tolerance of the nth component. More... | |
double | transient_atol (size_t n) |
Transient absolute tolerance of the nth component. More... | |
doublereal | upperBound (size_t n) const |
Upper bound on the nth component. More... | |
doublereal | lowerBound (size_t n) const |
Lower bound on the nth component. More... | |
void | initTimeInteg (doublereal dt, const doublereal *x0) |
Prepare to do time stepping with time step dt. More... | |
void | setSteadyMode () |
Prepare to solve the steady-state problem. More... | |
bool | steady () |
True if in steady-state mode. More... | |
bool | transient () |
True if not in steady-state mode. More... | |
void | needJacUpdate () |
size_t | index (size_t n, size_t j) const |
doublereal | value (const doublereal *x, size_t n, size_t j) const |
virtual void | setJac (MultiJac *jac) |
size_t | size () const |
void | locate () |
Find the index of the first grid point in this domain, and the start of its variables in the global solution vector. More... | |
virtual size_t | loc (size_t j=0) const |
Location of the start of the local solution vector in the global solution vector,. More... | |
size_t | firstPoint () const |
The index of the first (i.e., left-most) grid point belonging to this domain. More... | |
size_t | lastPoint () const |
The index of the last (i.e., right-most) grid point belonging to this domain. More... | |
void | linkLeft (Domain1D *left) |
Set the left neighbor to domain 'left. More... | |
void | linkRight (Domain1D *right) |
Set the right neighbor to domain 'right.'. More... | |
void | append (Domain1D *right) |
Append domain 'right' to this one, and update all links. More... | |
Domain1D * | left () const |
Return a pointer to the left neighbor. More... | |
Domain1D * | right () const |
Return a pointer to the right neighbor. More... | |
double | prevSoln (size_t n, size_t j) const |
Value of component n at point j in the previous solution. More... | |
void | setID (const std::string &s) |
Specify an identifying tag for this domain. More... | |
std::string | id () const |
virtual void | showSolution_s (std::ostream &s, const doublereal *x) |
doublereal | z (size_t jlocal) const |
doublereal | zmin () const |
doublereal | zmax () const |
void | setProfile (const std::string &name, double *values, double *soln) |
vector_fp & | grid () |
const vector_fp & | grid () const |
doublereal | grid (size_t point) const |
virtual doublereal | initialValue (size_t n, size_t j) |
Initial value of solution component n at grid point j. More... | |
void | forceFullUpdate (bool update) |
In some cases, for computational efficiency some properties (e.g. More... | |
Protected Member Functions | |
virtual void | evalResidual (double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) |
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 . More... | |
virtual void | updateDiffFluxes (const double *x, size_t j0, size_t j1) |
Update the diffusive mass fluxes. More... | |
virtual void | frozenIonMethod (const double *x, size_t j0, size_t j1) |
Solving phase one: the fluxes of charged species are turned off. More... | |
virtual void | electricFieldMethod (const double *x, size_t j0, size_t j1) |
Solving phase two: the electric field equation is added coupled by the electrical drift. More... | |
double | E (const double *x, size_t j) const |
electric field More... | |
double | dEdz (const double *x, size_t j) const |
double | ND (const double *x, size_t k, size_t j) const |
number density More... | |
double | rho_e (double *x, size_t j) const |
total charge density More... | |
Protected Member Functions inherited from StFlow | |
doublereal | wdot (size_t k, size_t j) const |
void | getWdot (doublereal *x, size_t j) |
Write the net production rates at point j into array m_wdot More... | |
virtual void | updateProperties (size_t jg, double *x, size_t jmin, size_t jmax) |
Update the properties (thermo, transport, and diffusion flux). More... | |
void | updateThermo (const doublereal *x, size_t j0, size_t j1) |
Update the thermodynamic properties from point j0 to point j1 (inclusive), based on solution x. More... | |
doublereal | shear (const doublereal *x, size_t j) const |
doublereal | divHeatFlux (const doublereal *x, size_t j) const |
size_t | mindex (size_t k, size_t j, size_t m) |
doublereal | T (const doublereal *x, size_t j) const |
doublereal & | T (doublereal *x, size_t j) |
doublereal | T_prev (size_t j) const |
doublereal | rho_u (const doublereal *x, size_t j) const |
doublereal | u (const doublereal *x, size_t j) const |
doublereal | V (const doublereal *x, size_t j) const |
doublereal | V_prev (size_t j) const |
doublereal | lambda (const doublereal *x, size_t j) const |
doublereal | Y (const doublereal *x, size_t k, size_t j) const |
doublereal & | Y (doublereal *x, size_t k, size_t j) |
doublereal | Y_prev (size_t k, size_t j) const |
doublereal | X (const doublereal *x, size_t k, size_t j) const |
doublereal | flux (size_t k, size_t j) const |
doublereal | dVdz (const doublereal *x, size_t j) const |
doublereal | dYdz (const doublereal *x, size_t k, size_t j) const |
doublereal | dTdz (const doublereal *x, size_t j) const |
Protected Attributes | |
std::vector< bool > | m_do_electric_field |
flag for solving electric field or not More... | |
bool | m_import_electron_transport |
flag for importing transport of electron More... | |
vector_int | m_speciesCharge |
electrical properties More... | |
std::vector< size_t > | m_kCharge |
index of species with charges More... | |
std::vector< size_t > | m_kNeutral |
index of neutral species More... | |
vector_fp | m_mobi_e_fix |
coefficients of polynomial fitting of fixed electron transport profile More... | |
vector_fp | m_diff_e_fix |
vector_fp | m_mobility |
mobility More... | |
int | m_stage |
solving stage More... | |
size_t | m_kElectron |
index of electron More... | |
Protected Attributes inherited from StFlow | |
doublereal | m_press |
vector_fp | m_dz |
vector_fp | m_rho |
vector_fp | m_wtm |
vector_fp | m_wt |
vector_fp | m_cp |
vector_fp | m_visc |
vector_fp | m_tcon |
vector_fp | m_diff |
vector_fp | m_multidiff |
Array2D | m_dthermal |
Array2D | m_flux |
Array2D | m_wdot |
size_t | m_nsp |
IdealGasPhase * | m_thermo |
Kinetics * | m_kin |
Transport * | m_trans |
doublereal | m_epsilon_left |
doublereal | m_epsilon_right |
std::vector< size_t > | m_kRadiating |
Indices within the ThermoPhase of the radiating species. More... | |
std::vector< bool > | m_do_energy |
bool | m_do_soret |
std::vector< bool > | m_do_species |
bool | m_do_multicomponent |
bool | m_do_radiation |
flag for the radiative heat loss More... | |
vector_fp | m_qdotRadiation |
radiative heat loss vector More... | |
vector_fp | m_fixedtemp |
vector_fp | m_zfix |
vector_fp | m_tfix |
size_t | m_kExcessLeft |
Index of species with a large mass fraction at each boundary, for which the mass fraction may be calculated as 1 minus the sum of the other mass fractions. More... | |
size_t | m_kExcessRight |
bool | m_dovisc |
Protected Attributes inherited from Domain1D | |
doublereal | m_rdt |
size_t | m_nv |
size_t | m_points |
vector_fp | m_slast |
vector_fp | m_max |
vector_fp | m_min |
vector_fp | m_rtol_ss |
vector_fp | m_rtol_ts |
vector_fp | m_atol_ss |
vector_fp | m_atol_ts |
vector_fp | m_z |
OneDim * | m_container |
size_t | m_index |
int | m_type |
size_t | m_iloc |
Starting location within the solution vector for unknowns that correspond to this domain. More... | |
size_t | m_jstart |
Domain1D * | m_left |
Domain1D * | m_right |
std::string | m_id |
Identity tag for the domain. More... | |
std::string | m_desc |
std::unique_ptr< Refiner > | m_refiner |
std::vector< std::string > | m_name |
int | m_bw |
bool | m_force_full_update |
Additional Inherited Members | |
Public Attributes inherited from StFlow | |
double | m_zfixed |
Location of the point where temperature is fixed. More... | |
double | m_tfixed |
Temperature at the point used to fix the flame location. More... | |
This class models the ion transportation in a flame.
There are three stages of the simulation.
The first stage turns off the diffusion of ions due to the fast diffusion rate of electron without internal electric forces (ambi- polar diffusion effect).
The second stage evaluates drift flux from electric field calculated from Poisson's equation, which is solved together with other equations. Poisson's equation is coupled because the total charge densities depends on the species' concentration. Reference: Pederson, Timothy, and R. C. Brown. "Simulation of electric field effects in premixed methane flames." Combustion and Flames 94.4(1993): 433-448.
|
virtual |
|
virtual |
Change the grid size. Called after grid refinement.
Reimplemented from StFlow.
Definition at line 59 of file IonFlow.cpp.
References IonFlow::m_do_electric_field, IonFlow::m_mobility, and StFlow::resize().
|
virtual |
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 StFlow.
Definition at line 272 of file IonFlow.cpp.
References StFlow::_finalize(), IonFlow::m_do_electric_field, and IonFlow::solveElectricField().
void solveElectricField | ( | size_t | j = npos | ) |
set to solve electric field on a point
Definition at line 205 of file IonFlow.cpp.
References IonFlow::m_do_electric_field, Domain1D::needJacUpdate(), and Cantera::npos.
Referenced by IonFlow::_finalize().
void fixElectricField | ( | size_t | j = npos | ) |
set to fix voltage on a point
Definition at line 230 of file IonFlow.cpp.
References IonFlow::m_do_electric_field, Domain1D::needJacUpdate(), and Cantera::npos.
Sometimes it is desired to carry out the simulation using a specified electron transport profile, rather than assuming it as a constant (0.4).
Reference: Bisetti, Fabrizio, and Mbark El Morsli. "Calculation and analysis of the mobility and diffusion coefficient of thermal electrons in methane/air premixed flames." Combustion and flame 159.12 (2012): 3518-3521. If in the future the class GasTranport is improved, this method may be discard. This method specifies this profile.
Definition at line 255 of file IonFlow.cpp.
References IonFlow::m_import_electron_transport, IonFlow::m_mobi_e_fix, and Cantera::polyfit().
|
protectedvirtual |
This function overloads the original function. The residual function of electric field is added.
Reimplemented from StFlow.
Definition at line 170 of file IonFlow.cpp.
References IonFlow::E(), Cantera::epsilon_0, StFlow::evalResidual(), IonFlow::m_kCharge, IonFlow::m_stage, and IonFlow::rho_e().
|
protectedvirtual |
Update the transport properties at grid points in the range from j0
to j1
, based on solution x
.
Reimplemented from StFlow.
Definition at line 66 of file IonFlow.cpp.
References Transport::getMobilities(), IonFlow::m_import_electron_transport, IonFlow::m_kElectron, IonFlow::m_mobi_e_fix, IonFlow::m_mobility, Cantera::poly5(), StFlow::setGasAtMidpoint(), Phase::temperature(), and StFlow::updateTransport().
|
protectedvirtual |
Update the diffusive mass fluxes.
Reimplemented from StFlow.
Definition at line 81 of file IonFlow.cpp.
References IonFlow::electricFieldMethod(), IonFlow::frozenIonMethod(), and IonFlow::m_stage.
|
protectedvirtual |
Solving phase one: the fluxes of charged species are turned off.
Definition at line 91 of file IonFlow.cpp.
References IonFlow::m_kCharge, and IonFlow::m_kNeutral.
Referenced by IonFlow::updateDiffFluxes().
|
protectedvirtual |
Solving phase two: the electric field equation is added coupled by the electrical drift.
Definition at line 118 of file IonFlow.cpp.
References IonFlow::E(), IonFlow::m_kCharge, IonFlow::m_kNeutral, IonFlow::m_mobility, and IonFlow::m_speciesCharge.
Referenced by IonFlow::updateDiffFluxes().
|
inlineprotected |
electric field
Definition at line 107 of file IonFlow.h.
Referenced by IonFlow::electricFieldMethod(), and IonFlow::evalResidual().
|
inlineprotected |
number density
Definition at line 116 of file IonFlow.h.
References Cantera::Avogadro.
Referenced by IonFlow::rho_e().
|
inlineprotected |
total charge density
Definition at line 121 of file IonFlow.h.
References IonFlow::m_kCharge, IonFlow::m_speciesCharge, and IonFlow::ND().
Referenced by IonFlow::evalResidual().
|
protected |
flag for solving electric field or not
Definition at line 79 of file IonFlow.h.
Referenced by IonFlow::_finalize(), IonFlow::fixElectricField(), IonFlow::resize(), and IonFlow::solveElectricField().
|
protected |
flag for importing transport of electron
Definition at line 82 of file IonFlow.h.
Referenced by IonFlow::setElectronTransport(), and IonFlow::updateTransport().
|
protected |
electrical properties
Definition at line 85 of file IonFlow.h.
Referenced by IonFlow::electricFieldMethod(), and IonFlow::rho_e().
|
protected |
index of species with charges
Definition at line 88 of file IonFlow.h.
Referenced by IonFlow::electricFieldMethod(), IonFlow::evalResidual(), IonFlow::frozenIonMethod(), and IonFlow::rho_e().
|
protected |
index of neutral species
Definition at line 91 of file IonFlow.h.
Referenced by IonFlow::electricFieldMethod(), and IonFlow::frozenIonMethod().
|
protected |
coefficients of polynomial fitting of fixed electron transport profile
Definition at line 94 of file IonFlow.h.
Referenced by IonFlow::setElectronTransport(), and IonFlow::updateTransport().
|
protected |
mobility
Definition at line 98 of file IonFlow.h.
Referenced by IonFlow::electricFieldMethod(), IonFlow::resize(), and IonFlow::updateTransport().
|
protected |
solving stage
Definition at line 101 of file IonFlow.h.
Referenced by IonFlow::evalResidual(), IonFlow::setSolvingStage(), and IonFlow::updateDiffFluxes().
|
protected |
index of electron
Definition at line 104 of file IonFlow.h.
Referenced by IonFlow::updateTransport().