43 #ifndef Panzer_Integrator_CurlBasisDotVector_impl_hpp
44 #define Panzer_Integrator_CurlBasisDotVector_impl_hpp
53 #include "Intrepid2_FunctionSpaceTools.hpp"
62 #include "Phalanx_KokkosDeviceTypes.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 useDescriptors_(false),
85 multiplier_(multiplier),
86 basisName_(basis.name())
96 using std::invalid_argument;
97 using std::logic_error;
103 "Integrator_CurlBasisDotVector called with an empty residual name.")
105 "Integrator_CurlBasisDotVector called with an empty value name.")
108 "Error: Integrator_CurlBasisDotVector: Basis of type \""
109 << tmpBasis->name() <<
"\" is not a vector basis.")
111 logic_error,
"Error: Integrator_CurlBasisDotVector: Basis of type \""
112 << tmpBasis->name() <<
"\" does not require orientations, though it " \
113 "should for its use in this Evaluator.")
115 "Error: Integrator_CurlBasisDotVector: Basis of type \""
116 << tmpBasis->name() <<
"\" does not support curl.")
118 (tmpBasis->dimension() == 3)), logic_error,
119 "Error: Integrator_CurlBasisDotVector requires either a two- or " \
120 "three-dimensional basis. The basis \"" << tmpBasis->name()
142 this->addContributedField(
field_);
144 this->addEvaluatedField(
field_);
149 kokkosFieldMults_ = View<Kokkos::View<const ScalarT**, typename PHX::DevLayout<ScalarT>::type, Kokkos::MemoryUnmanaged>*>(
150 "CurlBasisDotVector::KokkosFieldMultipliers", fmNames.size());
151 for (
const auto& name : fmNames)
158 string n(
"Integrator_CurlBasisDotVector (");
163 n +=
"): " +
field_.fieldTag().name();
172 template<
typename EvalT,
typename Traits>
179 p.get<std::string>(
"Residual Name"),
180 p.get<std::string>(
"Value Name"),
183 p.get<double>(
"Multiplier"),
184 p.isType<Teuchos::
RCP<
const std::vector<std::string>>>
185 (
"Field Multipliers") ?
186 (*p.get<Teuchos::
RCP<
const std::vector<std::string>>>
187 (
"Field Multipliers")) : std::vector<std::string>())
202 template<
typename EvalT,
typename TRAITS>
212 const std::vector<PHX::FieldTag>& multipliers
215 evalStyle_(evalStyle),
216 useDescriptors_(true),
219 multiplier_(multiplier),
224 using std::logic_error;
229 "Error: Integrator_CurlBasisDotVector: Basis of type \""
230 <<
bd_.
getType() <<
"\" does not support curl.")
232 logic_error,
"Error: Integrator_CurlBasisDotVector works on either " \
233 "two- or three-dimensional problems. You've provided spaceDim = "
252 this->addContributedField(
field_);
254 this->addEvaluatedField(
field_);
259 kokkosFieldMults_ = View<Kokkos::View<const ScalarT**, typename PHX::DevLayout<ScalarT>::type, Kokkos::MemoryUnmanaged>*>(
260 "CurlBasisDotVector::KokkosFieldMultipliers", multipliers.size());
261 for (
const auto& fm : multipliers)
268 string n(
"Integrator_CurlBasisDotVector (");
273 n +=
"): " +
field_.fieldTag().name();
282 template<
typename EvalT,
typename Traits>
294 auto kokkosFieldMults_h = Kokkos::create_mirror_view(kokkosFieldMults_);
295 for (
size_t i(0); i < fieldMults_.size(); ++i)
296 kokkosFieldMults_h(i) = fieldMults_[i].get_static_view();
297 Kokkos::deep_copy(kokkosFieldMults_, kokkosFieldMults_h);
300 if (not useDescriptors_)
306 fm.template getKokkosExtendedDataTypeDimensions<EvalT>(),
true);
309 "Integrator_CurlBasisDotVector: 2-D Result", vector2D_.extent(0),
310 vector2D_.extent(1));
313 "Integrator_CurlBasisDotVector: 3-D Result", vector3D_.extent(0),
314 vector3D_.extent(1), 3);
335 template<
typename ScalarT>
344 struct ScalarMultiplierTag
352 struct FieldMultiplierTag
367 KOKKOS_INLINE_FUNCTION
368 void operator()(
const ScalarMultiplierTag,
const unsigned cell)
const
370 int numQP(
result.extent(1));
371 for (
int qp(0); qp < numQP; ++qp)
387 KOKKOS_INLINE_FUNCTION
388 void operator()(
const FieldMultiplierTag,
const unsigned cell)
const
390 int numQP(
result.extent(1));
391 for (
int qp(0); qp < numQP; ++qp)
399 using execution_space = PHX::Device;
437 template<
typename ScalarT>
446 struct ScalarMultiplierTag
454 struct FieldMultiplierTag
469 KOKKOS_INLINE_FUNCTION
470 void operator()(
const ScalarMultiplierTag,
const unsigned cell)
const
473 for (
int qp(0); qp < numQP; ++qp)
474 for (
int dim(0); dim < numDim; ++dim)
490 KOKKOS_INLINE_FUNCTION
491 void operator()(
const FieldMultiplierTag,
const unsigned cell)
const
494 for (
int qp(0); qp < numQP; ++qp)
495 for (
int dim(0); dim < numDim; ++dim)
503 using execution_space = PHX::Device;
542 template<
typename ScalarT,
int spaceDim>
559 KOKKOS_INLINE_FUNCTION
560 void operator()(
const unsigned cell)
const
565 for (
int basis(0); basis < numBases; ++basis)
568 field(cell, basis) = 0.0;
569 for (
int qp(0); qp < numQP; ++qp)
570 for (
int dim(0); dim < spaceDim; ++dim)
580 using execution_space = PHX::Device;
617 template<
typename ScalarT>
633 KOKKOS_INLINE_FUNCTION
634 void operator()(
const unsigned cell)
const
639 for (
int basis(0); basis < numBases; ++basis)
642 field(cell, basis) = 0.0;
643 for (
int qp(0); qp < numQP; ++qp)
653 using execution_space = PHX::Device;
686 template<
typename EvalT,
typename Traits>
692 using Kokkos::parallel_for;
693 using Kokkos::RangePolicy;
702 *this->wda(workset).bases[basisIndex_];
709 using PreMultiply = PreMultiply2D<ScalarT>;
710 using ScalarMultiplierTag =
typename PreMultiply::ScalarMultiplierTag;
711 using FieldMultiplierTag =
typename PreMultiply::FieldMultiplierTag;
712 PreMultiply preMultiply;
713 preMultiply.result = result2D_;
714 preMultiply.multiplier = multiplier_;
715 preMultiply.vectorField = vector2D_;
719 parallel_for(RangePolicy<Device, ScalarMultiplierTag>(0,
724 for (
const auto&
field : fieldMults_)
726 preMultiply.fieldMult =
field;
727 parallel_for(RangePolicy<Device, FieldMultiplierTag>(0,
732 Integrate2D<ScalarT> integrate;
733 integrate.result = result2D_;
734 integrate.field = field_;
737 integrate.evalStyle = evalStyle_;
738 parallel_for(workset.
num_cells, integrate);
744 using PreMultiply = PreMultiply3D<ScalarT>;
745 using ScalarMultiplierTag =
typename PreMultiply::ScalarMultiplierTag;
746 using FieldMultiplierTag =
typename PreMultiply::FieldMultiplierTag;
747 PreMultiply preMultiply;
748 preMultiply.result = result3D_;
749 preMultiply.multiplier = multiplier_;
750 preMultiply.vectorField = vector3D_;
754 parallel_for(RangePolicy<Device, ScalarMultiplierTag>(0,
759 for (
const auto&
field : fieldMults_)
761 preMultiply.fieldMult =
field;
762 parallel_for(RangePolicy<Device, FieldMultiplierTag>(0,
767 Integrate3D<ScalarT, 3> integrate;
768 integrate.result = result3D_;
769 integrate.field = field_;
772 integrate.evalStyle = evalStyle_;
773 parallel_for(workset.
num_cells, integrate);
782 template<
typename EvalT,
typename TRAITS>
796 p->set<
string>(
"Residual Name",
"?");
797 p->set<
string>(
"Value Name",
"?");
799 p->set(
"Basis", basis);
802 p->set<
double>(
"Multiplier", 1.0);
804 p->set(
"Field Multipliers", fms);
810 #endif // Panzer_Integrator_CurlBasisDotVector_impl_hpp
Array_CellBasisIPDim weighted_curl_basis_vector
void evaluateFields(typename Traits::EvalData d)
Evaluate Fields.
int num_cells
DEPRECATED - use: numCells()
std::vector< PHX::MDField< const ScalarT, panzer::Cell, panzer::IP > > fieldMults_
The (possibly empty) list of fields that are multipliers out in front of the integral ( ...
PHX::MDField< const double, panzer::Cell, panzer::BASIS, panzer::IP, panzer::Dim > weightedCurlBasis
The vector basis information necessary for integration.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
PHX::MDField< const ScalarT, panzer::Cell, panzer::IP > vectorField
A field representing the vector-valued function we're integrating ( ).
Teuchos::RCP< Teuchos::ParameterList > getValidParameters() const
Get Valid Parameters.
EvaluatorStyle
An indication of how an Evaluator will behave.
PHX::MDField< Scalar, T0 > buildStaticArray(const std::string &str, int d0) const
ConstArray_CellBasisIPDim getCurlVectorBasis(const bool weighted, const bool cache=true, const bool force=false) const
Get the curl of a 3D vector basis evaluated at mesh points.
Teuchos::RCP< const PureBasis > getBasis() const
panzer::BasisDescriptor bd_
The BasisDescriptor for the basis to use.
PHX::MDField< ScalarT, panzer::Cell, panzer::IP > result
A field that will be used to build up the result of the integral we're performing.
double multiplier
The scalar multiplier out in front of the integral ( ).
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
const std::string & getType() const
Get type of basis.
Teuchos::RCP< PHX::DataLayout > dl_scalar
Data layout for scalar fields.
ConstArray_CellBasisIP getBasisValues(const bool weighted, const bool cache=true, const bool force=false) const
Get the basis values evaluated at mesh points.
panzer::EvaluatorStyle evalStyle
The EvaluatorStyle of the parent Integrator_CurlBasisDotVector object.
PHX::MDField< const ScalarT, panzer::Cell, panzer::IP, panzer::Dim > vector3D_
A field representing the vector-valued function we're integrating ( ) in a three-dimensional problem...
int spaceDim_
The spatial dimension of the vector-valued function we're integrating, either 2 or 3...
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.
Array_CellBasisIP weighted_curl_basis_scalar
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 sd, PHX::FieldManager< Traits > &fm)
Post-Registration Setup.
Teuchos::RCP< PHX::DataLayout > dl_vector
Data layout for vector fields.
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...
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...
PHX::View< Kokkos::View< const ScalarT **, typename PHX::DevLayout< ScalarT >::type, Kokkos::MemoryUnmanaged > * > kokkosFieldMults_
The PHX::View representation of the (possibly empty) list of fields that are multipliers out in front...
Integrator_CurlBasisDotVector(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.
const std::vector< std::pair< int, LocalOrdinal > > &pid_and_lid const
Description and data layouts associated with a particular basis.
PHX::MDField< const ScalarT, panzer::Cell, panzer::IP > fieldMult
One of the field multipliers ( , , etc.) out in front of the integral.
Teuchos::RCP< PHX::DataLayout > functional
<Cell,Basis>
PHX::MDField< const ScalarT, panzer::Cell, panzer::IP > vector2D_
A field representing the vector-valued function we're integrating ( ) in a two-dimensional problem...
ConstArray_CellBasisIP getCurl2DVectorBasis(const bool weighted, const bool cache=true, const bool force=false) const
Get the curl of a 2D vector basis evaluated at mesh points.
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_
typename EvalT::ScalarT ScalarT
The scalar type.