43 #ifndef PANZER_WEAKDIRICHLET_RESIDUAL_IMPL_HPP 44 #define PANZER_WEAKDIRICHLET_RESIDUAL_IMPL_HPP 51 #include "Intrepid2_FunctionSpaceTools.hpp" 59 std::string residual_name = p.get<std::string>(
"Residual Name");
60 std::string flux_name = p.get<std::string>(
"Flux Name");
61 std::string normal_name = p.get<std::string>(
"Normal Name");
62 std::string normal_dot_flux_name = normal_name +
" dot " + flux_name;
63 std::string dof_name = p.get<std::string>(
"DOF Name");
64 std::string value_name = p.get<std::string>(
"Value Name");
65 std::string sigma_name = p.get<std::string>(
"Sigma Name");
76 flux = PHX::MDField<ScalarT>(flux_name, ir->dl_vector);
77 normal = PHX::MDField<ScalarT>(normal_name, ir->dl_vector);
78 dof = PHX::MDField<ScalarT>(dof_name, ir->dl_scalar);
79 value = PHX::MDField<ScalarT>(value_name, ir->dl_scalar);
80 sigma = PHX::MDField<ScalarT>(sigma_name, ir->dl_scalar);
84 this->addDependentField(
normal);
85 this->addDependentField(
flux);
86 this->addDependentField(
dof);
87 this->addDependentField(
value);
88 this->addDependentField(
sigma);
92 std::string n =
"Weak Dirichlet Residual Evaluator";
99 this->utils.setFieldData(
residual,fm);
101 this->utils.setFieldData(
flux,fm);
102 this->utils.setFieldData(
normal,fm);
103 this->utils.setFieldData(
dof,fm);
104 this->utils.setFieldData(
value,fm);
105 this->utils.setFieldData(
sigma,fm);
119 for (index_t cell = 0; cell < workset.num_cells; ++cell) {
120 for (std::size_t ip = 0; ip <
num_ip; ++ip) {
122 for (std::size_t dim = 0; dim <
num_dim; ++dim) {
129 if(workset.num_cells>0)
130 Intrepid2::FunctionSpaceTools::
132 (this->wda(workset).bases[
basis_index])->weighted_basis_scalar,
133 Intrepid2::COMP_CPP);
PHX::MDField< ScalarT > residual
Evaluates a Dirichlet BC residual corresponding to a field value.
std::vector< std::string >::size_type getBasisIndex(std::string basis_name, panzer::Workset &workset, WorksetDetailsAccessor &wda)
Returns the index in the workset bases for a particular BasisIRLayout name.
PHX::MDField< ScalarT > sigma
PHX::MDField< ScalarT > flux
PHX::MDField< ScalarT > normal
Teuchos::RCP< panzer::BasisIRLayout > basisIRLayout(std::string basis_type, const int basis_order, const PointRule &pt_rule)
Nonmember constructor.
PHX::MDField< const ScalarT > dof
PHX_EVALUATOR_CTOR(BasisValues_Evaluator, p)
PHX_EVALUATE_FIELDS(BasisValues_Evaluator, workset)
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
#define TEUCHOS_ASSERT(assertion_test)
Teuchos::RCP< PHX::DataLayout > functional
<Cell,Basis> or <Cell,Basis>
PHX_POST_REGISTRATION_SETUP(BasisValues_Evaluator, sd, fm)
PHX::MDField< ScalarT > normal_dot_flux_plus_pen