43 #ifndef __Panzer_Integrator_DivBasisTimesScalar_impl_hpp__ 44 #define __Panzer_Integrator_DivBasisTimesScalar_impl_hpp__ 46 #include "Intrepid2_FunctionSpaceTools.hpp" 51 #include "Kokkos_ViewFactory.hpp" 56 template<
typename EvalT,
typename Traits>
59 const Teuchos::ParameterList& p) :
60 residual( p.get<
std::string>(
"Residual Name"),
65 p.validateParameters(*valid_params);
67 Teuchos::RCP<const PureBasis> basis
68 = p.get< Teuchos::RCP<BasisIRLayout> >(
"Basis")->getBasis();
71 TEUCHOS_TEST_FOR_EXCEPTION(!basis->supportsDiv(),std::logic_error,
72 "Integrator_DivBasisTimesScalar: Basis of type \"" << basis->name() <<
"\" does not support DIV.");
73 TEUCHOS_TEST_FOR_EXCEPTION(!basis->requiresOrientations(),std::logic_error,
74 "Integration_DivBasisTimesScalar: Basis of type \"" << basis->name() <<
"\" should require orientations. So we are throwing.");
76 scalar = PHX::MDField<const ScalarT,Cell,IP>( p.get<std::string>(
"Value Name"),
77 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_scalar );
80 this->addDependentField(
scalar);
83 if (p.isType<Teuchos::RCP<
const std::vector<std::string> > >(
"Field Multipliers")) {
85 const std::vector<std::string>& field_multiplier_names =
86 *(p.get<Teuchos::RCP<const std::vector<std::string> > >(
"Field Multipliers"));
88 for (
const std::string & name : field_multiplier_names) {
89 PHX::MDField<const ScalarT,Cell,IP> tmp_field(name, p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_scalar);
95 this->addDependentField(
field);
98 "Integrator_DivBasisTimesScalar: " +
residual.fieldTag().name();
104 template<
typename EvalT,
typename Traits>
111 this->utils.setFieldData(residual,fm);
112 this->utils.setFieldData(scalar,fm);
114 for (
auto &
field : field_multipliers)
115 this->utils.setFieldData(
field,fm);
117 num_nodes = residual.extent(1);
118 num_qp = scalar.extent(1);
126 template<
typename EvalT,
typename Traits>
133 residual.deep_copy(
ScalarT(0.0));
135 for (index_t cell = 0; cell < workset.
num_cells; ++cell) {
136 for (std::size_t qp = 0; qp < num_qp; ++qp) {
138 for (
const auto &
field : field_multipliers)
139 tmpVar = tmpVar *
field(cell,qp);
142 tmp(cell,qp) =
multiplier * tmpVar * scalar(cell,qp);
150 for (index_t cell = 0; cell < workset.
num_cells; ++cell)
151 for (std::size_t basis = 0; basis < num_nodes; ++basis) {
152 for (std::size_t qp = 0; qp < num_qp; ++qp)
168 template<
typename EvalT,
typename TRAITS>
169 Teuchos::RCP<Teuchos::ParameterList>
172 Teuchos::RCP<Teuchos::ParameterList> p = Teuchos::rcp(
new Teuchos::ParameterList);
173 p->set<std::string>(
"Residual Name",
"?");
174 p->set<std::string>(
"Value Name",
"?");
175 p->set<std::string>(
"Test Field Name",
"?");
176 Teuchos::RCP<panzer::BasisIRLayout> basis;
177 p->set(
"Basis", basis);
178 Teuchos::RCP<panzer::IntegrationRule> ir;
180 p->set<
double>(
"Multiplier", 1.0);
181 Teuchos::RCP<const std::vector<std::string> > fms;
182 p->set(
"Field Multipliers", fms);
void evaluateFields(typename Traits::EvalData d)
typename EvalT::ScalarT ScalarT
Integrator_DivBasisTimesScalar(const Teuchos::ParameterList &p)
double multiplier
The scalar multiplier out in front of the integral ( ).
Kokkos::DynRankView< typename InputArray::value_type, PHX::Device > createDynRankView(const InputArray &a, const std::string &name, const DimensionPack... dims)
Wrapper to simplify Panzer use of Sacado ViewFactory.
std::vector< PHX::MDField< const ScalarT, Cell, IP > > field_multipliers
PHX::MDField< const ScalarT, Cell, IP > scalar
std::vector< std::string >::size_type getBasisIndex(std::string basis_name, const panzer::Workset &workset, WorksetDetailsAccessor &wda)
Returns the index in the workset bases for a particular BasisIRLayout name.
Teuchos::RCP< Teuchos::ParameterList > getValidParameters() const
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
PHX::MDField< ScalarT, panzer::Cell, panzer::BASIS > field
A field to which we'll contribute, or in which we'll store, the result of computing this integral...
Array_CellBasisIP weighted_div_basis
PHX::MDField< ScalarT, Cell, BASIS > residual
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_