GRINS-0.6.0
Public Member Functions | Protected Attributes | Static Protected Attributes | Private Member Functions | List of all members
GRINS::AxisymmetricBoussinesqBuoyancy Class Reference

Adds Axisymmetric Boussinesq bouyancy source term. More...

#include <axisym_boussinesq_buoyancy.h>

Inheritance diagram for GRINS::AxisymmetricBoussinesqBuoyancy:
Inheritance graph
[legend]
Collaboration diagram for GRINS::AxisymmetricBoussinesqBuoyancy:
Collaboration graph
[legend]

Public Member Functions

 AxisymmetricBoussinesqBuoyancy (const std::string &physics_name, const GetPot &input)
 
 ~AxisymmetricBoussinesqBuoyancy ()
 
virtual void read_input_options (const GetPot &input)
 Read options from GetPot input file. More...
 
virtual void init_variables (libMesh::FEMSystem *system)
 Initialization of AxisymmetricBoussinesqBuoyancy variables. More...
 
virtual void element_time_derivative (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Source term contribution for AxisymmetricBoussinesqBuoyancy. More...
 
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 element_constraint (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Constraint part(s) of physics for element interiors. 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 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 nonlocal_mass_residual (bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
 Mass matrix part(s) for scalar variables. More...
 
void init_bcs (libMesh::FEMSystem *system)
 
void init_ics (libMesh::FEMSystem *system, libMesh::CompositeFunction< libMesh::Number > &all_ics)
 
void attach_neumann_bound_func (GRINS::NBCContainer &neumann_bcs)
 
void attach_dirichlet_bound_func (const GRINS::DBCContainer &dirichlet_bc)
 
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_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)
 
BCHandlingBaseget_bc_handler ()
 
ICHandlingBaseget_ic_handler ()
 
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 register_parameter (const std::string &param_name, libMesh::ParameterMultiPointer< libMesh::Number > &param_pointer) const
 Each subclass will register its copy of an independent. More...
 

Protected Attributes

unsigned int _dim
 Physical dimension of problem. More...
 
GRINSEnums::FEFamily _T_FE_family
 Element type, read from input. More...
 
GRINSEnums::FEFamily _V_FE_family
 
GRINSEnums::Order _T_order
 Temperature element order, read from input. More...
 
GRINSEnums::Order _V_order
 
VariableIndex _u_r_var
 Index for r-velocity field. More...
 
VariableIndex _u_z_var
 Index for z-velocity field. More...
 
VariableIndex _T_var
 Index for temperature field. More...
 
std::string _u_r_var_name
 Name of r-velocity. More...
 
std::string _u_z_var_name
 Name of z-velocity. More...
 
std::string _T_var_name
 Name of temperature. More...
 
libMesh::Number _rho_ref
 $ \rho_0 = $ reference 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...
 
const PhysicsName _physics_name
 Name of the physics object. Used for reading physics specific inputs. More...
 
GRINS::BCHandlingBase_bc_handler
 
GRINS::ICHandlingBase_ic_handler
 
std::set< libMesh::subdomain_id_type > _enabled_subdomains
 Subdomains on which the current Physics class is enabled. More...
 
bool _is_axisymmetric
 

Static Protected Attributes

static bool _is_steady = false
 Caches whether or not the solver that's being used is steady or not. More...
 

Private Member Functions

 AxisymmetricBoussinesqBuoyancy ()
 

Detailed Description

Adds Axisymmetric Boussinesq bouyancy source term.

This class implements the Axisymmetric Boussinesq approximation for thermal buoyancy. Namely: $ \mathbf{F} = -\rho_0 \beta_T \left( T - T_0 \right) \mathbf{g} $ where $ \rho_0 = $ reference density, $ T_0 = $ reference temperature, $ \beta_T = $ coefficient of thermal expansion, and $ \mathbf{g} = $ the gravitional vector. This source term is added to the governing flow equations through the element_time_derivative routine. This class requires an axisymmetric flow physics enabled and the AxisymmetricHeatTransfer physics class enabled.

Definition at line 55 of file axisym_boussinesq_buoyancy.h.

Constructor & Destructor Documentation

GRINS::AxisymmetricBoussinesqBuoyancy::AxisymmetricBoussinesqBuoyancy ( const std::string &  physics_name,
const GetPot &  input 
)

Definition at line 43 of file axisym_boussinesq_buoyancy.C.

References read_input_options().

45  : Physics(physics_name, input)
46  {
47  this->read_input_options(input);
48  return;
49  }
virtual void read_input_options(const GetPot &input)
Read options from GetPot input file.
GRINS::AxisymmetricBoussinesqBuoyancy::~AxisymmetricBoussinesqBuoyancy ( )

Definition at line 51 of file axisym_boussinesq_buoyancy.C.

52  {
53  return;
54  }
GRINS::AxisymmetricBoussinesqBuoyancy::AxisymmetricBoussinesqBuoyancy ( )
private

Member Function Documentation

void GRINS::Physics::attach_dirichlet_bound_func ( const GRINS::DBCContainer dirichlet_bc)
inherited

Definition at line 150 of file physics.C.

References GRINS::Physics::_bc_handler, and GRINS::BCHandlingBase::attach_dirichlet_bound_func().

151  {
152  _bc_handler->attach_dirichlet_bound_func( dirichlet_bc );
153  return;
154  }
GRINS::BCHandlingBase * _bc_handler
Definition: physics.h:256
void attach_dirichlet_bound_func(const GRINS::DBCContainer &dirichlet_bc)
void GRINS::Physics::attach_neumann_bound_func ( GRINS::NBCContainer neumann_bcs)
inherited

Definition at line 144 of file physics.C.

References GRINS::Physics::_bc_handler, and GRINS::BCHandlingBase::attach_neumann_bound_func().

145  {
146  _bc_handler->attach_neumann_bound_func( neumann_bcs );
147  return;
148  }
GRINS::BCHandlingBase * _bc_handler
Definition: physics.h:256
void attach_neumann_bound_func(GRINS::NBCContainer &neumann_bcs)
void GRINS::Physics::auxiliary_init ( MultiphysicsSystem system)
virtualinherited

Any auxillary initialization a Physics class may need.

This is called after all variables are added, so this method can safely query the MultiphysicsSystem about variable information.

Definition at line 113 of file physics.C.

114  {
115  return;
116  }
void GRINS::Physics::compute_element_constraint_cache ( const AssemblyContext context,
CachedValues cache 
)
virtualinherited

Definition at line 185 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::element_constraint().

187  {
188  return;
189  }
void GRINS::Physics::compute_element_time_derivative_cache ( const AssemblyContext context,
CachedValues cache 
)
virtualinherited
void GRINS::Physics::compute_mass_residual_cache ( const AssemblyContext context,
CachedValues cache 
)
virtualinherited

Definition at line 203 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::mass_residual().

205  {
206  return;
207  }
void GRINS::Physics::compute_nonlocal_constraint_cache ( const AssemblyContext context,
CachedValues cache 
)
virtualinherited

Definition at line 197 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::nonlocal_constraint().

199  {
200  return;
201  }
void GRINS::Physics::compute_nonlocal_mass_residual_cache ( const AssemblyContext context,
CachedValues cache 
)
virtualinherited

Definition at line 209 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::nonlocal_mass_residual().

211  {
212  return;
213  }
void GRINS::Physics::compute_nonlocal_time_derivative_cache ( const AssemblyContext context,
CachedValues cache 
)
virtualinherited

Definition at line 179 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::nonlocal_time_derivative().

181  {
182  return;
183  }
void GRINS::Physics::compute_postprocessed_quantity ( unsigned int  quantity_index,
const AssemblyContext context,
const libMesh::Point &  point,
libMesh::Real &  value 
)
virtualinherited
void GRINS::Physics::compute_side_constraint_cache ( const AssemblyContext context,
CachedValues cache 
)
virtualinherited

Definition at line 191 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::side_constraint().

193  {
194  return;
195  }
void GRINS::Physics::compute_side_time_derivative_cache ( const AssemblyContext context,
CachedValues cache 
)
virtualinherited

Reimplemented in GRINS::ReactingLowMachNavierStokes< Mixture, Evaluator >.

Definition at line 173 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::side_time_derivative().

175  {
176  return;
177  }
void GRINS::Physics::element_constraint ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtualinherited
void GRINS::AxisymmetricBoussinesqBuoyancy::element_time_derivative ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtual

Source term contribution for AxisymmetricBoussinesqBuoyancy.

This is the main part of the class. This will add the source term to the AxisymmetricIncompNavierStokes class.

Reimplemented from GRINS::Physics.

Definition at line 101 of file axisym_boussinesq_buoyancy.C.

References _beta_T, _g, _rho_ref, _T_ref, _T_var, _u_r_var, and _u_z_var.

104  {
105 #ifdef GRINS_USE_GRVY_TIMERS
106  this->_timer->BeginTimer("AxisymmetricBoussinesqBuoyancy::element_time_derivative");
107 #endif
108 
109  // The number of local degrees of freedom in each variable.
110  const unsigned int n_u_dofs = context.get_dof_indices(_u_r_var).size();
111  const unsigned int n_T_dofs = context.get_dof_indices(_T_var).size();
112 
113  // Element Jacobian * quadrature weights for interior integration.
114  const std::vector<libMesh::Real> &JxW =
115  context.get_element_fe(_u_r_var)->get_JxW();
116 
117  // The velocity shape functions at interior quadrature points.
118  const std::vector<std::vector<libMesh::Real> >& vel_phi =
119  context.get_element_fe(_u_r_var)->get_phi();
120 
121  // The temperature shape functions at interior quadrature points.
122  const std::vector<std::vector<libMesh::Real> >& T_phi =
123  context.get_element_fe(_T_var)->get_phi();
124 
125  // Physical location of the quadrature points
126  const std::vector<libMesh::Point>& u_qpoint =
127  context.get_element_fe(_u_r_var)->get_xyz();
128 
129  // Get residuals
130  libMesh::DenseSubVector<libMesh::Number> &Fr = context.get_elem_residual(_u_r_var); // R_{r}
131  libMesh::DenseSubVector<libMesh::Number> &Fz = context.get_elem_residual(_u_z_var); // R_{z}
132 
133  // Get Jacobians
134  libMesh::DenseSubMatrix<libMesh::Number> &KrT = context.get_elem_jacobian(_u_r_var, _T_var); // R_{r},{T}
135  libMesh::DenseSubMatrix<libMesh::Number> &KzT = context.get_elem_jacobian(_u_z_var, _T_var); // R_{z},{T}
136 
137  // Now we will build the element Jacobian and residual.
138  // Constructing the residual requires the solution and its
139  // gradient from the previous timestep. This must be
140  // calculated at each quadrature point by summing the
141  // solution degree-of-freedom values by the appropriate
142  // weight functions.
143  unsigned int n_qpoints = context.get_element_qrule().n_points();
144 
145  for (unsigned int qp=0; qp != n_qpoints; qp++)
146  {
147  const libMesh::Number r = u_qpoint[qp](0);
148 
149  // Compute the solution & its gradient at the old Newton iterate.
150  libMesh::Number T;
151  T = context.interior_value(_T_var, qp);
152 
153  // First, an i-loop over the velocity degrees of freedom.
154  // We know that n_u_dofs == n_v_dofs so we can compute contributions
155  // for both at the same time.
156  for (unsigned int i=0; i != n_u_dofs; i++)
157  {
158  Fr(i) += -_rho_ref*_beta_T*(T - _T_ref)*_g(0)*vel_phi[i][qp]*r*JxW[qp];
159  Fz(i) += -_rho_ref*_beta_T*(T - _T_ref)*_g(1)*vel_phi[i][qp]*r*JxW[qp];
160 
161  if (compute_jacobian && context.get_elem_solution_derivative())
162  {
163  for (unsigned int j=0; j != n_T_dofs; j++)
164  {
165  const libMesh::Number val =
166  -_rho_ref*_beta_T*vel_phi[i][qp]*T_phi[j][qp]*r*JxW[qp]
167  * context.get_elem_solution_derivative();
168  KrT(i,j) += val*_g(0);
169  KzT(i,j) += val*_g(1);
170  } // End j dof loop
171  } // End compute_jacobian check
172 
173  } // End i dof loop
174  } // End quadrature loop
175 
176 #ifdef GRINS_USE_GRVY_TIMERS
177  this->_timer->EndTimer("AxisymmetricBoussinesqBuoyancy::element_time_derivative");
178 #endif
179 
180  return;
181  }
libMesh::Number _rho_ref
reference density
VariableIndex _T_var
Index for temperature field.
VariableIndex _u_z_var
Index for z-velocity field.
VariableIndex _u_r_var
Index for r-velocity field.
libMesh::Point _g
Gravitational vector.
libMesh::Number _T_ref
reference temperature
libMesh::Number _beta_T
coefficient of thermal expansion
bool GRINS::Physics::enabled_on_elem ( const libMesh::Elem *  elem)
virtualinherited

Find if current physics is active on supplied element.

Definition at line 83 of file physics.C.

References GRINS::Physics::_enabled_subdomains.

84  {
85  // Check if enabled_subdomains flag has been set and if we're
86  // looking at a real element (rather than a nonlocal evaluation)
87  if( !elem || _enabled_subdomains.empty() )
88  return true;
89 
90  // Check if current physics is enabled on elem
91  if( _enabled_subdomains.find( elem->subdomain_id() ) == _enabled_subdomains.end() )
92  return false;
93 
94  return true;
95  }
std::set< libMesh::subdomain_id_type > _enabled_subdomains
Subdomains on which the current Physics class is enabled.
Definition: physics.h:261
BCHandlingBase * GRINS::Physics::get_bc_handler ( )
inlineinherited

Definition at line 282 of file physics.h.

References GRINS::Physics::_bc_handler.

283  {
284  return _bc_handler;
285  }
GRINS::BCHandlingBase * _bc_handler
Definition: physics.h:256
ICHandlingBase * GRINS::Physics::get_ic_handler ( )
inlineinherited

Definition at line 288 of file physics.h.

References GRINS::Physics::_ic_handler.

289  {
290  return _ic_handler;
291  }
GRINS::ICHandlingBase * _ic_handler
Definition: physics.h:258
void GRINS::Physics::init_bcs ( libMesh::FEMSystem *  system)
inherited

Definition at line 118 of file physics.C.

References GRINS::Physics::_bc_handler, GRINS::BCHandlingBase::init_bc_data(), GRINS::BCHandlingBase::init_dirichlet_bc_func_objs(), GRINS::BCHandlingBase::init_dirichlet_bcs(), and GRINS::BCHandlingBase::init_periodic_bcs().

119  {
120  // Only need to init BC's if the physics actually created a handler
121  if( _bc_handler )
122  {
123  _bc_handler->init_bc_data( *system );
124  _bc_handler->init_dirichlet_bcs( system );
126  _bc_handler->init_periodic_bcs( system );
127  }
128 
129  return;
130  }
GRINS::BCHandlingBase * _bc_handler
Definition: physics.h:256
virtual void init_dirichlet_bcs(libMesh::FEMSystem *system) const
virtual void init_periodic_bcs(libMesh::FEMSystem *system) const
virtual void init_dirichlet_bc_func_objs(libMesh::FEMSystem *system) const
virtual void init_bc_data(const libMesh::FEMSystem &system)
Override this method to initialize any system-dependent data.
void GRINS::Physics::init_context ( AssemblyContext context)
virtualinherited
void GRINS::Physics::init_ics ( libMesh::FEMSystem *  system,
libMesh::CompositeFunction< libMesh::Number > &  all_ics 
)
inherited

Definition at line 133 of file physics.C.

References GRINS::Physics::_ic_handler, and GRINS::ICHandlingBase::init_ic_data().

135  {
136  if( _ic_handler )
137  {
138  _ic_handler->init_ic_data( *system, all_ics );
139  }
140 
141  return;
142  }
GRINS::ICHandlingBase * _ic_handler
Definition: physics.h:258
virtual void init_ic_data(const libMesh::FEMSystem &system, libMesh::CompositeFunction< libMesh::Number > &all_ics)
Override this method to initialize any system-dependent data.
void GRINS::AxisymmetricBoussinesqBuoyancy::init_variables ( libMesh::FEMSystem *  system)
virtual

Initialization of AxisymmetricBoussinesqBuoyancy variables.

Implements GRINS::Physics.

Definition at line 91 of file axisym_boussinesq_buoyancy.C.

References _dim, _T_FE_family, _T_order, _T_var, _T_var_name, _u_r_var, _u_r_var_name, _u_z_var, _u_z_var_name, _V_FE_family, and _V_order.

92  {
93  this->_dim = system->get_mesh().mesh_dimension();
94 
95  _u_r_var = system->add_variable(_u_r_var_name, _V_order, _V_FE_family);
96  _u_z_var = system->add_variable(_u_z_var_name, _V_order, _V_FE_family);
97  _T_var = system->add_variable(_T_var_name, _T_order, _T_FE_family);
98  return;
99  }
GRINSEnums::Order _T_order
Temperature element order, read from input.
std::string _u_r_var_name
Name of r-velocity.
VariableIndex _T_var
Index for temperature field.
VariableIndex _u_z_var
Index for z-velocity field.
std::string _T_var_name
Name of temperature.
VariableIndex _u_r_var
Index for r-velocity field.
std::string _u_z_var_name
Name of z-velocity.
unsigned int _dim
Physical dimension of problem.
GRINSEnums::FEFamily _T_FE_family
Element type, read from input.
bool GRINS::Physics::is_steady ( ) const
inherited

Returns whether or not this physics is being solved with a steady solver.

Definition at line 103 of file physics.C.

References GRINS::Physics::_is_steady.

Referenced by GRINS::Physics::set_is_steady().

104  {
105  return _is_steady;
106  }
static bool _is_steady
Caches whether or not the solver that's being used is steady or not.
Definition: physics.h:266
void GRINS::Physics::mass_residual ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtualinherited
void GRINS::Physics::nonlocal_constraint ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtualinherited

Constraint part(s) of physics for scalar variables.

Reimplemented in GRINS::ScalarODE.

Definition at line 250 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::nonlocal_constraint().

253  {
254  return;
255  }
void GRINS::Physics::nonlocal_mass_residual ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtualinherited

Mass matrix part(s) for scalar variables.

Reimplemented in GRINS::ScalarODE, and GRINS::AveragedTurbine< Viscosity >.

Definition at line 264 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::nonlocal_mass_residual().

267  {
268  return;
269  }
void GRINS::Physics::nonlocal_time_derivative ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtualinherited

Time dependent part(s) of physics for scalar variables.

Reimplemented in GRINS::AveragedTurbine< Viscosity >, and GRINS::ScalarODE.

Definition at line 229 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::nonlocal_time_derivative().

232  {
233  return;
234  }
void GRINS::AxisymmetricBoussinesqBuoyancy::read_input_options ( const GetPot &  input)
virtual

Read options from GetPot input file.

Reimplemented from GRINS::Physics.

Definition at line 56 of file axisym_boussinesq_buoyancy.C.

References _beta_T, _g, _rho_ref, _T_FE_family, _T_order, _T_ref, _T_var_name, _u_r_var_name, _u_z_var_name, _V_FE_family, _V_order, GRINS::axisymmetric_boussinesq_buoyancy, GRINS::axisymmetric_heat_transfer, GRINS::incompressible_navier_stokes, GRINS::ParameterUser::set_parameter(), GRINS::T_var_name_default, GRINS::u_r_var_name_default, and GRINS::u_z_var_name_default.

Referenced by AxisymmetricBoussinesqBuoyancy().

57  {
58  this->_T_FE_family =
59  libMesh::Utility::string_to_enum<GRINSEnums::FEFamily>( input("Physics/"+axisymmetric_heat_transfer+"/FE_family", "LAGRANGE") );
60 
61  this->_T_order =
62  libMesh::Utility::string_to_enum<GRINSEnums::Order>( input("Physics/"+axisymmetric_heat_transfer+"/T_order", "SECOND") );
63 
64  this->_V_FE_family =
65  libMesh::Utility::string_to_enum<GRINSEnums::FEFamily>( input("Physics/"+incompressible_navier_stokes+"/FE_family", "LAGRANGE") );
66 
67  this->_V_order =
68  libMesh::Utility::string_to_enum<GRINSEnums::Order>( input("Physics/"+incompressible_navier_stokes+"/V_order", "SECOND") );
69 
70  // Read variable naming info
71  this->_u_r_var_name = input("Physics/VariableNames/r_velocity", u_r_var_name_default );
72  this->_u_z_var_name = input("Physics/VariableNames/z_velocity", u_z_var_name_default );
73  this->_T_var_name = input("Physics/VariableNames/Temperature", T_var_name_default );
74 
75  this->set_parameter
76  (_rho_ref, input,
77  "Physics/"+axisymmetric_boussinesq_buoyancy+"/rho_ref", 1.0);
78 
79  this->set_parameter
80  (_T_ref, input, "Physics/"+axisymmetric_boussinesq_buoyancy+"/T_ref", 1.0);
81 
82  this->set_parameter
83  (_beta_T, input, "Physics/"+axisymmetric_boussinesq_buoyancy+"/beta_T", 1.0);
84 
85  _g(0) = input("Physics/"+axisymmetric_boussinesq_buoyancy+"/g", 0.0, 0 );
86  _g(1) = input("Physics/"+axisymmetric_boussinesq_buoyancy+"/g", 0.0, 1 );
87 
88  return;
89  }
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.
const PhysicsName incompressible_navier_stokes
GRINSEnums::Order _T_order
Temperature element order, read from input.
const PhysicsName axisymmetric_heat_transfer
std::string _u_r_var_name
Name of r-velocity.
libMesh::Number _rho_ref
reference density
const std::string T_var_name_default
temperature
std::string _T_var_name
Name of temperature.
libMesh::Point _g
Gravitational vector.
const std::string u_z_var_name_default
z-velocity (axisymmetric case)
libMesh::Number _T_ref
reference temperature
const PhysicsName axisymmetric_boussinesq_buoyancy
std::string _u_z_var_name
Name of z-velocity.
libMesh::Number _beta_T
coefficient of thermal expansion
const std::string u_r_var_name_default
r-velocity (axisymmetric case)
GRINSEnums::FEFamily _T_FE_family
Element type, read from input.
void GRINS::ParameterUser::register_parameter ( const std::string &  param_name,
libMesh::ParameterMultiPointer< libMesh::Number > &  param_pointer 
) const
virtualinherited

Each subclass will register its copy of an independent.

Reimplemented in GRINS::AxisymmetricHeatTransfer< Conductivity >, GRINS::LowMachNavierStokesBase< Viscosity, SpecificHeat, ThermalConductivity >, GRINS::IncompressibleNavierStokesBase< Viscosity >, GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >, GRINS::HeatConduction< Conductivity >, GRINS::HeatTransferBase< Conductivity >, and GRINS::BoussinesqBuoyancyAdjointStabilization< Viscosity >.

Definition at line 50 of file parameter_user.C.

Referenced by GRINS::BoussinesqBuoyancyAdjointStabilization< Viscosity >::register_parameter(), GRINS::HeatTransferBase< Conductivity >::register_parameter(), GRINS::HeatConduction< Conductivity >::register_parameter(), GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::register_parameter(), GRINS::IncompressibleNavierStokesBase< Viscosity >::register_parameter(), GRINS::LowMachNavierStokesBase< Viscosity, SpecificHeat, ThermalConductivity >::register_parameter(), and GRINS::AxisymmetricHeatTransfer< Conductivity >::register_parameter().

53  {
54  std::map<std::string, libMesh::Number*>::const_iterator it =
55  _my_parameters.find(param_name);
56 
57  if (it != _my_parameters.end())
58  {
59  std::cout << _my_name << " uses parameter " << param_name
60  << std::endl;
61  param_pointer.push_back(it->second);
62  }
63  }
std::map< std::string, libMesh::Number * > _my_parameters
void GRINS::Physics::register_postprocessing_vars ( const GetPot &  input,
PostProcessedQuantities< libMesh::Real > &  postprocessing 
)
virtualinherited

Register name of postprocessed quantity with PostProcessedQuantities.

Each Physics class will need to cache an unsigned int corresponding to each postprocessed quantity. This will be used in computing the values and putting them in the CachedVariables object.

Reimplemented in GRINS::ParsedVelocitySource< Viscosity >, GRINS::VelocityPenalty< Viscosity >, GRINS::IncompressibleNavierStokes< Viscosity >, GRINS::LowMachNavierStokes< Viscosity, SpecificHeat, ThermalConductivity >, GRINS::HeatTransfer< Conductivity >, GRINS::ElasticCable< StressStrainLaw >, GRINS::ElasticMembrane< StressStrainLaw >, and GRINS::ReactingLowMachNavierStokes< Mixture, Evaluator >.

Definition at line 161 of file physics.C.

163  {
164  return;
165  }
void GRINS::Physics::set_is_steady ( bool  is_steady)
inherited

Sets whether this physics is to be solved with a steady solver or not.

Since the member variable is static, only needs to be called on a single physics.

Definition at line 97 of file physics.C.

References GRINS::Physics::_is_steady, and GRINS::Physics::is_steady().

98  {
100  return;
101  }
bool is_steady() const
Returns whether or not this physics is being solved with a steady solver.
Definition: physics.C:103
static bool _is_steady
Caches whether or not the solver that's being used is steady or not.
Definition: physics.h:266
void GRINS::ParameterUser::set_parameter ( libMesh::Number &  param_variable,
const GetPot &  input,
const std::string &  param_name,
libMesh::Number  param_default 
)
virtualinherited

Each subclass can simultaneously read a parameter value from.

Definition at line 35 of file parameter_user.C.

References GRINS::ParameterUser::_my_name, and GRINS::ParameterUser::_my_parameters.

Referenced by GRINS::AveragedFanAdjointStabilization< Viscosity >::AveragedFanAdjointStabilization(), GRINS::AveragedTurbineAdjointStabilization< Viscosity >::AveragedTurbineAdjointStabilization(), GRINS::BoussinesqBuoyancyAdjointStabilization< Viscosity >::BoussinesqBuoyancyAdjointStabilization(), GRINS::BoussinesqBuoyancyBase::BoussinesqBuoyancyBase(), GRINS::BoussinesqBuoyancySPGSMStabilization< Viscosity >::BoussinesqBuoyancySPGSMStabilization(), GRINS::ConstantConductivity::ConstantConductivity(), GRINS::ConstantPrandtlConductivity::ConstantPrandtlConductivity(), GRINS::ConstantSourceFunction::ConstantSourceFunction(), GRINS::ConstantSourceTerm::ConstantSourceTerm(), GRINS::ConstantSpecificHeat::ConstantSpecificHeat(), GRINS::ConstantViscosity::ConstantViscosity(), GRINS::ElasticCable< StressStrainLaw >::ElasticCable(), GRINS::ElasticCableConstantGravity::ElasticCableConstantGravity(), GRINS::ElasticMembrane< StressStrainLaw >::ElasticMembrane(), GRINS::ElasticMembraneConstantPressure::ElasticMembraneConstantPressure(), GRINS::HeatConduction< Conductivity >::HeatConduction(), GRINS::HeatTransferBase< Conductivity >::HeatTransferBase(), GRINS::IncompressibleNavierStokesBase< Viscosity >::IncompressibleNavierStokesBase(), GRINS::AverageNusseltNumber::init(), GRINS::MooneyRivlin::MooneyRivlin(), GRINS::ReactingLowMachNavierStokesBase< Mixture, Evaluator >::ReactingLowMachNavierStokesBase(), GRINS::HookesLaw1D::read_input_options(), GRINS::HookesLaw::read_input_options(), read_input_options(), and GRINS::VelocityDragAdjointStabilization< Viscosity >::VelocityDragAdjointStabilization().

39  {
40  param_variable = input(param_name, param_default);
41 
42  libmesh_assert_msg(!_my_parameters.count(param_name),
43  "ERROR: " << _my_name << " double-registered parameter " <<
44  param_name);
45 
46  _my_parameters[param_name] = &param_variable;
47  }
std::map< std::string, libMesh::Number * > _my_parameters
void GRINS::Physics::set_time_evolving_vars ( libMesh::FEMSystem *  system)
virtualinherited

Set which variables are time evolving.

Set those variables which evolve in time (as opposed to variables that behave like constraints). This is done separately from init_variables() because the MultiphysicsSystem must initialize its base class before time_evolving variables can be set. Default implementation is no time evolving variables.

Reimplemented in GRINS::SpalartAllmaras< Viscosity >, GRINS::ScalarODE, GRINS::AxisymmetricHeatTransfer< Conductivity >, GRINS::IncompressibleNavierStokesBase< Viscosity >, GRINS::AveragedTurbineBase< Viscosity >, GRINS::LowMachNavierStokesBase< Viscosity, SpecificHeat, ThermalConductivity >, GRINS::HeatTransferBase< Conductivity >, GRINS::ParsedVelocitySourceBase< Viscosity >, GRINS::ReactingLowMachNavierStokesBase< Mixture, Evaluator >, GRINS::ElasticCableBase, GRINS::ElasticMembraneBase, and GRINS::HeatConduction< Conductivity >.

Definition at line 108 of file physics.C.

109  {
110  return;
111  }
void GRINS::Physics::side_constraint ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtualinherited

Constraint part(s) of physics for boundaries of elements on the domain boundary.

Reimplemented in GRINS::IncompressibleNavierStokes< Viscosity >, GRINS::LowMachNavierStokes< Viscosity, SpecificHeat, ThermalConductivity >, and GRINS::ReactingLowMachNavierStokes< Mixture, Evaluator >.

Definition at line 243 of file physics.C.

Referenced by GRINS::MultiphysicsSystem::side_constraint().

246  {
247  return;
248  }
void GRINS::Physics::side_time_derivative ( bool  compute_jacobian,
AssemblyContext context,
CachedValues cache 
)
virtualinherited

Member Data Documentation

GRINS::BCHandlingBase* GRINS::Physics::_bc_handler
protectedinherited
libMesh::Number GRINS::AxisymmetricBoussinesqBuoyancy::_beta_T
protected

$ \beta_T = $ coefficient of thermal expansion

Definition at line 116 of file axisym_boussinesq_buoyancy.h.

Referenced by element_time_derivative(), and read_input_options().

unsigned int GRINS::AxisymmetricBoussinesqBuoyancy::_dim
protected

Physical dimension of problem.

Definition at line 80 of file axisym_boussinesq_buoyancy.h.

Referenced by init_variables().

std::set<libMesh::subdomain_id_type> GRINS::Physics::_enabled_subdomains
protectedinherited

Subdomains on which the current Physics class is enabled.

Definition at line 261 of file physics.h.

Referenced by GRINS::Physics::enabled_on_elem(), and GRINS::Physics::read_input_options().

libMesh::Point GRINS::AxisymmetricBoussinesqBuoyancy::_g
protected

Gravitational vector.

Definition at line 120 of file axisym_boussinesq_buoyancy.h.

Referenced by element_time_derivative(), and read_input_options().

GRINS::ICHandlingBase* GRINS::Physics::_ic_handler
protectedinherited
bool GRINS::Physics::_is_axisymmetric
protectedinherited
bool GRINS::Physics::_is_steady = false
staticprotectedinherited

Caches whether or not the solver that's being used is steady or not.

This is need, for example, in flow stabilization as the tau terms change depending on whether the solver is steady or unsteady.

Definition at line 266 of file physics.h.

Referenced by GRINS::Physics::is_steady(), and GRINS::Physics::set_is_steady().

const PhysicsName GRINS::Physics::_physics_name
protectedinherited

Name of the physics object. Used for reading physics specific inputs.

We use a reference because the physics names are const global objects in GRINS namespace

No, we use a copy, because otherwise as soon as the memory in std::set<std::string> requested_physics gets overwritten we get in trouble.

Definition at line 254 of file physics.h.

Referenced by GRINS::SourceTermBase::parse_var_info(), and GRINS::Physics::read_input_options().

libMesh::Number GRINS::AxisymmetricBoussinesqBuoyancy::_rho_ref
protected

$ \rho_0 = $ reference density

Definition at line 110 of file axisym_boussinesq_buoyancy.h.

Referenced by element_time_derivative(), and read_input_options().

GRINSEnums::FEFamily GRINS::AxisymmetricBoussinesqBuoyancy::_T_FE_family
protected

Element type, read from input.

Definition at line 83 of file axisym_boussinesq_buoyancy.h.

Referenced by init_variables(), and read_input_options().

GRINSEnums::Order GRINS::AxisymmetricBoussinesqBuoyancy::_T_order
protected

Temperature element order, read from input.

Definition at line 86 of file axisym_boussinesq_buoyancy.h.

Referenced by init_variables(), and read_input_options().

libMesh::Number GRINS::AxisymmetricBoussinesqBuoyancy::_T_ref
protected

$ T_0 = $ reference temperature

Definition at line 113 of file axisym_boussinesq_buoyancy.h.

Referenced by element_time_derivative(), and read_input_options().

VariableIndex GRINS::AxisymmetricBoussinesqBuoyancy::_T_var
protected

Index for temperature field.

Definition at line 96 of file axisym_boussinesq_buoyancy.h.

Referenced by element_time_derivative(), and init_variables().

std::string GRINS::AxisymmetricBoussinesqBuoyancy::_T_var_name
protected

Name of temperature.

Definition at line 107 of file axisym_boussinesq_buoyancy.h.

Referenced by init_variables(), and read_input_options().

VariableIndex GRINS::AxisymmetricBoussinesqBuoyancy::_u_r_var
protected

Index for r-velocity field.

Definition at line 90 of file axisym_boussinesq_buoyancy.h.

Referenced by element_time_derivative(), and init_variables().

std::string GRINS::AxisymmetricBoussinesqBuoyancy::_u_r_var_name
protected

Name of r-velocity.

Definition at line 101 of file axisym_boussinesq_buoyancy.h.

Referenced by init_variables(), and read_input_options().

VariableIndex GRINS::AxisymmetricBoussinesqBuoyancy::_u_z_var
protected

Index for z-velocity field.

Definition at line 93 of file axisym_boussinesq_buoyancy.h.

Referenced by element_time_derivative(), and init_variables().

std::string GRINS::AxisymmetricBoussinesqBuoyancy::_u_z_var_name
protected

Name of z-velocity.

Definition at line 104 of file axisym_boussinesq_buoyancy.h.

Referenced by init_variables(), and read_input_options().

GRINSEnums::FEFamily GRINS::AxisymmetricBoussinesqBuoyancy::_V_FE_family
protected

Definition at line 83 of file axisym_boussinesq_buoyancy.h.

Referenced by init_variables(), and read_input_options().

GRINSEnums::Order GRINS::AxisymmetricBoussinesqBuoyancy::_V_order
protected

Definition at line 86 of file axisym_boussinesq_buoyancy.h.

Referenced by init_variables(), and read_input_options().


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

Generated on Mon Jun 22 2015 21:32:22 for GRINS-0.6.0 by  doxygen 1.8.9.1