34 #include "libmesh/utility.h"
35 #include "libmesh/string_to_enum.h"
36 #include "libmesh/getpot.h"
37 #include "libmesh/fem_system.h"
38 #include "libmesh/quadrature.h"
59 libMesh::Utility::string_to_enum<GRINSEnums::FEFamily>( input(
"Physics/"+
axisymmetric_heat_transfer+
"/FE_family",
"LAGRANGE") );
93 this->
_dim = system->get_mesh().mesh_dimension();
105 #ifdef GRINS_USE_GRVY_TIMERS
106 this->_timer->BeginTimer(
"AxisymmetricBoussinesqBuoyancy::element_time_derivative");
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();
114 const std::vector<libMesh::Real> &JxW =
115 context.get_element_fe(
_u_r_var)->get_JxW();
118 const std::vector<std::vector<libMesh::Real> >& vel_phi =
119 context.get_element_fe(
_u_r_var)->get_phi();
122 const std::vector<std::vector<libMesh::Real> >& T_phi =
123 context.get_element_fe(
_T_var)->get_phi();
126 const std::vector<libMesh::Point>& u_qpoint =
127 context.get_element_fe(
_u_r_var)->get_xyz();
130 libMesh::DenseSubVector<libMesh::Number> &Fr = context.get_elem_residual(
_u_r_var);
131 libMesh::DenseSubVector<libMesh::Number> &Fz = context.get_elem_residual(
_u_z_var);
134 libMesh::DenseSubMatrix<libMesh::Number> &KrT = context.get_elem_jacobian(
_u_r_var,
_T_var);
135 libMesh::DenseSubMatrix<libMesh::Number> &KzT = context.get_elem_jacobian(
_u_z_var,
_T_var);
143 unsigned int n_qpoints = context.get_element_qrule().n_points();
145 for (
unsigned int qp=0; qp != n_qpoints; qp++)
147 const libMesh::Number r = u_qpoint[qp](0);
151 T = context.interior_value(
_T_var, qp);
156 for (
unsigned int i=0; i != n_u_dofs; i++)
161 if (compute_jacobian && context.get_elem_solution_derivative())
163 for (
unsigned int j=0; j != n_T_dofs; j++)
165 const libMesh::Number val =
167 * context.get_elem_solution_derivative();
168 KrT(i,j) += val*
_g(0);
169 KzT(i,j) += val*
_g(1);
176 #ifdef GRINS_USE_GRVY_TIMERS
177 this->_timer->EndTimer(
"AxisymmetricBoussinesqBuoyancy::element_time_derivative");
virtual void set_parameter(libMesh::Number ¶m_variable, const GetPot &input, const std::string ¶m_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.
Physics abstract base class. Defines API for physics to be added to MultiphysicsSystem.
const PhysicsName axisymmetric_heat_transfer
GRINSEnums::FEFamily _V_FE_family
std::string _u_r_var_name
Name of r-velocity.
libMesh::Number _rho_ref
reference density
virtual void read_input_options(const GetPot &input)
Read options from GetPot input file.
GRINSEnums::Order _V_order
~AxisymmetricBoussinesqBuoyancy()
virtual void element_time_derivative(bool compute_jacobian, AssemblyContext &context, CachedValues &cache)
Source term contribution for AxisymmetricBoussinesqBuoyancy.
VariableIndex _T_var
Index for temperature field.
const std::string T_var_name_default
temperature
VariableIndex _u_z_var
Index for z-velocity field.
virtual void init_variables(libMesh::FEMSystem *system)
Initialization of AxisymmetricBoussinesqBuoyancy variables.
std::string _T_var_name
Name of temperature.
VariableIndex _u_r_var
Index for r-velocity field.
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
AxisymmetricBoussinesqBuoyancy()
std::string _u_z_var_name
Name of z-velocity.
libMesh::Number _beta_T
coefficient of thermal expansion
unsigned int _dim
Physical dimension of problem.
const std::string u_r_var_name_default
r-velocity (axisymmetric case)
GRINSEnums::FEFamily _T_FE_family
Element type, read from input.