8#include "cantera/oneD/refine.h"
18IonFlow::IonFlow(ThermoPhase* ph,
size_t nsp,
size_t points) :
19 Flow1D(ph, nsp, points)
22 for (
size_t k = 0; k < m_nsp; k++) {
23 m_speciesCharge.push_back(m_thermo->charge(k));
27 for (
size_t k = 0; k < m_nsp; k++){
28 if (m_speciesCharge[k] != 0){
29 m_kCharge.push_back(k);
31 m_kNeutral.push_back(k);
36 if (m_thermo->speciesIndex(
"E") != npos ) {
37 m_kElectron = m_thermo->speciesIndex(
"E");
41 setBounds(c_offset_E, -1.0e20, 1.0e20);
46 for (
size_t k : m_kCharge) {
47 setBounds(c_offset_Y + k, -1e-14, 1.0);
49 setBounds(c_offset_Y + m_kElectron, -1e-18, 1.0);
51 m_refiner->setActive(c_offset_E,
false);
52 m_mobility.resize(m_nsp*m_points);
53 m_do_electric_field.resize(m_points,
false);
56IonFlow::IonFlow(shared_ptr<Solution> sol,
const string&
id,
size_t points)
57 :
IonFlow(sol->thermo().get(), sol->thermo()->nSpecies(), points)
65 "An appropriate transport model\nshould be set when instantiating the "
66 "Solution ('gas') object.");
68 m_solution->registerChangedCallback(
this, [
this]() {
76 return "free-ion-flow";
79 return "axisymmetric-ion-flow";
81 return "unstrained-ion-flow";
87 m_do_species.resize(
m_nsp,
true);
103 for (
size_t j = j0; j < j1; j++) {
110 double rho = m_thermo->
density();
129 for (
size_t j = j0; j < j1; j++) {
130 double dz = z(j+1) - z(j);
134 m_flux(k,j) *= (X(x,k,j) - X(x,k,j+1))/dz;
140 m_flux(k,j) += sum*Y(x,k,j);
154 for (
size_t j = j0; j < j1; j++) {
155 double rho = density(j);
156 double dz = z(j+1) - z(j);
159 for (
size_t k = 0; k <
m_nsp; k++) {
161 m_flux(k,j) *= (X(x,k,j) - X(x,k,j+1))/dz;
165 double E_ambi =
E(x,j);
167 double Yav = 0.5 * (Y(x,k,j) + Y(x,k,j+1));
168 double drift = rho * Yav * E_ambi
170 m_flux(k,j) += drift;
174 double sum_flux = 0.0;
175 for (
size_t k = 0; k <
m_nsp; k++) {
176 sum_flux -= m_flux(k,j);
178 double sum_ion = 0.0;
184 m_flux(k,j) += Y(x,k,j) / (1-sum_ion) * sum_flux;
191 if (stage == 1 || stage == 2) {
195 "solution stage must be set to: "
196 "1) frozenIonMethod, "
197 "2) electricFieldEqnMethod");
203 double rdt,
size_t jmin,
size_t jmax)
219 size_t j0 = std::max<size_t>(jmin, 1);
220 size_t j1 = std::min(jmax,
m_points - 2);
221 for (
size_t j = j0; j <= j1; j++) {
228 double rdt,
size_t jmin,
size_t jmax)
240 rsd[index(
c_offset_Y + k, jmin)] = Y(x,k,jmin) - Y(x,k,jmin + 1);
247 bool changed =
false;
249 for (
size_t i = 0; i <
m_points; i++) {
272 bool changed =
false;
274 for (
size_t i = 0; i <
m_points; i++) {
296 vector<double>& mobi_e)
300 size_t n = tfix.size();
302 for (
size_t i = 0; i < n; i++) {
303 tlog.push_back(log(tfix[i]));
305 vector<double> w(n, -1.0);
306 m_diff_e_fix.resize(degree + 1);
308 polyfit(n, degree, tlog.data(), diff_e.data(), w.data(), m_diff_e_fix.data());
Headers for the Transport object, which is the virtual base class for all transport property evaluato...
Base class for exceptions thrown by Cantera classes.
shared_ptr< Solution > m_solution
Composite thermo/kinetics/transport handler.
void setSolution(shared_ptr< Solution > sol)
Set the solution manager.
size_t m_points
Number of grid points.
string m_id
Identity tag for the domain.
void needJacUpdate()
Set this if something has changed in the governing equations (for example, the value of a constant ha...
void setTransport(shared_ptr< Transport > trans) override
Set transport model to existing instance.
void setKinetics(shared_ptr< Kinetics > kin) override
Set the kinetics manager.
void resize(size_t components, size_t points) override
Change the grid size. Called after grid refinement.
vector< double > m_diff
Array of size m_nsp by m_points for saving density times diffusion coefficient times species molar ma...
virtual bool componentActive(size_t n) const
Returns true if the specified component is an active part of the solver state.
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.
void _finalize(const double *x) override
In some cases, a domain may need to set parameters that depend on the initial solution estimate.
size_t m_nsp
Number of species in the mechanism.
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.
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.
This class models the ion transportation in a flame.
vector< size_t > m_kCharge
index of species with charges
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.
double E(const double *x, size_t j) const
electric field
vector< bool > m_do_electric_field
flag for solving electric field or not
size_t m_kElectron
index of electron
void frozenIonMethod(const double *x, size_t j0, size_t j1)
Solving phase one: the fluxes of charged species are turned off.
void resize(size_t components, size_t points) override
Change the grid size. Called after grid refinement.
void setElectronTransport(vector< double > &tfix, vector< double > &diff_e, vector< double > &mobi_e)
Sometimes it is desired to carry out the simulation using a specified electron transport profile,...
void evalElectricField(double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) override
Evaluate the electric field equation residual by Gauss's law.
double rho_e(double *x, size_t j) const
total charge density
void updateTransport(double *x, size_t j0, size_t j1) override
Update the transport properties at grid points in the range from j0 to j1, based on solution x.
vector< double > m_mobility
mobility
void updateDiffFluxes(const double *x, size_t j0, size_t j1) override
Update the diffusive mass fluxes.
void evalSpecies(double *x, double *rsd, int *diag, double rdt, size_t jmin, size_t jmax) override
Evaluate the species equations' residual.
size_t m_stage
solving stage
void _finalize(const double *x) override
In some cases, a domain may need to set parameters that depend on the initial solution estimate.
void setSolvingStage(const size_t stage) override
Solving stage mode for handling ionized species (used by IonFlow specialization)
bool m_import_electron_transport
flag for importing transport of electron
void solveElectricField(size_t j=npos) override
Set to solve electric field in a point (used by IonFlow specialization)
void fixElectricField(size_t j=npos) override
Set to fix voltage in a point (used by IonFlow specialization)
string domainType() const override
Domain type flag.
vector< size_t > m_kNeutral
index of neutral species
vector< double > m_mobi_e_fix
coefficients of polynomial fitting of fixed electron transport profile
bool componentActive(size_t n) const override
Returns true if the specified component is an active part of the solver state.
vector< double > m_speciesCharge
electrical properties
double temperature() const
Temperature (K).
double meanMolecularWeight() const
The mean molecular weight. Units: (kg/kmol)
virtual double density() const
Density (kg/m^3).
virtual string transportModel() const
Identifies the model represented by this Transport object.
virtual void getMobilities(double *const mobil_e)
Get the Electrical mobilities (m^2/V/s).
Header for a file containing miscellaneous numerical functions.
This file contains definitions for utility functions and text for modules, inputfiles and logging,...
R poly5(D x, R *c)
Templated evaluation of a polynomial of order 5.
double polyfit(size_t n, size_t deg, const double *xp, const double *yp, const double *wp, double *pp)
Fits a polynomial function to a set of data points.
const double epsilon_0
Permittivity of free space [F/m].
Namespace for the Cantera kernel.
const size_t npos
index returned by functions to indicate "no position"
@ c_offset_U
axial velocity [m/s]
@ c_offset_E
electric field
@ c_offset_Y
mass fractions
@ c_offset_T
temperature [kelvin]
Various templated functions that carry out common vector and polynomial operations (see Templated Arr...