Cantera  3.1.0a1
StFlow.h
Go to the documentation of this file.
1 //! @file StFlow.h
2 
3 // This file is part of Cantera. See License.txt in the top-level directory or
4 // at https://cantera.org/license.txt for license and copyright information.
5 
6 #ifndef CT_STFLOW_H
7 #define CT_STFLOW_H
8 
9 #include "Domain1D.h"
10 #include "cantera/base/Array.h"
11 #include "cantera/base/Solution.h"
14 
15 namespace Cantera
16 {
17 
18 //------------------------------------------
19 // constants
20 //------------------------------------------
21 
22 //! Offsets of solution components in the 1D solution array.
23 enum offset
24 {
25  c_offset_U //! axial velocity
26  , c_offset_V //! strain rate
27  , c_offset_T //! temperature
28  , c_offset_L //! (1/r)dP/dr
29  , c_offset_E //! electric field equation
30  , c_offset_Y //! mass fractions
31 };
32 
33 class Transport;
34 
35 //! @defgroup flowGroup Flow Domains
36 //! One-dimensional flow domains.
37 //! @ingroup onedGroup
38 
39 /**
40  * This class represents 1D flow domains that satisfy the one-dimensional
41  * similarity solution for chemically-reacting, axisymmetric flows.
42  * @ingroup flowGroup
43  */
44 class StFlow : public Domain1D
45 {
46 public:
47  //--------------------------------
48  // construction and destruction
49  //--------------------------------
50 
51  //! Create a new flow domain.
52  //! @param ph Object representing the gas phase. This object will be used
53  //! to evaluate all thermodynamic, kinetic, and transport properties.
54  //! @param nsp Number of species.
55  //! @param points Initial number of grid points
56  StFlow(ThermoPhase* ph = 0, size_t nsp = 1, size_t points = 1);
57 
58  //! Delegating constructor
59  StFlow(shared_ptr<ThermoPhase> th, size_t nsp = 1, size_t points = 1);
60 
61  //! Create a new flow domain.
62  //! @param sol Solution object used to evaluate all thermodynamic, kinetic, and
63  //! transport properties
64  //! @param id name of flow domain
65  //! @param points initial number of grid points
66  StFlow(shared_ptr<Solution> sol, const string& id="", size_t points=1);
67 
68  ~StFlow();
69 
70  string domainType() const override;
71 
72  //! @name Problem Specification
73  //! @{
74 
75  void setupGrid(size_t n, const double* z) override;
76 
77  void resetBadValues(double* xg) override;
78 
79  ThermoPhase& phase() {
80  return *m_thermo;
81  }
82 
83  Kinetics& kinetics() {
84  return *m_kin;
85  }
86 
87  void setKinetics(shared_ptr<Kinetics> kin) override;
88 
89  void setTransport(shared_ptr<Transport> trans) override;
90 
91  //! Set the transport model
92  //! @since New in %Cantera 3.0.
93  void setTransportModel(const string& trans);
94 
95  //! Retrieve transport model
96  //! @since New in %Cantera 3.0.
97  string transportModel() const;
98 
99  //! Enable thermal diffusion, also known as Soret diffusion.
100  //! Requires that multicomponent transport properties be
101  //! enabled to carry out calculations.
102  void enableSoret(bool withSoret) {
103  m_do_soret = withSoret;
104  }
105  bool withSoret() const {
106  return m_do_soret;
107  }
108 
109  //! Set the pressure. Since the flow equations are for the limit of small
110  //! Mach number, the pressure is very nearly constant throughout the flow.
111  void setPressure(double p) {
112  m_press = p;
113  }
114 
115  //! The current pressure [Pa].
116  double pressure() const {
117  return m_press;
118  }
119 
120  //! Write the initial solution estimate into array x.
121  void _getInitialSoln(double* x) override;
122 
123  void _finalize(const double* x) override;
124 
125  //! Sometimes it is desired to carry out the simulation using a specified
126  //! temperature profile, rather than computing it by solving the energy
127  //! equation. This method specifies this profile.
128  void setFixedTempProfile(vector<double>& zfixed, vector<double>& tfixed) {
129  m_zfix = zfixed;
130  m_tfix = tfixed;
131  }
132 
133  /**
134  * Set the temperature fixed point at grid point j, and disable the energy
135  * equation so that the solution will be held to this value.
136  */
137  void setTemperature(size_t j, double t) {
138  m_fixedtemp[j] = t;
139  m_do_energy[j] = false;
140  }
141 
142  //! The fixed temperature value at point j.
143  double T_fixed(size_t j) const {
144  return m_fixedtemp[j];
145  }
146 
147  //! @}
148 
149  string componentName(size_t n) const override;
150 
151  size_t componentIndex(const string& name) const override;
152 
153  //! Returns true if the specified component is an active part of the solver state
154  virtual bool componentActive(size_t n) const;
155 
156  //! Print the solution.
157  void show(const double* x) override;
158 
159  shared_ptr<SolutionArray> asArray(const double* soln) const override;
160  void fromArray(SolutionArray& arr, double* soln) override;
161 
162  //! Set flow configuration for freely-propagating flames, using an internal point
163  //! with a fixed temperature as the condition to determine the inlet mass flux.
164  void setFreeFlow() {
165  m_dovisc = false;
166  m_isFree = true;
167  m_usesLambda = false;
168  }
169 
170  //! Set flow configuration for axisymmetric counterflow flames, using specified
171  //! inlet mass fluxes.
173  m_dovisc = true;
174  m_isFree = false;
175  m_usesLambda = true;
176  }
177 
178  //! Set flow configuration for burner-stabilized flames, using specified inlet mass
179  //! fluxes.
181  m_dovisc = false;
182  m_isFree = false;
183  m_usesLambda = false;
184  }
185 
186  void solveEnergyEqn(size_t j=npos);
187 
188  //! Get the solving stage (used by IonFlow specialization)
189  //! @since New in %Cantera 3.0
190  virtual size_t getSolvingStage() const;
191 
192  //! Solving stage mode for handling ionized species (used by IonFlow specialization)
193  //! - @c stage=1: the fluxes of charged species are set to zero
194  //! - @c stage=2: the electric field equation is solved, and the drift flux for
195  //! ionized species is evaluated
196  virtual void setSolvingStage(const size_t stage);
197 
198  //! Set to solve electric field in a point (used by IonFlow specialization)
199  virtual void solveElectricField(size_t j=npos);
200 
201  //! Set to fix voltage in a point (used by IonFlow specialization)
202  virtual void fixElectricField(size_t j=npos);
203 
204  //! Retrieve flag indicating whether electric field is solved or not (used by
205  //! IonFlow specialization)
206  virtual bool doElectricField(size_t j) const;
207 
208  //! Turn radiation on / off.
209  void enableRadiation(bool doRadiation) {
210  m_do_radiation = doRadiation;
211  }
212 
213  //! Returns `true` if the radiation term in the energy equation is enabled
214  bool radiationEnabled() const {
215  return m_do_radiation;
216  }
217 
218  //! Return radiative heat loss at grid point j
219  double radiativeHeatLoss(size_t j) const {
220  return m_qdotRadiation[j];
221  }
222 
223  //! Set the emissivities for the boundary values
224  /*!
225  * Reads the emissivities for the left and right boundary values in the
226  * radiative term and writes them into the variables, which are used for the
227  * calculation.
228  */
229  void setBoundaryEmissivities(double e_left, double e_right);
230 
231  //! Return emissivity at left boundary
232  double leftEmissivity() const { return m_epsilon_left; }
233 
234  //! Return emissivity at right boundary
235  double rightEmissivity() const { return m_epsilon_right; }
236 
237  void fixTemperature(size_t j=npos);
238 
239  bool doEnergy(size_t j) {
240  return m_do_energy[j];
241  }
242 
243  //! Change the grid size. Called after grid refinement.
244  void resize(size_t components, size_t points) override;
245 
246  //! Set the gas object state to be consistent with the solution at point j.
247  void setGas(const double* x, size_t j);
248 
249  //! Set the gas state to be consistent with the solution at the midpoint
250  //! between j and j + 1.
251  void setGasAtMidpoint(const double* x, size_t j);
252 
253  double density(size_t j) const {
254  return m_rho[j];
255  }
256 
257  /**
258  * Retrieve flag indicating whether flow is freely propagating.
259  * The flow is unstrained and the axial mass flow rate is not specified.
260  * For free flame propagation, the axial velocity is determined by the solver.
261  * @since New in %Cantera 3.0
262  */
263  bool isFree() const {
264  return m_isFree;
265  }
266 
267  /**
268  * Retrieve flag indicating whether flow uses radial momentum.
269  * If `true`, radial momentum equation for @f$ V @f$ as well as
270  * @f$ d\Lambda/dz = 0 @f$ are solved; if `false`, @f$ \Lambda(z) = 0 @f$ and
271  * @f$ V(z) = 0 @f$ by definition.
272  * @since New in %Cantera 3.0
273  */
274  bool isStrained() const {
275  return m_usesLambda;
276  }
277 
278  void setViscosityFlag(bool dovisc) {
279  m_dovisc = dovisc;
280  }
281 
282  /**
283  * Evaluate the residual functions for axisymmetric stagnation flow.
284  * If jGlobal == npos, the residual function is evaluated at all grid points.
285  * Otherwise, the residual function is only evaluated at grid points j-1, j,
286  * and j+1. This option is used to efficiently evaluate the Jacobian numerically.
287  *
288  * These residuals at all the boundary grid points are evaluated using a default
289  * boundary condition that may be modified by a boundary object that is attached
290  * to the domain. The boundary object connected will modify these equations by
291  * subtracting the boundary object's values for V, T, mdot, etc. As a result,
292  * these residual equations will force the solution variables to the values of
293  * the connected boundary object.
294  *
295  * @param jGlobal Global grid point at which to update the residual
296  * @param[in] xGlobal Global state vector
297  * @param[out] rsdGlobal Global residual vector
298  * @param[out] diagGlobal Global boolean mask indicating whether each solution
299  * component has a time derivative (1) or not (0).
300  * @param[in] rdt Reciprocal of the timestep (`rdt=0` implies steady-state.)
301  */
302  void eval(size_t jGlobal, double* xGlobal, double* rsdGlobal,
303  integer* diagGlobal, double rdt) override;
304 
305  //! Index of the species on the left boundary with the largest mass fraction
306  size_t leftExcessSpecies() const {
307  return m_kExcessLeft;
308  }
309 
310  //! Index of the species on the right boundary with the largest mass fraction
311  size_t rightExcessSpecies() const {
312  return m_kExcessRight;
313  }
314 
315 protected:
316  AnyMap getMeta() const override;
317  void setMeta(const AnyMap& state) override;
318 
319  double wdot(size_t k, size_t j) const {
320  return m_wdot(k,j);
321  }
322 
323  //! Update the properties (thermo, transport, and diffusion flux).
324  //! This function is called in eval after the points which need
325  //! to be updated are defined.
326  virtual void updateProperties(size_t jg, double* x, size_t jmin, size_t jmax);
327 
328  /**
329  * Computes the radiative heat loss vector over points jmin to jmax and stores
330  * the data in the qdotRadiation variable.
331  *
332  * The simple radiation model used was established by Liu and Rogg
333  * @cite liu1991. This model considers the radiation of CO2 and H2O.
334  *
335  * This model uses the optically thin limit and the gray-gas approximation to
336  * simply calculate a volume specified heat flux out of the Planck absorption
337  * coefficients, the boundary emissivities and the temperature. Polynomial lines
338  * calculate the species Planck coefficients for H2O and CO2. The data for the
339  * lines are taken from the RADCAL program @cite RADCAL.
340  * The coefficients for the polynomials are taken from
341  * [TNF Workshop](https://tnfworkshop.org/radiation/) material.
342  */
343  void computeRadiation(double* x, size_t jmin, size_t jmax);
344 
345  /**
346  * Evaluate the continuity equation residual.
347  *
348  * This function calculates the residual of the continuity equation
349  * @f[
350  * \frac{d(\rho u)}{dz} + 2\rho V = 0
351  * @f]
352  *
353  * Axisymmetric flame:
354  * The continuity equation propagates information from right-to-left.
355  * The @f$ \rho u @f$ at point 0 is dependent on @f$ \rho u @f$ at point 1,
356  * but not on @f$ \dot{m} @f$ from the inlet.
357  *
358  * Freely-propagating flame:
359  * The continuity equation propagates information away from a fixed temperature
360  * point that is set in the domain.
361  *
362  * Unstrained flame:
363  * A specified mass flux; the main example being burner-stabilized flames.
364  *
365  * The default boundary condition for the continuity equation is
366  * (@f$ u = 0 @f$) at the left and right boundary.
367  *
368  * @param[in] x Local domain state vector, includes variables like temperature,
369  * density, etc.
370  * @param[out] rsd Local domain residual vector that stores the continuity
371  * equation residuals.
372  * @param[out] diag Local domain diagonal matrix that controls whether an entry
373  * has a time-derivative (used by the solver).
374  * @param[in] rdt Reciprocal of the timestep.
375  * @param[in] jmin The index for the starting point in the local domain grid.
376  * @param[in] jmax The index for the ending point in the local domain grid.
377  */
378  virtual void evalContinuity(double* x, double* rsd, int* diag,
379  double rdt, size_t jmin, size_t jmax);
380 
381  /**
382  * Evaluate the momentum equation residual.
383  *
384  * The function calculates the radial momentum equation defined as
385  * @f[
386  * \rho u \frac{dV}{dz} + \rho V^2 =
387  * \frac{d}{dz}\left( \mu \frac{dV}{dz} \right) - \Lambda
388  * @f]
389  *
390  * The radial momentum equation is used for axisymmetric flows, and incorporates
391  * terms for time and spatial variations of radial velocity (@f$ V @f$). The
392  * default boundary condition is zero radial velocity (@f$ V @f$) at the left
393  * and right boundary.
394  *
395  * For argument explanation, see evalContinuity().
396  */
397  virtual void evalMomentum(double* x, double* rsd, int* diag,
398  double rdt, size_t jmin, size_t jmax);
399 
400  /**
401  * Evaluate the lambda equation residual.
402  *
403  * The function calculates the lambda equation as
404  * @f[
405  * \frac{d\Lambda}{dz} = 0
406  * @f]
407  *
408  * The lambda equation serves as an eigenvalue that allows the momentum
409  * equation and continuity equations to be simultaneously satisfied in
410  * axisymmetric flows. The lambda equation propagates information from
411  * left-to-right. The default boundary condition is @f$ \Lambda = 0 @f$
412  * at the left and zero flux at the right boundary.
413  *
414  * For argument explanation, see evalContinuity().
415  */
416  virtual void evalLambda(double* x, double* rsd, int* diag,
417  double rdt, size_t jmin, size_t jmax);
418 
419  /**
420  * Evaluate the energy equation residual.
421  *
422  * The function calculates the energy equation:
423  * @f[
424  * \rho c_p u \frac{dT}{dz} =
425  * \frac{d}{dz}\left( \lambda \frac{dT}{dz} \right)
426  * - \sum_k h_kW_k\dot{\omega}_k
427  * - \sum_k j_k \frac{dh_k}{dz}
428  * @f]
429  *
430  * The energy equation includes contributions from
431  * chemical reactions and diffusion. Default is zero temperature (@f$ T @f$)
432  * at the left and right boundaries. These boundary values are updated by the
433  * specific boundary object connected to the domain.
434  *
435  * For argument explanation, see evalContinuity().
436  */
437  virtual void evalEnergy(double* x, double* rsd, int* diag,
438  double rdt, size_t jmin, size_t jmax);
439 
440  /**
441  * Evaluate the species equations' residuals.
442  *
443  * The function calculates the species equations as
444  * @f[
445  * \rho u \frac{dY_k}{dz} + \frac{dj_k}{dz} = W_k\dot{\omega}_k
446  * @f]
447  *
448  * The species equations include terms for temporal and spatial variations
449  * of species mass fractions (@f$ Y_k @f$). The default boundary condition is zero
450  * flux for species at the left and right boundary.
451  *
452  * For argument explanation, see evalContinuity().
453  */
454  virtual void evalSpecies(double* x, double* rsd, int* diag,
455  double rdt, size_t jmin, size_t jmax);
456 
457  /**
458  * Evaluate the electric field equation residual to be zero everywhere.
459  *
460  * The electric field equation is implemented in the IonFlow class. The default
461  * boundary condition is zero electric field (@f$ E @f$) at the boundary,
462  * and @f$ E @f$ is zero within the domain.
463  *
464  * For argument explanation, see evalContinuity().
465  */
466  virtual void evalElectricField(double* x, double* rsd, int* diag,
467  double rdt, size_t jmin, size_t jmax);
468 
469  /**
470  * Update the thermodynamic properties from point j0 to point j1
471  * (inclusive), based on solution x.
472  *
473  * The gas state is set to be consistent with the solution at the
474  * points from j0 to j1.
475  *
476  * Properties that are computed and cached are:
477  * * #m_rho (density)
478  * * #m_wtm (mean molecular weight)
479  * * #m_cp (specific heat capacity)
480  * * #m_hk (species specific enthalpies)
481  * * #m_wdot (species production rates)
482  */
483  void updateThermo(const double* x, size_t j0, size_t j1) {
484  for (size_t j = j0; j <= j1; j++) {
485  setGas(x,j);
486  m_rho[j] = m_thermo->density();
487  m_wtm[j] = m_thermo->meanMolecularWeight();
488  m_cp[j] = m_thermo->cp_mass();
489  m_thermo->getPartialMolarEnthalpies(&m_hk(0, j));
490  m_kin->getNetProductionRates(&m_wdot(0, j));
491  }
492  }
493 
494  //! @name Solution components
495  //! @{
496 
497  double T(const double* x, size_t j) const {
498  return x[index(c_offset_T, j)];
499  }
500  double& T(double* x, size_t j) {
501  return x[index(c_offset_T, j)];
502  }
503  double T_prev(size_t j) const {
504  return prevSoln(c_offset_T, j);
505  }
506 
507  double rho_u(const double* x, size_t j) const {
508  return m_rho[j]*x[index(c_offset_U, j)];
509  }
510 
511  double u(const double* x, size_t j) const {
512  return x[index(c_offset_U, j)];
513  }
514 
515  double V(const double* x, size_t j) const {
516  return x[index(c_offset_V, j)];
517  }
518  double V_prev(size_t j) const {
519  return prevSoln(c_offset_V, j);
520  }
521 
522  double lambda(const double* x, size_t j) const {
523  return x[index(c_offset_L, j)];
524  }
525 
526  double Y(const double* x, size_t k, size_t j) const {
527  return x[index(c_offset_Y + k, j)];
528  }
529 
530  double& Y(double* x, size_t k, size_t j) {
531  return x[index(c_offset_Y + k, j)];
532  }
533 
534  double Y_prev(size_t k, size_t j) const {
535  return prevSoln(c_offset_Y + k, j);
536  }
537 
538  double X(const double* x, size_t k, size_t j) const {
539  return m_wtm[j]*Y(x,k,j)/m_wt[k];
540  }
541 
542  double flux(size_t k, size_t j) const {
543  return m_flux(k, j);
544  }
545  //! @}
546 
547  //! @name convective spatial derivatives.
548  //!
549  //! These use upwind differencing, assuming u(z) is negative
550  //! @{
551  double dVdz(const double* x, size_t j) const {
552  size_t jloc = (u(x,j) > 0.0 ? j : j + 1);
553  return (V(x,jloc) - V(x,jloc-1))/m_dz[jloc-1];
554  }
555 
556  double dYdz(const double* x, size_t k, size_t j) const {
557  size_t jloc = (u(x,j) > 0.0 ? j : j + 1);
558  return (Y(x,k,jloc) - Y(x,k,jloc-1))/m_dz[jloc-1];
559  }
560 
561  double dTdz(const double* x, size_t j) const {
562  size_t jloc = (u(x,j) > 0.0 ? j : j + 1);
563  return (T(x,jloc) - T(x,jloc-1))/m_dz[jloc-1];
564  }
565  //! @}
566 
567  double shear(const double* x, size_t j) const {
568  double c1 = m_visc[j-1]*(V(x,j) - V(x,j-1));
569  double c2 = m_visc[j]*(V(x,j+1) - V(x,j));
570  return 2.0*(c2/(z(j+1) - z(j)) - c1/(z(j) - z(j-1)))/(z(j+1) - z(j-1));
571  }
572 
573  double divHeatFlux(const double* x, size_t j) const {
574  double c1 = m_tcon[j-1]*(T(x,j) - T(x,j-1));
575  double c2 = m_tcon[j]*(T(x,j+1) - T(x,j));
576  return -2.0*(c2/(z(j+1) - z(j)) - c1/(z(j) - z(j-1)))/(z(j+1) - z(j-1));
577  }
578 
579  size_t mindex(size_t k, size_t j, size_t m) {
580  return m*m_nsp*m_nsp + m_nsp*j + k;
581  }
582 
583  //! Update the diffusive mass fluxes.
584  virtual void updateDiffFluxes(const double* x, size_t j0, size_t j1);
585 
586  //! Get the gradient of species specific molar enthalpies
587  virtual void grad_hk(const double* x, size_t j);
588 
589  //---------------------------------------------------------
590  // member data
591  //---------------------------------------------------------
592 
593  double m_press = -1.0; // pressure
594 
595  // grid parameters
596  vector<double> m_dz;
597 
598  // mixture thermo properties
599  vector<double> m_rho;
600  vector<double> m_wtm;
601 
602  // species thermo properties
603  vector<double> m_wt;
604  vector<double> m_cp;
605 
606  // transport properties
607  vector<double> m_visc;
608  vector<double> m_tcon;
609  //! Array of size #m_nsp by #m_points for saving density times diffusion
610  //! coefficient times species molar mass divided by mean molecular weight
611  vector<double> m_diff;
612  vector<double> m_multidiff;
613  Array2D m_dthermal;
614  Array2D m_flux;
615 
616  //! Array of size #m_nsp by #m_points for saving molar enthalpies
618 
619  //! Array of size #m_nsp by #m_points-1 for saving enthalpy fluxes
621 
622  // production rates
623  Array2D m_wdot;
624 
625  size_t m_nsp; //!< Number of species in the mechanism
626 
627  ThermoPhase* m_thermo = nullptr;
628  Kinetics* m_kin = nullptr;
629  Transport* m_trans = nullptr;
630 
631  // boundary emissivities for the radiation calculations
632  double m_epsilon_left = 0.0;
633  double m_epsilon_right = 0.0;
634 
635  //! Indices within the ThermoPhase of the radiating species. First index is
636  //! for CO2, second is for H2O.
637  vector<size_t> m_kRadiating;
638 
639  // flags
640  vector<bool> m_do_energy;
641  bool m_do_soret = false;
642  vector<bool> m_do_species;
643  bool m_do_multicomponent = false;
644 
645  //! flag for the radiative heat loss
646  bool m_do_radiation = false;
647 
648  //! radiative heat loss vector
649  vector<double> m_qdotRadiation;
650 
651  // fixed T and Y values
652  vector<double> m_fixedtemp;
653  vector<double> m_zfix;
654  vector<double> m_tfix;
655 
656  //! Index of species with a large mass fraction at each boundary, for which
657  //! the mass fraction may be calculated as 1 minus the sum of the other mass
658  //! fractions
659  size_t m_kExcessLeft = 0;
660  size_t m_kExcessRight = 0;
661 
662  bool m_dovisc;
663  bool m_isFree;
664  bool m_usesLambda;
665 
666  //! Update the transport properties at grid points in the range from `j0`
667  //! to `j1`, based on solution `x`.
668  virtual void updateTransport(double* x, size_t j0, size_t j1);
669 
670 public:
671  //! Location of the point where temperature is fixed
672  double m_zfixed = Undef;
673 
674  //! Temperature at the point used to fix the flame location
675  double m_tfixed = -1.0;
676 
677 private:
678  vector<double> m_ybar;
679 };
680 
681 }
682 
683 #endif
Header file for class Cantera::Array2D.
Base class for kinetics managers and also contains the kineticsmgr module documentation (see Kinetics...
Header file for class ThermoPhase, the base class for phases with thermodynamic properties,...
A map of string keys to values whose type can vary at runtime.
Definition: AnyMap.h:427
A class for 2D arrays stored in column-major (Fortran-compatible) form.
Definition: Array.h:32
Base class for one-dimensional domains.
Definition: Domain1D.h:28
double prevSoln(size_t n, size_t j) const
Value of component n at point j in the previous solution.
Definition: Domain1D.h:441
Public interface for kinetics managers.
Definition: Kinetics.h:125
virtual void getNetProductionRates(double *wdot)
Species net production rates [kmol/m^3/s or kmol/m^2/s].
Definition: Kinetics.cpp:363
double meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
Definition: Phase.h:655
virtual double density() const
Density (kg/m^3).
Definition: Phase.h:587
A container class holding arrays of state information.
Definition: SolutionArray.h:33
This class represents 1D flow domains that satisfy the one-dimensional similarity solution for chemic...
Definition: StFlow.h:45
void eval(size_t jGlobal, double *xGlobal, double *rsdGlobal, integer *diagGlobal, double rdt) override
Evaluate the residual functions for axisymmetric stagnation flow.
Definition: StFlow.cpp:296
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...
Definition: StFlow.h:137
size_t m_kExcessLeft
Index of species with a large mass fraction at each boundary, for which the mass fraction may be calc...
Definition: StFlow.h:659
void setMeta(const AnyMap &state) override
Retrieve meta data.
Definition: StFlow.cpp:872
void setTransportModel(const string &trans)
Set the transport model.
Definition: StFlow.cpp:211
double leftEmissivity() const
Return emissivity at left boundary.
Definition: StFlow.h:232
void setTransport(shared_ptr< Transport > trans) override
Set transport model to existing instance.
Definition: StFlow.cpp:140
void setUnstrainedFlow()
Set flow configuration for burner-stabilized flames, using specified inlet mass fluxes.
Definition: StFlow.h:180
void setKinetics(shared_ptr< Kinetics > kin) override
Set the kinetics manager.
Definition: StFlow.cpp:134
void resetBadValues(double *xg) override
When called, this function should reset "bad" values in the state vector such as negative species con...
Definition: StFlow.cpp:201
size_t rightExcessSpecies() const
Index of the species on the right boundary with the largest mass fraction.
Definition: StFlow.h:311
vector< double > m_qdotRadiation
radiative heat loss vector
Definition: StFlow.h:649
virtual void evalMomentum(double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax)
Evaluate the momentum equation residual.
Definition: StFlow.cpp:462
double pressure() const
The current pressure [Pa].
Definition: StFlow.h:116
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.
Definition: StFlow.h:483
void resize(size_t components, size_t points) override
Change the grid size. Called after grid refinement.
Definition: StFlow.cpp:160
StFlow(ThermoPhase *ph=0, size_t nsp=1, size_t points=1)
Create a new flow domain.
Definition: StFlow.cpp:19
void enableSoret(bool withSoret)
Enable thermal diffusion, also known as Soret diffusion.
Definition: StFlow.h:102
void setFixedTempProfile(vector< double > &zfixed, vector< double > &tfixed)
Sometimes it is desired to carry out the simulation using a specified temperature profile,...
Definition: StFlow.h:128
virtual void evalContinuity(double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax)
Evaluate the continuity equation residual.
Definition: StFlow.cpp:405
void setBoundaryEmissivities(double e_left, double e_right)
Set the emissivities for the boundary values.
Definition: StFlow.cpp:982
virtual void evalEnergy(double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax)
Evaluate the energy equation residual.
Definition: StFlow.cpp:522
void enableRadiation(bool doRadiation)
Turn radiation on / off.
Definition: StFlow.h:209
vector< double > m_diff
Array of size m_nsp by m_points for saving density times diffusion coefficient times species molar ma...
Definition: StFlow.h:611
shared_ptr< SolutionArray > asArray(const double *soln) const override
Save the state of this domain as a SolutionArray.
Definition: StFlow.cpp:808
size_t componentIndex(const string &name) const override
index of component with name name.
Definition: StFlow.cpp:722
double rightEmissivity() const
Return emissivity at right boundary.
Definition: StFlow.h:235
void setGas(const double *x, size_t j)
Set the gas object state to be consistent with the solution at point j.
Definition: StFlow.cpp:229
double m_tfixed
Temperature at the point used to fix the flame location.
Definition: StFlow.h:675
bool radiationEnabled() const
Returns true if the radiation term in the energy equation is enabled.
Definition: StFlow.h:214
virtual bool componentActive(size_t n) const
Returns true if the specified component is an active part of the solver state.
Definition: StFlow.cpp:745
Array2D m_hk
Array of size m_nsp by m_points for saving molar enthalpies.
Definition: StFlow.h:617
void setFreeFlow()
Set flow configuration for freely-propagating flames, using an internal point with a fixed temperatur...
Definition: StFlow.h:164
virtual bool doElectricField(size_t j) const
Retrieve flag indicating whether electric field is solved or not (used by IonFlow specialization)
Definition: StFlow.cpp:976
virtual void evalSpecies(double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax)
Evaluate the species equations' residuals.
Definition: StFlow.cpp:560
void setupGrid(size_t n, const double *z) override
called to set up initial grid, and after grid refinement
Definition: StFlow.cpp:186
size_t leftExcessSpecies() const
Index of the species on the left boundary with the largest mass fraction.
Definition: StFlow.h:306
Array2D m_dhk_dz
Array of size m_nsp by m_points-1 for saving enthalpy fluxes.
Definition: StFlow.h:620
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.
Definition: StFlow.cpp:599
double radiativeHeatLoss(size_t j) const
Return radiative heat loss at grid point j.
Definition: StFlow.h:219
double m_zfixed
Location of the point where temperature is fixed.
Definition: StFlow.h:672
void _finalize(const double *x) override
In some cases, a domain may need to set parameters that depend on the initial solution estimate.
Definition: StFlow.cpp:249
virtual size_t getSolvingStage() const
Get the solving stage (used by IonFlow specialization)
Definition: StFlow.cpp:952
size_t m_nsp
Number of species in the mechanism.
Definition: StFlow.h:625
virtual void evalLambda(double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax)
Evaluate the lambda equation residual.
Definition: StFlow.cpp:494
void fromArray(SolutionArray &arr, double *soln) override
Restore the solution for this domain from a SolutionArray.
Definition: StFlow.cpp:842
AnyMap getMeta() const override
Retrieve meta data.
Definition: StFlow.cpp:759
virtual void updateDiffFluxes(const double *x, size_t j0, size_t j1)
Update the diffusive mass fluxes.
Definition: StFlow.cpp:660
string componentName(size_t n) const override
Name of the nth component. May be overloaded.
Definition: StFlow.cpp:700
bool isFree() const
Retrieve flag indicating whether flow is freely propagating.
Definition: StFlow.h:263
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: StFlow.cpp:237
virtual void grad_hk(const double *x, size_t j)
Get the gradient of species specific molar enthalpies.
Definition: StFlow.cpp:1020
bool isStrained() const
Retrieve flag indicating whether flow uses radial momentum.
Definition: StFlow.h:274
string transportModel() const
Retrieve transport model.
Definition: StFlow.cpp:216
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 qdotRadia...
Definition: StFlow.cpp:358
virtual void updateProperties(size_t jg, double *x, size_t jmin, size_t jmax)
Update the properties (thermo, transport, and diffusion flux).
Definition: StFlow.cpp:334
string domainType() const override
Domain type flag.
Definition: StFlow.cpp:124
void show(const double *x) override
Print the solution.
Definition: StFlow.cpp:643
virtual void setSolvingStage(const size_t stage)
Solving stage mode for handling ionized species (used by IonFlow specialization)
Definition: StFlow.cpp:958
void setPressure(double p)
Set the pressure.
Definition: StFlow.h:111
virtual void fixElectricField(size_t j=npos)
Set to fix voltage in a point (used by IonFlow specialization)
Definition: StFlow.cpp:970
void setAxisymmetricFlow()
Set flow configuration for axisymmetric counterflow flames, using specified inlet mass fluxes.
Definition: StFlow.h:172
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.
Definition: StFlow.cpp:608
virtual void solveElectricField(size_t j=npos)
Set to solve electric field in a point (used by IonFlow specialization)
Definition: StFlow.cpp:964
void _getInitialSoln(double *x) override
Write the initial solution estimate into array x.
Definition: StFlow.cpp:220
vector< size_t > m_kRadiating
Indices within the ThermoPhase of the radiating species.
Definition: StFlow.h:637
double T_fixed(size_t j) const
The fixed temperature value at point j.
Definition: StFlow.h:143
bool m_do_radiation
flag for the radiative heat loss
Definition: StFlow.h:646
Base class for a phase with thermodynamic properties.
Definition: ThermoPhase.h:390
virtual void getPartialMolarEnthalpies(double *hbar) const
Returns an array of partial molar enthalpies for the species in the mixture.
Definition: ThermoPhase.h:801
double cp_mass() const
Specific heat at constant pressure. Units: J/kg/K.
Definition: ThermoPhase.h:1048
Base class for transport property managers.
Definition: Transport.h:72
Namespace for the Cantera kernel.
Definition: AnyMap.cpp:564
const size_t npos
index returned by functions to indicate "no position"
Definition: ct_defs.h:180
const double Undef
Fairly random number to be used to initialize variables against to see if they are subsequently defin...
Definition: ct_defs.h:164
offset
Offsets of solution components in the 1D solution array.
Definition: StFlow.h:24
@ c_offset_U
axial velocity
Definition: StFlow.h:25
@ c_offset_L
(1/r)dP/dr
Definition: StFlow.h:28
@ c_offset_V
strain rate
Definition: StFlow.h:26
@ c_offset_E
electric field equation
Definition: StFlow.h:29
@ c_offset_Y
mass fractions
Definition: StFlow.h:30
@ c_offset_T
temperature
Definition: StFlow.h:27