43 #ifndef __Panzer_Integrator_DivBasisTimesScalar_impl_hpp__
44 #define __Panzer_Integrator_DivBasisTimesScalar_impl_hpp__
53 #include "Intrepid2_FunctionSpaceTools.hpp"
56 #include "Kokkos_ViewFactory.hpp"
71 template<
typename EvalT,
typename Traits>
75 const std::string& resName,
76 const std::string& valName,
80 const std::vector<std::string>& fmNames
83 evalStyle_(evalStyle),
84 multiplier_(multiplier),
85 basisName_(basis.name()),
86 use_shared_memory(false)
95 using std::invalid_argument;
96 using std::logic_error;
102 "Integrator_DivBasisTimesScalar called with an empty residual name.")
104 "Integrator_DivBasisTimesScalar called with an empty value name.")
107 "Error: Integrator_DivBasisTimesScalar: Basis of type \""
108 << tmpBasis->name() <<
"\" does not support DIV.")
110 "Error: Integration_DivBasisTimesScalar: Basis of type \""
111 << tmpBasis->name() <<
"\" should require orientations.")
115 this->addDependentField(
scalar_);
121 this->addContributedField(
field_);
123 this->addEvaluatedField(
field_);
129 View<PHX::UnmanagedView<const ScalarT**>*>(
"DivBasisTimesScalar::KokkosFieldMultipliers",
131 for (
const auto& name : fmNames)
138 string n(
"Integrator_DivBasisTimesScalar (");
143 n +=
"): " +
field_.fieldTag().name();
152 template<
typename EvalT,
typename Traits>
159 p.get<std::string>(
"Residual Name"),
160 p.get<std::string>(
"Value Name"),
163 p.get<double>(
"Multiplier"),
164 p.isType<Teuchos::
RCP<
const std::vector<std::string>>>
165 (
"Field Multipliers") ?
166 (*p.get<Teuchos::
RCP<
const std::vector<std::string>>>
167 (
"Field Multipliers")) : std::vector<std::string>())
182 template<
typename EvalT,
typename Traits>
192 auto kokkosFieldMults_h = Kokkos::create_mirror_view(kokkosFieldMults_);
195 for (
size_t i(0); i < fieldMults_.size(); ++i)
196 kokkosFieldMults_h(i) = fieldMults_[i].get_static_view();
198 Kokkos::deep_copy(kokkosFieldMults_, kokkosFieldMults_h);
202 if (!use_shared_memory) {
203 if (Sacado::IsADType<ScalarT>::value) {
204 const auto fadSize = Kokkos::dimension_scalar(field_.get_view());
205 tmp_ = PHX::View<ScalarT*>(
"panzer::Integrator::DivBasisTimesScalar::tmp_",field_.extent(0),fadSize);
207 tmp_ = PHX::View<ScalarT*>(
"panzer::Integrator::DivBasisTimesScalar::tmp_",field_.extent(0));
220 template<
typename EvalT,
typename Traits>
221 template<
int NUM_FIELD_MULT>
222 KOKKOS_INLINE_FUNCTION
227 const Kokkos::TeamPolicy<PHX::exec_space>::member_type& team)
const
230 const int cell = team.league_rank();
233 const int numQP(scalar_.extent(1)), numBases(basis_.extent(1));
235 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases), [&] (
const int basis) {
236 field_(cell, basis) = 0.0;
242 if (NUM_FIELD_MULT == 0)
247 for (
int qp(0); qp < numQP; ++qp)
249 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases), [&] (
const int basis) {
250 field_(cell, basis) += basis_(cell, basis, qp) * multiplier_ * scalar_(cell, qp);
254 else if (NUM_FIELD_MULT == 1)
259 for (
int qp(0); qp < numQP; ++qp)
261 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases), [&] (
const int basis) {
262 field_(cell, basis) += basis_(cell, basis, qp) * multiplier_ * scalar_(cell, qp) * kokkosFieldMults_(0)(cell, qp);
272 const int numFieldMults(kokkosFieldMults_.extent(0));
273 for (
int qp(0); qp < numQP; ++qp)
277 for (
int fm(0); fm < numFieldMults; ++fm)
278 tmp_(cell) *= kokkosFieldMults_(fm)(cell, qp);
279 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases),[&] (
const int& basis) {
280 field_(cell, basis) += basis_(cell, basis, qp) * multiplier_ * scalar_(cell, qp) * tmp_(cell);
291 template<
typename EvalT,
typename Traits>
292 template<
int NUM_FIELD_MULT>
293 KOKKOS_INLINE_FUNCTION
298 const Kokkos::TeamPolicy<PHX::exec_space>::member_type& team)
const
301 const int cell = team.league_rank();
302 const int numQP = scalar_.extent(1);
303 const int numBases = basis_.extent(1);
304 const int fadSize = Kokkos::dimension_scalar(field_.get_view());
308 if (Sacado::IsADType<ScalarT>::value) {
310 tmp_field =
scratch_view(team.team_shmem(),numBases,fadSize);
318 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases),[&] (
const int& basis) {
319 tmp_field(basis) = 0.0;
324 if (NUM_FIELD_MULT == 0)
329 for (
int qp(0); qp < numQP; ++qp)
332 tmp(0) = multiplier_ * scalar_(cell, qp);
333 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases),[&] (
const int& basis) {
334 tmp_field(basis) += basis_(cell, basis, qp) * tmp(0);
338 else if (NUM_FIELD_MULT == 1)
343 for (
int qp(0); qp < numQP; ++qp)
346 tmp(0) = multiplier_ * scalar_(cell, qp) * kokkosFieldMults_(0)(cell, qp);
347 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases),[&] (
const int& basis) {
348 tmp_field(basis) += basis_(cell, basis, qp) * tmp(0);
358 const int numFieldMults = kokkosFieldMults_.extent(0);
359 for (
int qp(0); qp < numQP; ++qp)
363 for (
int fm(0); fm < numFieldMults; ++fm)
364 fieldMultsTotal *= kokkosFieldMults_(fm)(cell, qp);
365 tmp(0) = multiplier_ * scalar_(cell, qp) * fieldMultsTotal;
366 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases),[&] (
const int& basis) {
367 tmp_field(basis) += basis_(cell, basis, qp) * tmp(0);
374 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases),[&] (
const int& basis) {
375 field_(cell,basis) = tmp_field(basis);
379 Kokkos::parallel_for(Kokkos::TeamThreadRange(team,0,numBases),[&] (
const int& basis) {
380 field_(cell,basis) += tmp_field(basis);
391 template<
typename EvalT,
typename Traits>
397 using Kokkos::parallel_for;
398 using Kokkos::TeamPolicy;
401 basis_ = this->wda(workset).bases[basisIndex_]->weighted_div_basis;
405 if (use_shared_memory) {
407 if (Sacado::IsADType<ScalarT>::value) {
408 const int fadSize = Kokkos::dimension_scalar(field_.get_view());
409 bytes = scratch_view::shmem_size(1,fadSize) + scratch_view::shmem_size(basis_.extent(1),fadSize);
412 bytes = scratch_view::shmem_size(1) + scratch_view::shmem_size(basis_.extent(1));
417 if (fieldMults_.size() == 0) {
419 parallel_for(this->getName(), policy, *
this);
420 }
else if (fieldMults_.size() == 1) {
422 parallel_for(this->getName(), policy, *
this);
425 parallel_for(this->getName(), policy, *
this);
432 if (fieldMults_.size() == 0) {
434 parallel_for(this->getName(), policy, *
this);
435 }
else if (fieldMults_.size() == 1) {
437 parallel_for(this->getName(), policy, *
this);
440 parallel_for(this->getName(), policy, *
this);
450 template<
typename EvalT,
typename TRAITS>
464 p->set<
string>(
"Residual Name",
"?");
465 p->set<
string>(
"Value Name",
"?");
466 p->set<
string>(
"Test Field Name",
"?");
468 p->set(
"Basis", basis);
471 p->set<
double>(
"Multiplier", 1.0);
473 p->set(
"Field Multipliers", fms);
479 #endif // __Panzer_Integrator_DivBasisTimesScalar_impl_hpp__
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.
int num_cells
DEPRECATED - use: numCells()
void evaluateFields(typename Traits::EvalData d)
Evaluate Fields.
Kokkos::TeamPolicy< TeamPolicyProperties...> teamPolicy(const int &league_size)
Returns a TeamPolicy for hierarchic parallelism.
typename EvalT::ScalarT ScalarT
The scalar type.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
EvaluatorStyle
An indication of how an Evaluator will behave.
Teuchos::RCP< const PureBasis > getBasis() const
const panzer::EvaluatorStyle evalStyle_
An enum determining the behavior of this Evaluator.
double multiplier
The scalar multiplier out in front of the integral ( ).
PHX::MDField< const ScalarT, Cell, IP > scalar_
A field representing the scalar function we're integrating ( ).
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< PHX::DataLayout > dl_scalar
Data layout for scalar fields.
panzer::EvaluatorStyle evalStyle
The EvaluatorStyle of the parent Integrator_CurlBasisDotVector object.
This empty struct allows us to optimize operator()() depending on the number of field multipliers...
KOKKOS_INLINE_FUNCTION void operator()(const FieldMultTag< NUM_FIELD_MULT > &tag, const Kokkos::TeamPolicy< PHX::exec_space >::member_type &team) const
Perform the integration. Main memory version.
bool useSharedMemory() const
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.
PHX::MDField< ScalarT, Cell, BASIS > field_
A field to which we'll contribute, or in which we'll store, the result of computing this integral...
Teuchos::RCP< Teuchos::ParameterList > getValidParameters() const
Get Valid Parameters.
void validateParameters(ParameterList const &validParamList, int const depth=1000, EValidateUsed const validateUsed=VALIDATE_USED_ENABLED, EValidateDefaults const validateDefaults=VALIDATE_DEFAULTS_ENABLED) const
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &fm)
Post-Registration Setup.
Integrator_DivBasisTimesScalar(const panzer::EvaluatorStyle &evalStyle, const std::string &resName, const std::string &valName, const panzer::BasisIRLayout &basis, const panzer::IntegrationRule &ir, const double &multiplier=1, const std::vector< std::string > &fmNames=std::vector< std::string >())
Main Constructor.
static HP & inst()
Private ctor.
const std::vector< std::pair< int, LocalOrdinal > > &pid_and_lid const
Description and data layouts associated with a particular basis.
PHX::View< PHX::UnmanagedView< const ScalarT ** > * > kokkosFieldMults_
The PHX::View representation of the (possibly empty) list of fields that are multipliers out in front...
This empty struct allows us to optimize operator()() depending on the number of field multipliers...
Teuchos::RCP< PHX::DataLayout > functional
<Cell,Basis>
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_
std::vector< PHX::MDField< const ScalarT, Cell, IP > > fieldMults_
The (possibly empty) list of fields that are multipliers out in front of the integral ( ...