GRINS-0.7.0
List of all members | Public Member Functions | Protected Attributes | Private Member Functions
GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity > Class Template Reference

Adds Boussinesq bouyancy SPGSM stabilization source term. More...

#include <boussinesq_buoyancy_spgsm_stab.h>

Inheritance diagram for GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >:
Inheritance graph
[legend]
Collaboration diagram for GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >:
Collaboration graph
[legend]

Public Member Functions

 BoussinesqBuoyancySPGSMStabilization (const std::string &physics_name, const GetPot &input)
 
 ~BoussinesqBuoyancySPGSMStabilization ()
 
virtual void element_time_derivative (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Time dependent part(s) of physics for element interiors. More...
 
virtual void element_constraint (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Constraint part(s) of physics for element interiors. More...
 
virtual void mass_residual (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Mass matrix part(s) for element interiors. All boundary terms lie within the time_derivative part. More...
 
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::BoussinesqBuoyancyBase
 BoussinesqBuoyancyBase (const std::string &physics_name, const GetPot &input)
 
 ~BoussinesqBuoyancyBase ()
 
virtual void init_variables (libMesh::FEMSystem *system)
 Initialization of BoussinesqBuoyancy variables. More...
 
- Public Member Functions inherited from GRINS::Physics
 Physics (const GRINS::PhysicsName &physics_name, const GetPot &input)
 
virtual ~Physics ()
 
virtual bool enabled_on_elem (const libMesh::Elem *elem)
 Find if current physics is active on supplied element. More...
 
void set_is_steady (bool is_steady)
 Sets whether this physics is to be solved with a steady solver or not. More...
 
bool is_steady () const
 Returns whether or not this physics is being solved with a steady solver. More...
 
virtual void set_time_evolving_vars (libMesh::FEMSystem *system)
 Set which variables are time evolving. More...
 
virtual void auxiliary_init (MultiphysicsSystem &system)
 Any auxillary initialization a Physics class may need. More...
 
virtual void register_postprocessing_vars (const GetPot &input, PostProcessedQuantities< libMesh::Real > &postprocessing)
 Register name of postprocessed quantity with PostProcessedQuantities. More...
 
virtual void init_context (AssemblyContext &context)
 Initialize context for added physics variables. More...
 
virtual void side_time_derivative (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Time dependent part(s) of physics for boundaries of elements on the domain boundary. More...
 
virtual void nonlocal_time_derivative (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Time dependent part(s) of physics for scalar variables. More...
 
virtual void side_constraint (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Constraint part(s) of physics for boundaries of elements on the domain boundary. More...
 
virtual void nonlocal_constraint (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Constraint part(s) of physics for scalar variables. More...
 
virtual void damping_residual (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Damping matrix part(s) for element interiors. All boundary terms lie within the time_derivative part. More...
 
virtual void nonlocal_mass_residual (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Mass matrix part(s) for scalar variables. More...
 
void init_ics (libMesh::FEMSystem *system, libMesh::CompositeFunction< libMesh::Number > &all_ics)
 
virtual void compute_element_time_derivative_cache (const AssemblyContext &context, CachedValues &cache)
 
virtual void compute_side_time_derivative_cache (const AssemblyContext &context, CachedValues &cache)
 
virtual void compute_nonlocal_time_derivative_cache (const AssemblyContext &context, CachedValues &cache)
 
virtual void compute_element_constraint_cache (const AssemblyContext &context, CachedValues &cache)
 
virtual void compute_side_constraint_cache (const AssemblyContext &context, CachedValues &cache)
 
virtual void compute_nonlocal_constraint_cache (const AssemblyContext &context, CachedValues &cache)
 
virtual void compute_damping_residual_cache (const AssemblyContext &context, CachedValues &cache)
 
virtual void compute_mass_residual_cache (const AssemblyContext &context, CachedValues &cache)
 
virtual void compute_nonlocal_mass_residual_cache (const AssemblyContext &context, CachedValues &cache)
 
virtual void compute_postprocessed_quantity (unsigned int quantity_index, const AssemblyContext &context, const libMesh::Point &point, libMesh::Real &value)
 
ICHandlingBaseget_ic_handler ()
 
- 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

IncompressibleNavierStokesStabilizationHelper _flow_stab_helper
 
HeatTransferStabilizationHelper _temp_stab_helper
 
libMesh::Number _Cp
 
libMesh::Number _k
 
Viscosity _mu
 Viscosity object. More...
 
- Protected Attributes inherited from GRINS::BoussinesqBuoyancyBase
VelocityFEVariables _flow_vars
 
PressureFEVariable _press_var
 
PrimitiveTempFEVariables _temp_vars
 
libMesh::Number _rho
 $ \rho = $ density More...
 
libMesh::Number _T_ref
 $ T_0 = $ reference temperature More...
 
libMesh::Number _beta_T
 $ \beta_T = $ coefficient of thermal expansion More...
 
libMesh::Point _g
 Gravitational vector. More...
 
unsigned int _dim
 Physical dimension of problem. More...
 
- Protected Attributes inherited from GRINS::Physics
const PhysicsName _physics_name
 Name of the physics object. Used for reading physics specific inputs. More...
 
GRINS::ICHandlingBase_ic_handler
 
std::set< libMesh::subdomain_id_type > _enabled_subdomains
 Subdomains on which the current Physics class is enabled. More...
 

Private Member Functions

 BoussinesqBuoyancySPGSMStabilization ()
 

Additional Inherited Members

- Static Public Member Functions inherited from GRINS::Physics
static void set_is_axisymmetric (bool is_axisymmetric)
 Set whether we should treat the problem as axisymmetric. More...
 
static bool is_axisymmetric ()
 
- 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::BoussinesqBuoyancyBase
void read_property (const GetPot &input, const std::string &old_option, const std::string &property, libMesh::Real &value)
 Helper function for parsing/maintaing backward compatibility. More...
 
void no_input_warning (const GetPot &input, const std::string &old_option, const std::string &material, const std::string &property)
 Helper function for parsing/maintaing backward compatibility. More...
 
- Protected Member Functions inherited from GRINS::Physics
libMesh::UniquePtr< libMesh::FEGenericBase< libMesh::Real > > build_new_fe (const libMesh::Elem *elem, const libMesh::FEGenericBase< libMesh::Real > *fe, const libMesh::Point p)
 
void parse_enabled_subdomains (const GetPot &input, const std::string &physics_name)
 
- Static Protected Attributes inherited from GRINS::Physics
static bool _is_steady = false
 Caches whether or not the solver that's being used is steady or not. More...
 
static bool _is_axisymmetric = false
 Caches whether we are solving an axisymmetric problem or not. More...
 

Detailed Description

template<class Viscosity>
class GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >

Adds Boussinesq bouyancy SPGSM stabilization source term.

This class implements the SPGSM stabilization term for the BoussinesqBuoyancy Physics. Intended to be used with IncompressibleNavierStokesSPGSMStabilization and HeatTransferSPGSMStabilization. This, essentially, is just adding SUPG terms to the momentum equation and PSPG terms to the pressure equation.

Definition at line 45 of file boussinesq_buoyancy_spgsm_stab.h.

Constructor & Destructor Documentation

template<class Mu >
GRINS::BoussinesqBuoyancySPGSMStabilization< Mu >::BoussinesqBuoyancySPGSMStabilization ( const std::string &  physics_name,
const GetPot &  input 
)

Definition at line 42 of file boussinesq_buoyancy_spgsm_stab.C.

References GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::_Cp, GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::_k, GRINS::PhysicsNaming::boussinesq_buoyancy(), GRINS::PhysicsNaming::heat_transfer(), GRINS::MaterialsParsing::read_specific_heat(), and GRINS::ParameterUser::set_parameter().

43  : BoussinesqBuoyancyBase(physics_name,input),
44  _flow_stab_helper(physics_name+"FlowStabHelper", input),
45  _temp_stab_helper(physics_name+"TempStabHelper", input),
46  _Cp(0.0),
47  _k(1.0),
49  {
51 
52  this->set_parameter
53  (_k, input, "Physics/"+PhysicsNaming::heat_transfer()+"/k", _k);
54  }
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.
static PhysicsName heat_transfer()
IncompressibleNavierStokesStabilizationHelper _flow_stab_helper
static void read_specific_heat(const std::string &core_physics_name, const GetPot &input, ParameterUser &params, libMesh::Real &cp)
Helper function to reading scalar specific heat from input.
static PhysicsName boussinesq_buoyancy()
static std::string material_name(const GetPot &input, const std::string &physics)
Get the name of the material in the Physics/physics section.

Definition at line 57 of file boussinesq_buoyancy_spgsm_stab.C.

58  {
59  return;
60  }
template<class Viscosity >
GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::BoussinesqBuoyancySPGSMStabilization ( )
private

Member Function Documentation

template<class Mu >
void GRINS::BoussinesqBuoyancySPGSMStabilization< Mu >::element_constraint ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtual

Constraint part(s) of physics for element interiors.

Reimplemented from GRINS::Physics.

Definition at line 162 of file boussinesq_buoyancy_spgsm_stab.C.

165  {
166 #ifdef GRINS_USE_GRVY_TIMERS
167  this->_timer->BeginTimer("BoussinesqBuoyancySPGSMStabilization::element_constraint");
168 #endif
169 
170  // The number of local degrees of freedom in each variable.
171  const unsigned int n_p_dofs = context.get_dof_indices(_press_var.p()).size();
172 
173  // Element Jacobian * quadrature weights for interior integration.
174  const std::vector<libMesh::Real> &JxW =
175  context.get_element_fe(_flow_vars.u())->get_JxW();
176 
177  const std::vector<std::vector<libMesh::RealGradient> >& p_dphi =
178  context.get_element_fe(this->_press_var.p())->get_dphi();
179 
180  libMesh::DenseSubVector<libMesh::Number> &Fp = context.get_elem_residual(this->_press_var.p()); // R_{p}
181 
182  // Now we will build the element Jacobian and residual.
183  // Constructing the residual requires the solution and its
184  // gradient from the previous timestep. This must be
185  // calculated at each quadrature point by summing the
186  // solution degree-of-freedom values by the appropriate
187  // weight functions.
188  unsigned int n_qpoints = context.get_element_qrule().n_points();
189 
190  libMesh::FEBase* fe = context.get_element_fe(this->_flow_vars.u());
191 
192  for (unsigned int qp=0; qp != n_qpoints; qp++)
193  {
194  libMesh::RealGradient g = this->_flow_stab_helper.compute_g( fe, context, qp );
195  libMesh::RealTensor G = this->_flow_stab_helper.compute_G( fe, context, qp );
196 
197  libMesh::RealGradient U( context.interior_value( this->_flow_vars.u(), qp ),
198  context.interior_value( this->_flow_vars.v(), qp ) );
199  if( this->_dim == 3 )
200  {
201  U(2) = context.interior_value( this->_flow_vars.w(), qp );
202  }
203 
204  // Compute the viscosity at this qp
205  libMesh::Real mu_qp = this->_mu(context, qp);
206 
207  libMesh::Real tau_M = this->_flow_stab_helper.compute_tau_momentum( context, qp, g, G, this->_rho, U, mu_qp, this->_is_steady );
208 
209  // Compute the solution & its gradient at the old Newton iterate.
210  libMesh::Number T;
211  T = context.interior_value(_temp_vars.T(), qp);
212 
213  libMesh::RealGradient residual = _rho*_beta_T*(T-_T_ref)*_g;
214 
215  // First, an i-loop over the velocity degrees of freedom.
216  // We know that n_u_dofs == n_v_dofs so we can compute contributions
217  // for both at the same time.
218  for (unsigned int i=0; i != n_p_dofs; i++)
219  {
220  Fp(i) += -tau_M*residual*p_dphi[i][qp]*JxW[qp];
221 
222  if (compute_jacobian)
223  {
224  libmesh_not_implemented();
225  } // End compute_jacobian check
226 
227  } // End i dof loop
228 
229  } // End quadrature loop
230 
231 #ifdef GRINS_USE_GRVY_TIMERS
232  this->_timer->EndTimer("BoussinesqBuoyancySPGSMStabilization::element_constraint");
233 #endif
234 
235  return;
236  }
libMesh::Point _g
Gravitational vector.
libMesh::RealTensor compute_G(libMesh::FEBase *fe, AssemblyContext &c, unsigned int qp) const
Definition: stab_helper.C:64
unsigned int _dim
Physical dimension of problem.
libMesh::RealGradient compute_g(libMesh::FEBase *fe, AssemblyContext &c, unsigned int qp) const
Definition: stab_helper.C:47
VariableIndex p() const
IncompressibleNavierStokesStabilizationHelper _flow_stab_helper
libMesh::Number _T_ref
reference temperature
static bool _is_steady
Caches whether or not the solver that's being used is steady or not.
Definition: physics.h:277
libMesh::Number _beta_T
coefficient of thermal expansion
libMesh::Real compute_tau_momentum(AssemblyContext &c, unsigned int qp, libMesh::RealGradient &g, libMesh::RealTensor &G, libMesh::Real rho, libMesh::Gradient U, libMesh::Real mu, bool is_steady) const
PrimitiveTempFEVariables _temp_vars
template<class Mu >
void GRINS::BoussinesqBuoyancySPGSMStabilization< Mu >::element_time_derivative ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtual

Time dependent part(s) of physics for element interiors.

Reimplemented from GRINS::Physics.

Definition at line 63 of file boussinesq_buoyancy_spgsm_stab.C.

66  {
67 #ifdef GRINS_USE_GRVY_TIMERS
68  this->_timer->BeginTimer("BoussinesqBuoyancySPGSMStabilization::element_time_derivative");
69 #endif
70 
71  // The number of local degrees of freedom in each variable.
72  const unsigned int n_u_dofs = context.get_dof_indices(_flow_vars.u()).size();
73 
74  // Element Jacobian * quadrature weights for interior integration.
75  const std::vector<libMesh::Real> &JxW =
76  context.get_element_fe(_flow_vars.u())->get_JxW();
77 
78  /*
79  const std::vector<std::vector<libMesh::Real> >& u_phi =
80  context.get_element_fe(this->_flow_vars.u())->get_phi();
81  */
82 
83  const std::vector<std::vector<libMesh::RealGradient> >& u_gradphi =
84  context.get_element_fe(this->_flow_vars.u())->get_dphi();
85 
86  // Get residuals
87  libMesh::DenseSubVector<libMesh::Number> &Fu = context.get_elem_residual(_flow_vars.u()); // R_{u}
88  libMesh::DenseSubVector<libMesh::Number> &Fv = context.get_elem_residual(_flow_vars.v()); // R_{v}
89  libMesh::DenseSubVector<libMesh::Number> *Fw = NULL;
90  if(this->_dim == 3)
91  {
92  Fw = &context.get_elem_residual(this->_flow_vars.w()); // R_{w}
93  }
94 
95  // Now we will build the element Jacobian and residual.
96  // Constructing the residual requires the solution and its
97  // gradient from the previous timestep. This must be
98  // calculated at each quadrature point by summing the
99  // solution degree-of-freedom values by the appropriate
100  // weight functions.
101  unsigned int n_qpoints = context.get_element_qrule().n_points();
102 
103  libMesh::FEBase* fe = context.get_element_fe(this->_flow_vars.u());
104 
105  for (unsigned int qp=0; qp != n_qpoints; qp++)
106  {
107  libMesh::RealGradient g = this->_flow_stab_helper.compute_g( fe, context, qp );
108  libMesh::RealTensor G = this->_flow_stab_helper.compute_G( fe, context, qp );
109 
110  libMesh::RealGradient U( context.interior_value( this->_flow_vars.u(), qp ),
111  context.interior_value( this->_flow_vars.v(), qp ) );
112  if( this->_dim == 3 )
113  {
114  U(2) = context.interior_value( this->_flow_vars.w(), qp );
115  }
116 
117  // Compute the viscosity at this qp
118  libMesh::Real mu_qp = this->_mu(context, qp);
119 
120  libMesh::Real tau_M = this->_flow_stab_helper.compute_tau_momentum( context, qp, g, G, this->_rho, U, mu_qp, this->_is_steady );
121 
122  //libMesh::Real tau_E = this->_temp_stab_helper.compute_tau_energy( context, G, _rho, _Cp, _k, U, this->_is_steady );
123 
124  //libMesh::Real RE = this->_temp_stab_helper.compute_res_energy_steady( context, qp, _rho, _Cp, _k );
125 
126  // Compute the solution & its gradient at the old Newton iterate.
127  libMesh::Number T;
128  T = context.interior_value(_temp_vars.T(), qp);
129 
130  libMesh::RealGradient residual = _rho*_beta_T*(T-_T_ref)*_g;
131 
132  for (unsigned int i=0; i != n_u_dofs; i++)
133  {
134  Fu(i) += ( -tau_M*residual(0)*_rho*U*u_gradphi[i][qp] )*JxW[qp];
135  // + _rho*_beta_T*tau_E*RE*_g(0)*u_phi[i][qp] )*JxW[qp];
136 
137  Fv(i) += ( -tau_M*residual(1)*_rho*U*u_gradphi[i][qp] )*JxW[qp];
138  // + _rho*_beta_T*tau_E*RE*_g(1)*u_phi[i][qp] )*JxW[qp];
139 
140  if (_dim == 3)
141  {
142  (*Fw)(i) += ( -tau_M*residual(2)*_rho*U*u_gradphi[i][qp] )*JxW[qp];
143  // + _rho*_beta_T*tau_E*RE*_g(2)*u_phi[i][qp] )*JxW[qp];
144  }
145 
146  if (compute_jacobian)
147  {
148  libmesh_not_implemented();
149  } // End compute_jacobian check
150 
151  } // End i dof loop
152  } // End quadrature loop
153 
154 #ifdef GRINS_USE_GRVY_TIMERS
155  this->_timer->EndTimer("BoussinesqBuoyancySPGSMStabilization::element_time_derivative");
156 #endif
157 
158  return;
159  }
libMesh::Point _g
Gravitational vector.
libMesh::RealTensor compute_G(libMesh::FEBase *fe, AssemblyContext &c, unsigned int qp) const
Definition: stab_helper.C:64
unsigned int _dim
Physical dimension of problem.
libMesh::RealGradient compute_g(libMesh::FEBase *fe, AssemblyContext &c, unsigned int qp) const
Definition: stab_helper.C:47
IncompressibleNavierStokesStabilizationHelper _flow_stab_helper
libMesh::Number _T_ref
reference temperature
static bool _is_steady
Caches whether or not the solver that's being used is steady or not.
Definition: physics.h:277
libMesh::Number _beta_T
coefficient of thermal expansion
libMesh::Real compute_tau_momentum(AssemblyContext &c, unsigned int qp, libMesh::RealGradient &g, libMesh::RealTensor &G, libMesh::Real rho, libMesh::Gradient U, libMesh::Real mu, bool is_steady) const
PrimitiveTempFEVariables _temp_vars
template<class Mu >
void GRINS::BoussinesqBuoyancySPGSMStabilization< Mu >::mass_residual ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtual

Mass matrix part(s) for element interiors. All boundary terms lie within the time_derivative part.

Reimplemented from GRINS::Physics.

Definition at line 239 of file boussinesq_buoyancy_spgsm_stab.C.

242  {
243  /*
244 #ifdef GRINS_USE_GRVY_TIMERS
245  this->_timer->BeginTimer("BoussinesqBuoyancySPGSMStabilization::mass_residual");
246 #endif
247 
248  // The number of local degrees of freedom in each variable.
249  const unsigned int n_u_dofs = context.get_dof_indices(_flow_vars.u()).size();
250 
251  // Element Jacobian * quadrature weights for interior integration.
252  const std::vector<libMesh::Real> &JxW =
253  context.get_element_fe(_flow_vars.u())->get_JxW();
254 
255  const std::vector<std::vector<libMesh::Real> >& u_phi =
256  context.get_element_fe(this->_flow_vars.u())->get_phi();
257 
258  // Get residuals
259  libMesh::DenseSubVector<libMesh::Number> &Fu = context.get_elem_residual(_flow_vars.u()); // R_{u}
260  libMesh::DenseSubVector<libMesh::Number> &Fv = context.get_elem_residual(_flow_vars.v()); // R_{v}
261  libMesh::DenseSubVector<libMesh::Number> *Fw = NULL;
262  if(this->_dim == 3)
263  {
264  Fw = &context.get_elem_residual(this->_flow_vars.w()); // R_{w}
265  }
266 
267  // Now we will build the element Jacobian and residual.
268  // Constructing the residual requires the solution and its
269  // gradient from the previous timestep. This must be
270  // calculated at each quadrature point by summing the
271  // solution degree-of-freedom values by the appropriate
272  // weight functions.
273  unsigned int n_qpoints = context.get_element_qrule().n_points();
274 
275  libMesh::FEBase* fe = context.get_element_fe(this->_flow_vars.u());
276 
277  for (unsigned int qp=0; qp != n_qpoints; qp++)
278  {
279  libMesh::RealGradient g = this->_flow_stab_helper.compute_g( fe, context, qp );
280  libMesh::RealTensor G = this->_flow_stab_helper.compute_G( fe, context, qp );
281 
282  libMesh::RealGradient U( context.fixed_interior_value( this->_flow_vars.u(), qp ),
283  context.fixed_interior_value( this->_flow_vars.v(), qp ) );
284  if( this->_dim == 3 )
285  {
286  U(2) = context.fixed_interior_value( this->_flow_vars.w(), qp );
287  }
288 
289  libMesh::Real tau_E = this->_temp_stab_helper.compute_tau_energy( context, G, _rho, _Cp, _k, U, false );
290 
291  libMesh::Real RE = this->_temp_stab_helper.compute_res_energy_transient( context, qp, _rho, _Cp );
292 
293 
294  for (unsigned int i=0; i != n_u_dofs; i++)
295  {
296  Fu(i) += -_rho*_beta_T*tau_E*RE*_g(0)*u_phi[i][qp]*JxW[qp];
297 
298  Fv(i) += -_rho*_beta_T*tau_E*RE*_g(1)*u_phi[i][qp]*JxW[qp];
299 
300  if (_dim == 3)
301  {
302  (*Fw)(i) += -_rho*_beta_T*tau_E*RE*_g(2)*u_phi[i][qp]*JxW[qp];
303  }
304 
305  if (compute_jacobian)
306  {
307  libmesh_not_implemented();
308  } // End compute_jacobian check
309 
310  } // End i dof loop
311  } // End quadrature loop
312 
313 #ifdef GRINS_USE_GRVY_TIMERS
314  this->_timer->EndTimer("BoussinesqBuoyancySPGSMStabilization::mass_residual");
315 #endif
316  */
317 
318  return;
319  }
template<class Mu >
void GRINS::BoussinesqBuoyancySPGSMStabilization< Mu >::register_parameter ( const std::string &  param_name,
libMesh::ParameterMultiAccessor< libMesh::Number > &  param_pointer 
) const
virtual

Each subclass will register its copy of an independent.

Reimplemented from GRINS::ParameterUser.

Definition at line 323 of file boussinesq_buoyancy_spgsm_stab.C.

References GRINS::ParameterUser::register_parameter().

326  {
327  ParameterUser::register_parameter(param_name, param_pointer);
328  _mu.register_parameter(param_name, param_pointer);
329  }
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.

Member Data Documentation

template<class Viscosity >
libMesh::Number GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::_Cp
protected
template<class Viscosity >
IncompressibleNavierStokesStabilizationHelper GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::_flow_stab_helper
protected

Definition at line 76 of file boussinesq_buoyancy_spgsm_stab.h.

template<class Viscosity >
libMesh::Number GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::_k
protected
template<class Viscosity >
Viscosity GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::_mu
protected

Viscosity object.

Definition at line 82 of file boussinesq_buoyancy_spgsm_stab.h.

template<class Viscosity >
HeatTransferStabilizationHelper GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::_temp_stab_helper
protected

Definition at line 77 of file boussinesq_buoyancy_spgsm_stab.h.


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

Generated on Thu Jun 2 2016 21:52:30 for GRINS-0.7.0 by  doxygen 1.8.10