GRINS-0.8.0
List of all members | Public Member Functions | Protected Attributes
GRINS::GasRecombinationCatalyticWall< Chemistry > Class Template Reference

#include <gas_recombination_catalytic_wall.h>

Inheritance diagram for GRINS::GasRecombinationCatalyticWall< Chemistry >:
Inheritance graph
[legend]
Collaboration diagram for GRINS::GasRecombinationCatalyticWall< Chemistry >:
Collaboration graph
[legend]

Public Member Functions

 GasRecombinationCatalyticWall (SharedPtr< Chemistry > &chem, SharedPtr< CatalycityBase > &gamma, const std::vector< VariableIndex > &species_vars, VariableIndex T_var, libMesh::Real p0, unsigned int reactant_species_idx, unsigned int product_species_idx)
 
 GasRecombinationCatalyticWall (const Chemistry &chem_mixture, CatalycityBase &gamma, const unsigned int reactant_species_idx, const unsigned int product_species_idx)
 Deprecated. More...
 
virtual ~GasRecombinationCatalyticWall ()
 
virtual void init (const libMesh::FEMSystem &system)
 Deprecated. More...
 
virtual void apply_fluxes (AssemblyContext &context, const CachedValues &cache, const bool request_jacobian)
 Deprecated. More...
 
virtual bool eval_flux (bool compute_jacobian, AssemblyContext &context, libMesh::Real sign, bool is_axisymmetric)
 
libMesh::Real compute_reactant_mass_flux (const libMesh::Real rho, const libMesh::Real Y_r, const libMesh::Real T)
 
libMesh::Real compute_product_mass_flux (const libMesh::Real rho, const libMesh::Real Y_r, const libMesh::Real T)
 
- Public Member Functions inherited from GRINS::CatalyticWallBase< Chemistry >
 CatalyticWallBase (SharedPtr< Chemistry > &chem, SharedPtr< CatalycityBase > &gamma, const std::vector< VariableIndex > &species_vars, VariableIndex T_var, libMesh::Real p0, unsigned int reactant_species_idx)
 
 CatalyticWallBase (const Chemistry &chem_mixture, CatalycityBase &gamma, const unsigned int reactant_species_idx)
 Deprecated. More...
 
virtual ~CatalyticWallBase ()
 
libMesh::Real rho (libMesh::Real T, libMesh::Real p0, libMesh::Real R_mix) const
 
libMesh::Real omega_dot (const libMesh::Real rho_s, const libMesh::Real T) const
 $ \rho_s \gamma \sqrt{ \frac{R_s T}{2\pi} } $ More...
 
libMesh::Real domega_dot_dws (const libMesh::Real rho_s, const libMesh::Real w_s, const libMesh::Real T, const libMesh::Real R) const
 
libMesh::Real domega_dot_dT (const libMesh::Real rho_s, const libMesh::Real T) const
 
void set_catalycity_params (const std::vector< libMesh::Real > &params)
 
virtual void register_parameter (const std::string &param_name, libMesh::ParameterMultiAccessor< libMesh::Number > &param_pointer) const
 Each subclass will register its copy of an independent. More...
 
- Public Member Functions inherited from GRINS::NeumannBCAbstract
 NeumannBCAbstract ()
 
virtual ~NeumannBCAbstract ()
 
- Public Member Functions inherited from GRINS::ParameterUser
 ParameterUser (const std::string &user_name)
 
virtual ~ParameterUser ()
 
virtual void set_parameter (libMesh::Number &param_variable, const GetPot &input, const std::string &param_name, libMesh::Number param_default)
 Each subclass can simultaneously read a parameter value from. More...
 
virtual void set_parameter (libMesh::ParsedFunction< libMesh::Number, libMesh::Gradient > &func, const GetPot &input, const std::string &func_param_name, const std::string &param_default)
 Each subclass can simultaneously read a parsed function from. More...
 
virtual void set_parameter (libMesh::ParsedFEMFunction< libMesh::Number > &func, const GetPot &input, const std::string &func_param_name, const std::string &param_default)
 Each subclass can simultaneously read a parsed function from. More...
 
virtual void move_parameter (const libMesh::Number &old_parameter, libMesh::Number &new_parameter)
 When cloning an object, we need to update parameter pointers. More...
 
virtual void move_parameter (const libMesh::ParsedFunction< libMesh::Number, libMesh::Gradient > &old_func, libMesh::ParsedFunction< libMesh::Number, libMesh::Gradient > &new_func)
 When cloning an object, we need to update parameter pointers. More...
 
virtual void move_parameter (const libMesh::ParsedFEMFunction< libMesh::Number > &old_func, libMesh::ParsedFEMFunction< libMesh::Number > &new_func)
 When cloning an object, we need to update parameter pointers. More...
 

Protected Attributes

const unsigned int _reactant_species_idx
 
VariableIndex _reactant_var_idx
 
const unsigned int _product_species_idx
 
VariableIndex _product_var_idx
 
- Protected Attributes inherited from GRINS::CatalyticWallBase< Chemistry >
SharedPtr< Chemistry > _chem_ptr
 
const Chemistry & _chemistry
 Deprecated. More...
 
SharedPtr< CatalycityBase_gamma_ptr
 
libMesh::UniquePtr< CatalycityBase_gamma_s
 Deprecated. More...
 
const libMesh::Real _C
 $ \sqrt{ \frac{R_s}{2\pi} } $ More...
 
std::vector< VariableIndex_species_vars
 
VariableIndex _T_var
 
libMesh::Real _p0
 Thermodynamic pressure. More...
 

Additional Inherited Members

- Static Public Attributes inherited from GRINS::ParameterUser
static std::string zero_vector_function = std::string("{0}")
 A parseable function string with LIBMESH_DIM components, all 0. More...
 
- Protected Member Functions inherited from GRINS::CatalyticWallBase< Chemistry >
libMesh::Real eval_gamma (libMesh::Real T) const
 Temporary helper to deal with intermediate refactoring. More...
 
libMesh::Real eval_gamma_dT (libMesh::Real T) const
 Temporary helper to deal with intermediate refactoring. More...
 

Detailed Description

template<typename Chemistry>
class GRINS::GasRecombinationCatalyticWall< Chemistry >

Definition at line 36 of file gas_recombination_catalytic_wall.h.

Constructor & Destructor Documentation

template<typename Chemistry >
GRINS::GasRecombinationCatalyticWall< Chemistry >::GasRecombinationCatalyticWall ( SharedPtr< Chemistry > &  chem,
SharedPtr< CatalycityBase > &  gamma,
const std::vector< VariableIndex > &  species_vars,
VariableIndex  T_var,
libMesh::Real  p0,
unsigned int  reactant_species_idx,
unsigned int  product_species_idx 
)

Definition at line 37 of file gas_recombination_catalytic_wall.C.

44  : CatalyticWallBase<Chemistry>(chem,gamma,species_vars,T_var,p0,reactant_species_idx),
45  _reactant_species_idx(reactant_species_idx),
46  _reactant_var_idx(species_vars[reactant_species_idx]),
47  _product_species_idx(product_species_idx),
48  _product_var_idx(species_vars[product_species_idx])
49  {}
template<typename Chemistry >
GRINS::GasRecombinationCatalyticWall< Chemistry >::GasRecombinationCatalyticWall ( const Chemistry &  chem_mixture,
CatalycityBase gamma,
const unsigned int  reactant_species_idx,
const unsigned int  product_species_idx 
)

Deprecated.

Definition at line 52 of file gas_recombination_catalytic_wall.C.

56  : CatalyticWallBase<Chemistry>(chem_mixture,gamma,reactant_species_idx),
57  _reactant_species_idx(reactant_species_idx),
58  _product_species_idx(product_species_idx)
59  {
60  libmesh_deprecated();
61  }
template<typename Chemistry>
virtual GRINS::GasRecombinationCatalyticWall< Chemistry >::~GasRecombinationCatalyticWall ( )
inlinevirtual

Definition at line 54 of file gas_recombination_catalytic_wall.h.

54 {};

Member Function Documentation

template<typename Chemistry >
void GRINS::GasRecombinationCatalyticWall< Chemistry >::apply_fluxes ( AssemblyContext context,
const CachedValues cache,
const bool  request_jacobian 
)
virtual

Deprecated.

Implements GRINS::CatalyticWallBase< Chemistry >.

Definition at line 150 of file gas_recombination_catalytic_wall.C.

References GRINS::CachedValues::get_cached_values(), GRINS::CachedValues::get_cached_vector_values(), GRINS::Physics::is_axisymmetric(), GRINS::Cache::MASS_FRACTIONS, GRINS::Cache::MIXTURE_DENSITY, and GRINS::Cache::TEMPERATURE.

153  {
154  libmesh_do_once(libmesh_deprecated());
155 
156  libMesh::FEGenericBase<libMesh::Real>* side_fe = NULL;
157  context.get_side_fe( _reactant_var_idx, side_fe );
158 
159  // The number of local degrees of freedom in each variable.
160  const unsigned int n_var_dofs = context.get_dof_indices(_reactant_var_idx).size();
161 
162  libmesh_assert_equal_to( n_var_dofs, context.get_dof_indices(_product_var_idx).size() );
163 
164  // Element Jacobian * quadrature weight for side integration.
165  const std::vector<libMesh::Real> &JxW_side = side_fe->get_JxW();
166 
167  // The var shape functions at side quadrature points.
168  const std::vector<std::vector<libMesh::Real> >& var_phi_side = side_fe->get_phi();
169 
170  // Physical location of the quadrature points
171  const std::vector<libMesh::Point>& var_qpoint = side_fe->get_xyz();
172 
173  // reactant residual
174  libMesh::DenseSubVector<libMesh::Number> &F_r_var = context.get_elem_residual(_reactant_var_idx);
175 
176  // product residual
177  libMesh::DenseSubVector<libMesh::Number> &F_p_var = context.get_elem_residual(_product_var_idx);
178 
179  unsigned int n_qpoints = context.get_side_qrule().n_points();
180 
181  for (unsigned int qp=0; qp != n_qpoints; qp++)
182  {
183  libMesh::Real jac = JxW_side[qp];
184 
186  {
187  const libMesh::Number r = var_qpoint[qp](0);
188  jac *= r;
189  }
190 
191  const libMesh::Real rho = cache.get_cached_values(Cache::MIXTURE_DENSITY)[qp];
192 
193  const libMesh::Real Y_r = cache.get_cached_vector_values(Cache::MASS_FRACTIONS)[qp][this->_reactant_species_idx];
194 
195  const libMesh::Real T = cache.get_cached_values(Cache::TEMPERATURE)[qp];
196 
197  const libMesh::Real r_value = this->compute_reactant_mass_flux(rho, Y_r, T);
198 
199  const libMesh::Real p_value = -r_value;
200 
201  for (unsigned int i=0; i != n_var_dofs; i++)
202  {
203  F_r_var(i) += r_value*var_phi_side[i][qp]*jac;
204 
205  F_p_var(i) += p_value*var_phi_side[i][qp]*jac;
206 
207  if( request_jacobian )
208  {
209  libmesh_not_implemented();
210  }
211  }
212  }
213  }
static bool is_axisymmetric()
Definition: physics.h:132
libMesh::Real compute_reactant_mass_flux(const libMesh::Real rho, const libMesh::Real Y_r, const libMesh::Real T)
libMesh::Real rho(libMesh::Real T, libMesh::Real p0, libMesh::Real R_mix) const
template<typename Chemistry >
libMesh::Real GRINS::GasRecombinationCatalyticWall< Chemistry >::compute_product_mass_flux ( const libMesh::Real  rho,
const libMesh::Real  Y_r,
const libMesh::Real  T 
)
inline

Definition at line 102 of file gas_recombination_catalytic_wall.h.

105  {
106  return -(this->compute_reactant_mass_flux(rho,Y_r,T));
107  }
libMesh::Real compute_reactant_mass_flux(const libMesh::Real rho, const libMesh::Real Y_r, const libMesh::Real T)
libMesh::Real rho(libMesh::Real T, libMesh::Real p0, libMesh::Real R_mix) const
template<typename Chemistry >
libMesh::Real GRINS::GasRecombinationCatalyticWall< Chemistry >::compute_reactant_mass_flux ( const libMesh::Real  rho,
const libMesh::Real  Y_r,
const libMesh::Real  T 
)
inline

Definition at line 89 of file gas_recombination_catalytic_wall.h.

92  {
93  const libMesh::Real rho_r = rho*Y_r;
94 
95  const libMesh::Real omega_dot = this->omega_dot( rho_r, T );
96 
97  return -omega_dot;
98  }
libMesh::Real omega_dot(const libMesh::Real rho_s, const libMesh::Real T) const
libMesh::Real rho(libMesh::Real T, libMesh::Real p0, libMesh::Real R_mix) const
template<typename Chemistry >
bool GRINS::GasRecombinationCatalyticWall< Chemistry >::eval_flux ( bool  compute_jacobian,
AssemblyContext context,
libMesh::Real  sign,
bool  is_axisymmetric 
)
virtual

Implements GRINS::NeumannBCAbstract.

Definition at line 81 of file gas_recombination_catalytic_wall.C.

85  {
86  libMesh::FEGenericBase<libMesh::Real>* side_fe = NULL;
87  context.get_side_fe( _reactant_var_idx, side_fe );
88 
89  // The number of local degrees of freedom in each variable.
90  const unsigned int n_var_dofs = context.get_dof_indices(_reactant_var_idx).size();
91 
92  libmesh_assert_equal_to( n_var_dofs, context.get_dof_indices(_product_var_idx).size() );
93 
94  // Element Jacobian * quadrature weight for side integration.
95  const std::vector<libMesh::Real> &JxW_side = side_fe->get_JxW();
96 
97  // The var shape functions at side quadrature points.
98  const std::vector<std::vector<libMesh::Real> >& var_phi_side = side_fe->get_phi();
99 
100  // Physical location of the quadrature points
101  const std::vector<libMesh::Point>& var_qpoint = side_fe->get_xyz();
102 
103  // reactant residual
104  libMesh::DenseSubVector<libMesh::Number> &F_r_var = context.get_elem_residual(_reactant_var_idx);
105 
106  // product residual
107  libMesh::DenseSubVector<libMesh::Number> &F_p_var = context.get_elem_residual(_product_var_idx);
108 
109  unsigned int n_qpoints = context.get_side_qrule().n_points();
110 
111  for (unsigned int qp=0; qp != n_qpoints; qp++)
112  {
113  libMesh::Real jac = JxW_side[qp];
114 
115  if(is_axisymmetric)
116  {
117  const libMesh::Number r = var_qpoint[qp](0);
118  jac *= r;
119  }
120 
121  std::vector<libMesh::Real> mass_fractions(this->_chem_ptr->n_species());
122  for( unsigned int s = 0; s < this->_chem_ptr->n_species(); s++ )
123  mass_fractions[s] = context.side_value(this->_species_vars[s], qp);
124 
125  libMesh::Real Y_r = mass_fractions[this->_reactant_species_idx];
126  libMesh::Real T = context.side_value(this->_T_var, qp);
127  libMesh::Real R_mix = this->_chem_ptr->R_mix(mass_fractions);
128  libMesh::Real rho = this->rho( T, this->_p0, R_mix );
129 
130  const libMesh::Real r_value = this->compute_reactant_mass_flux(rho, Y_r, T);
131 
132  const libMesh::Real p_value = -r_value;
133 
134  for (unsigned int i=0; i != n_var_dofs; i++)
135  {
136  F_r_var(i) += sign*r_value*var_phi_side[i][qp]*jac;
137 
138  F_p_var(i) += sign*p_value*var_phi_side[i][qp]*jac;
139 
140  if( compute_jacobian )
141  libmesh_not_implemented();
142  }
143  }
144 
145  // We're not computing the Jacobian yet
146  return false;
147  }
libMesh::Real _p0
Thermodynamic pressure.
SharedPtr< Chemistry > _chem_ptr
libMesh::Real compute_reactant_mass_flux(const libMesh::Real rho, const libMesh::Real Y_r, const libMesh::Real T)
libMesh::Real rho(libMesh::Real T, libMesh::Real p0, libMesh::Real R_mix) const
template<typename Chemistry >
void GRINS::GasRecombinationCatalyticWall< Chemistry >::init ( const libMesh::FEMSystem &  system)
virtual

Deprecated.

Reimplemented from GRINS::CatalyticWallBase< Chemistry >.

Definition at line 64 of file gas_recombination_catalytic_wall.C.

65  {
66  libmesh_do_once(libmesh_deprecated());
67 
68  const std::string r_var_name = std::string("w_"+this->_chemistry.species_name( this->_reactant_species_idx ) );
69 
70  const std::string p_var_name = std::string("w_"+this->_chemistry.species_name( this->_product_species_idx ) );
71 
72  libmesh_assert( system.has_variable( r_var_name ) );
73  libmesh_assert( system.has_variable( p_var_name ) );
74 
75  this->_reactant_var_idx = system.variable_number( r_var_name );
76 
77  this->_product_var_idx = system.variable_number( p_var_name );
78  }
const Chemistry & _chemistry
Deprecated.

Member Data Documentation

template<typename Chemistry>
const unsigned int GRINS::GasRecombinationCatalyticWall< Chemistry >::_product_species_idx
protected

Definition at line 82 of file gas_recombination_catalytic_wall.h.

template<typename Chemistry>
VariableIndex GRINS::GasRecombinationCatalyticWall< Chemistry >::_product_var_idx
protected

Definition at line 83 of file gas_recombination_catalytic_wall.h.

template<typename Chemistry>
const unsigned int GRINS::GasRecombinationCatalyticWall< Chemistry >::_reactant_species_idx
protected

Definition at line 79 of file gas_recombination_catalytic_wall.h.

template<typename Chemistry>
VariableIndex GRINS::GasRecombinationCatalyticWall< Chemistry >::_reactant_var_idx
protected

Definition at line 80 of file gas_recombination_catalytic_wall.h.


The documentation for this class was generated from the following files:

Generated on Tue Dec 19 2017 12:47:30 for GRINS-0.8.0 by  doxygen 1.8.9.1