35 #include "libmesh/getpot.h"
36 #include "libmesh/fem_system.h"
37 #include "libmesh/quadrature.h"
65 this->parse_thermal_conductivity(input);
68 ( _scaling, input,
"QoI/NusseltNumber/scaling", 1.0 );
72 std::cerr <<
"Error: thermal conductivity for AverageNusseltNumber must be positive." << std::endl
73 <<
"Found k = " << _k << std::endl;
78 int num_bcs = input.vector_variable_size(
"QoI/NusseltNumber/bc_ids");
82 std::cerr <<
"Error: Must specify at least one boundary id to compute"
83 <<
" average Nusselt number." << std::endl
84 <<
"Found: " << num_bcs << std::endl;
88 for(
int i = 0; i < num_bcs; i++ )
90 _bc_ids.insert( input(
"QoI/NusseltNumber/bc_ids", -1, i ) );
94 std::string T_var_name = input(
"Physics/VariableNames/Temperature",
97 this->_T_var = system.variable_number(T_var_name);
104 libMesh::FEBase* T_fe;
106 context.get_side_fe<libMesh::Real>(this->
_T_var, T_fe);
115 const unsigned int qoi_index )
117 bool on_correct_side =
false;
119 for (std::set<libMesh::boundary_id_type>::const_iterator
id =
121 if( context.has_side_boundary_id( (*
id) ) )
123 on_correct_side =
true;
127 if (!on_correct_side)
130 libMesh::FEBase* side_fe;
131 context.get_side_fe<libMesh::Real>(this->
_T_var, side_fe);
133 const std::vector<libMesh::Real> &JxW = side_fe->get_JxW();
135 const std::vector<libMesh::Point>& normals = side_fe->get_normals();
137 unsigned int n_qpoints = context.get_side_qrule().n_points();
139 libMesh::Number& qoi = context.get_qois()[qoi_index];
143 for (
unsigned int qp = 0; qp != n_qpoints; qp++)
146 libMesh::Gradient grad_T = 0.0;
147 context.side_gradient(this->_T_var, qp, grad_T);
150 qoi += (this->
_scaling)*(this->
_k)*(grad_T*normals[qp])*JxW[qp];
156 const unsigned int qoi_index )
159 for( std::set<libMesh::boundary_id_type>::const_iterator
id =
_bc_ids.begin();
162 if( context.has_side_boundary_id( (*
id) ) )
164 libMesh::FEBase* T_side_fe;
165 context.get_side_fe<libMesh::Real>(this->
_T_var, T_side_fe);
167 const std::vector<libMesh::Real> &JxW = T_side_fe->get_JxW();
169 const std::vector<libMesh::Point>& normals = T_side_fe->get_normals();
171 unsigned int n_qpoints = context.get_side_qrule().n_points();
173 const unsigned int n_T_dofs = context.get_dof_indices(_T_var).size();
175 const std::vector<std::vector<libMesh::Gradient> >& T_gradphi = T_side_fe->get_dphi();
177 libMesh::DenseSubVector<libMesh::Number>& dQ_dT =
178 context.get_qoi_derivatives(qoi_index, _T_var);
181 for (
unsigned int qp = 0; qp != n_qpoints; qp++)
184 libMesh::Gradient grad_T = 0.0;
185 context.side_gradient(this->_T_var, qp, grad_T);
190 for(
unsigned int i = 0; i != n_T_dofs; i++ )
192 dQ_dT(i) +=
_scaling*
_k*T_gradphi[i][qp]*normals[qp]*JxW[qp];
206 std::string material = input(
"QoI/NusseltNumber/material",
"NoMaterial!");
209 "Materials/"+material+
"/ThermalConductivity/value",
210 "QoI/NusseltNumber/thermal_conductivity");
213 if( input.have_variable(
"QoI/NusseltNumber/thermal_conductivity") )
216 "ThermalConductivity/value" );
219 (
_k, input,
"QoI/NusseltNumber/thermal_conductivity", -1.0 );
222 else if( input.have_variable(
"Materials/"+material+
"/ThermalConductivity/value") )
225 if( input(
"Materials/"+material+
"/ThermalConductivity/model",
"DIE!")
226 != std::string(
"constant") )
228 libmesh_error_msg(
"ERROR: Only constant ThermalConductivity model supported in NusseltNumber!");
232 (
_k, input,
"Materials/"+material+
"/ThermalConductivity/value", -1.0 );
235 libmesh_error_msg(
"ERROR: Could not find valid thermal conducitivity value!");
239 libmesh_error_msg(
"ERROR: Detected non-positive thermal conductivity!");
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.
std::set< libMesh::boundary_id_type > _bc_ids
List of boundary ids for which we want to compute this QoI.
virtual void init_context(AssemblyContext &context)
virtual ~AverageNusseltNumber()
libMesh::Real _scaling
Scaling constant.
virtual void init(const GetPot &input, const MultiphysicsSystem &system, unsigned int qoi_num)
Method to allow QoI to cache any system information needed for QoI calculation, for example...
const std::string T_var_name_default
temperature
virtual QoIBase * clone() const
Clone this QoI.
virtual void side_qoi_derivative(AssemblyContext &context, const unsigned int qoi_index)
Compute the qoi derivative with respect to the solution on the domain boundary.
libMesh::Real _k
Thermal conductivity.
virtual void side_qoi(AssemblyContext &context, const unsigned int qoi_index)
Compute the qoi value on the domain boundary.
VariableIndex _T_var
Temperature variable index.
Interface with libMesh for solving Multiphysics problems.
static void dep_input_warning(const std::string &old_option, const std::string &property)
Helper function for parsing/maintaing backward compatibility.
static void duplicate_input_test(const GetPot &input, const std::string &option1, const std::string &option2)
Helper function for parsing/maintaing backward compatibility.
virtual void move_parameter(const libMesh::Number &old_parameter, libMesh::Number &new_parameter)
When cloning an object, we need to update parameter pointers.
void parse_thermal_conductivity(const GetPot &input)