GRINS-0.8.0
gas_solid_catalytic_wall.C
Go to the documentation of this file.
1 //-----------------------------------------------------------------------bl-
2 //--------------------------------------------------------------------------
3 //
4 // GRINS - General Reacting Incompressible Navier-Stokes
5 //
6 // Copyright (C) 2014-2017 Paul T. Bauman, Roy H. Stogner
7 // Copyright (C) 2010-2013 The PECOS Development Team
8 //
9 // This library is free software; you can redistribute it and/or
10 // modify it under the terms of the Version 2.1 GNU Lesser General
11 // Public License as published by the Free Software Foundation.
12 //
13 // This library is distributed in the hope that it will be useful,
14 // but WITHOUT ANY WARRANTY; without even the implied warranty of
15 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16 // Lesser General Public License for more details.
17 //
18 // You should have received a copy of the GNU Lesser General Public
19 // License along with this library; if not, write to the Free Software
20 // Foundation, Inc. 51 Franklin Street, Fifth Floor,
21 // Boston, MA 02110-1301 USA
22 //
23 //-----------------------------------------------------------------------el-
24 
25 // This class
27 
28 // GRINS
29 #include "grins/physics.h"
30 
31 // libMesh
32 #include "libmesh/quadrature.h"
33 
34 namespace GRINS
35 {
36  template<typename Chemistry>
38  SharedPtr<CatalycityBase>& gamma,
39  const std::vector<VariableIndex>& species_vars,
40  VariableIndex T_var,
41  libMesh::Real p0,
42  unsigned int reactant_gas_species_idx,
43  unsigned int reactant_solid_species_idx,
44  unsigned int product_species_idx)
45  : CatalyticWallBase<Chemistry>(chem,gamma,species_vars,T_var,p0,reactant_gas_species_idx),
46  _reactant_gas_species_idx(reactant_gas_species_idx),
47  _reactant_gas_var_idx(species_vars[reactant_gas_species_idx]),
48  _reactant_solid_species_idx(reactant_solid_species_idx),
49  _product_species_idx(product_species_idx),
50  _product_var_idx(species_vars[product_species_idx])
51  {}
52 
53  template<typename Chemistry>
55  CatalycityBase& gamma,
56  const unsigned int reactant_gas_species_idx,
57  const unsigned int reactant_solid_species_idx,
58  const unsigned int product_species_idx )
59  : CatalyticWallBase<Chemistry>(chem_mixture,gamma,reactant_gas_species_idx),
60  _reactant_gas_species_idx(reactant_gas_species_idx),
61  _reactant_solid_species_idx(reactant_solid_species_idx),
62  _product_species_idx(product_species_idx)
63  {
64  libmesh_deprecated();
65  }
66 
67  template<typename Chemistry>
68  void GasSolidCatalyticWall<Chemistry>::init( const libMesh::FEMSystem& system )
69  {
70  libmesh_do_once(libmesh_deprecated());
71 
72  const std::string r_var_name = std::string("w_"+this->_chemistry.species_name( this->_reactant_gas_species_idx ) );
73 
74  const std::string p_var_name = std::string("w_"+this->_chemistry.species_name( this->_product_species_idx ) );
75 
76  libmesh_assert( system.has_variable( r_var_name ) );
77  libmesh_assert( system.has_variable( p_var_name ) );
78 
79  this->_reactant_gas_var_idx = system.variable_number( r_var_name );
80 
81  this->_product_var_idx = system.variable_number( p_var_name );
82 
83  return;
84  }
85 
86  template<typename Chemistry>
87  bool GasSolidCatalyticWall<Chemistry>::eval_flux( bool compute_jacobian,
88  AssemblyContext& context,
89  libMesh::Real sign,
90  bool is_axisymmetric )
91  {
92  libMesh::FEGenericBase<libMesh::Real>* side_fe = NULL;
93  context.get_side_fe( _reactant_gas_var_idx, side_fe );
94 
95  // The number of local degrees of freedom in each variable.
96  const unsigned int n_var_dofs = context.get_dof_indices(_reactant_gas_var_idx).size();
97 
98  libmesh_assert_equal_to( n_var_dofs, context.get_dof_indices(_product_var_idx).size() );
99 
100  // Element Jacobian * quadrature weight for side integration.
101  const std::vector<libMesh::Real> &JxW_side = side_fe->get_JxW();
102 
103  // The var shape functions at side quadrature points.
104  const std::vector<std::vector<libMesh::Real> >& var_phi_side = side_fe->get_phi();
105 
106  // Physical location of the quadrature points
107  const std::vector<libMesh::Point>& var_qpoint = side_fe->get_xyz();
108 
109  // reactant residual
110  libMesh::DenseSubVector<libMesh::Number> &F_r_var = context.get_elem_residual(_reactant_gas_var_idx);
111 
112  // product residual
113  libMesh::DenseSubVector<libMesh::Number> &F_p_var = context.get_elem_residual(_product_var_idx);
114 
115  unsigned int n_qpoints = context.get_side_qrule().n_points();
116 
117  for (unsigned int qp=0; qp != n_qpoints; qp++)
118  {
119  libMesh::Real jac = JxW_side[qp];
120 
121  if(is_axisymmetric)
122  {
123  const libMesh::Number r = var_qpoint[qp](0);
124  jac *= r;
125  }
126 
127  std::vector<libMesh::Real> mass_fractions(this->_chem_ptr->n_species());
128  for( unsigned int s = 0; s < this->_chem_ptr->n_species(); s++ )
129  mass_fractions[s] = context.side_value(this->_species_vars[s], qp);
130 
131  libMesh::Real T = context.side_value(this->_T_var, qp);
132  libMesh::Real R_mix = this->_chem_ptr->R_mix(mass_fractions);
133  libMesh::Real rho = this->rho( T, this->_p0, R_mix );
134 
135  libMesh::Real Y_r = mass_fractions[this->_reactant_gas_species_idx];
136 
137  const libMesh::Real r_value = this->compute_reactant_gas_mass_flux(rho, Y_r, T);
138 
139  const libMesh::Real p_value = this->compute_product_mass_flux(rho, Y_r, T);
140 
141  for (unsigned int i=0; i != n_var_dofs; i++)
142  {
143  F_r_var(i) += sign*r_value*var_phi_side[i][qp]*jac;
144 
145  F_p_var(i) += sign*p_value*var_phi_side[i][qp]*jac;
146 
147  if( compute_jacobian )
148  {
149  libmesh_not_implemented();
150  }
151  }
152  }
153 
154  // We're not computing the Jacobian yet
155  return false;
156  }
157 
158  template<typename Chemistry>
160  const CachedValues& cache,
161  const bool request_jacobian )
162  {
163  libmesh_do_once(libmesh_deprecated());
164 
165  libMesh::FEGenericBase<libMesh::Real>* side_fe = NULL;
166  context.get_side_fe( _reactant_gas_var_idx, side_fe );
167 
168  // The number of local degrees of freedom in each variable.
169  const unsigned int n_var_dofs = context.get_dof_indices(_reactant_gas_var_idx).size();
170 
171  libmesh_assert_equal_to( n_var_dofs, context.get_dof_indices(_product_var_idx).size() );
172 
173  // Element Jacobian * quadrature weight for side integration.
174  const std::vector<libMesh::Real> &JxW_side = side_fe->get_JxW();
175 
176  // The var shape functions at side quadrature points.
177  const std::vector<std::vector<libMesh::Real> >& var_phi_side = side_fe->get_phi();
178 
179  // Physical location of the quadrature points
180  const std::vector<libMesh::Point>& var_qpoint = side_fe->get_xyz();
181 
182  // reactant residual
183  libMesh::DenseSubVector<libMesh::Number> &F_r_var = context.get_elem_residual(_reactant_gas_var_idx);
184 
185  // product residual
186  libMesh::DenseSubVector<libMesh::Number> &F_p_var = context.get_elem_residual(_product_var_idx);
187 
188  unsigned int n_qpoints = context.get_side_qrule().n_points();
189 
190  for (unsigned int qp=0; qp != n_qpoints; qp++)
191  {
192  libMesh::Real jac = JxW_side[qp];
193 
195  {
196  const libMesh::Number r = var_qpoint[qp](0);
197  jac *= r;
198  }
199 
200  const libMesh::Real rho = cache.get_cached_values(Cache::MIXTURE_DENSITY)[qp];
201 
202  const libMesh::Real Y_r = cache.get_cached_vector_values(Cache::MASS_FRACTIONS)[qp][this->_reactant_gas_species_idx];
203 
204  const libMesh::Real T = cache.get_cached_values(Cache::TEMPERATURE)[qp];
205 
206  const libMesh::Real r_value = this->compute_reactant_gas_mass_flux(rho, Y_r, T);
207 
208  const libMesh::Real p_value = this->compute_product_mass_flux(rho, Y_r, T);
209 
210  for (unsigned int i=0; i != n_var_dofs; i++)
211  {
212  F_r_var(i) += r_value*var_phi_side[i][qp]*jac;
213 
214  F_p_var(i) += p_value*var_phi_side[i][qp]*jac;
215 
216  if( request_jacobian )
217  {
218  libmesh_not_implemented();
219  }
220  }
221  }
222 
223  return;
224  }
225 
226 } // end namespace GRINS
static bool is_axisymmetric()
Definition: physics.h:132
unsigned int VariableIndex
More descriptive name of the type used for variable indices.
Definition: var_typedefs.h:42
virtual void init(const libMesh::FEMSystem &system)
Deprecated.
GRINS namespace.
const std::vector< std::vector< libMesh::Number > > & get_cached_vector_values(unsigned int quantity) const
virtual bool eval_flux(bool compute_jacobian, AssemblyContext &context, libMesh::Real sign, bool is_axisymmetric)
GasSolidCatalyticWall(SharedPtr< Chemistry > &chem, SharedPtr< CatalycityBase > &gamma, const std::vector< VariableIndex > &species_vars, VariableIndex T_var, libMesh::Real p0, unsigned int reactant_gas_species_idx, unsigned int reactant_solid_species_idx, unsigned int product_species_idx)
virtual void apply_fluxes(AssemblyContext &context, const CachedValues &cache, const bool request_jacobian)
Deprecated.
const std::vector< libMesh::Number > & get_cached_values(unsigned int quantity) const
Definition: cached_values.C:99

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