43#ifndef PANZER_DOF_CURL_IMPL_HPP
44#define PANZER_DOF_CURL_IMPL_HPP
49#include "Intrepid2_FunctionSpaceTools.hpp"
50#include "Phalanx_KokkosDeviceTypes.hpp"
57template <
typename ScalarT,
typename Array,
int spaceDim>
58class EvaluateCurlWithSens_Vector {
67 typedef typename PHX::Device execution_space;
69 EvaluateCurlWithSens_Vector(PHX::MDField<const ScalarT,Cell,Point> in_dof_value,
70 PHX::MDField<ScalarT,Cell,Point,Dim> in_dof_curl,
77 KOKKOS_INLINE_FUNCTION
78 void operator()(
const unsigned int cell)
const
81 for (
int d=0; d<spaceDim; d++) {
92template <
typename ScalarT,
typename ArrayT>
93void evaluateCurl_withSens_vector(
int numCells,
94 PHX::MDField<ScalarT,Cell,Point,Dim> &
dof_curl,
95 PHX::MDField<const ScalarT,Cell,Point> &
dof_value,
104 for (
int cell=0; cell<numCells; cell++) {
106 for (
int d=0; d<spaceDim; d++) {
119template <
typename ScalarT,
typename Array>
120class EvaluateCurlWithSens_Scalar {
121 PHX::MDField<const ScalarT,Cell,Point>
dof_value;
122 PHX::MDField<ScalarT,Cell,Point>
dof_curl;
129 typedef typename PHX::Device execution_space;
131 EvaluateCurlWithSens_Scalar(PHX::MDField<const ScalarT,Cell,Point> in_dof_value,
132 PHX::MDField<ScalarT,Cell,Point> in_dof_curl,
139 KOKKOS_INLINE_FUNCTION
140 void operator()(
const unsigned int cell)
const
152template <
typename ScalarT,
typename ArrayT>
153void evaluateCurl_withSens_scalar(
int numCells,
154 PHX::MDField<ScalarT,Cell,Point> &
dof_curl,
155 PHX::MDField<const ScalarT,Cell,Point> &
dof_value,
163 for (
int cell=0; cell<numCells; cell++) {
176template <
typename ScalarT,
typename Array,
int spaceDim>
177class EvaluateCurlFastSens_Vector {
178 PHX::MDField<const ScalarT,Cell,Point>
dof_value;
179 PHX::MDField<ScalarT,Cell,Point,Dim>
dof_curl;
187 typedef typename PHX::Device execution_space;
189 EvaluateCurlFastSens_Vector(PHX::MDField<const ScalarT,Cell,Point> in_dof_value,
190 PHX::MDField<ScalarT,Cell,Point,Dim> in_dof_curl,
191 PHX::View<const int*> in_offsets,
198 KOKKOS_INLINE_FUNCTION
199 void operator()(
const unsigned int cell)
const
202 for (
int d=0; d<spaceDim; d++) {
215template <
typename ScalarT,
typename ArrayT>
216void evaluateCurl_fastSens_vector(
int numCells,
217 PHX::MDField<ScalarT,Cell,Point,Dim> &
dof_curl,
218 PHX::MDField<const ScalarT,Cell,Point> &
dof_value,
219 const std::vector<int> &
offsets,
227 for (
int cell=0; cell<numCells; cell++) {
229 for (
int d=0; d<spaceDim; d++) {
245template <
typename ScalarT,
typename Array>
246class EvaluateCurlFastSens_Scalar {
247 PHX::MDField<const ScalarT,Cell,Point>
dof_value;
248 PHX::MDField<ScalarT,Cell,Point>
dof_curl;
256 typedef typename PHX::Device execution_space;
258 EvaluateCurlFastSens_Scalar(PHX::MDField<const ScalarT,Cell,Point> in_dof_value,
259 PHX::MDField<ScalarT,Cell,Point> in_dof_curl,
260 PHX::View<const int*> in_offsets,
267 KOKKOS_INLINE_FUNCTION
268 void operator()(
const unsigned int cell)
const
282template <
typename ScalarT,
typename ArrayT>
283void evaluateCurl_fastSens_scalar(
int numCells,
284 PHX::MDField<ScalarT,Cell,Point> &
dof_curl,
285 PHX::MDField<const ScalarT,Cell,Point> &
dof_value,
286 const std::vector<int> &
offsets,
293 for (
int cell=0; cell<numCells; cell++) {
317template<
typename EvalT,
typename TRAITS>
319DOFCurl(
const Teuchos::ParameterList & p) :
320 use_descriptors_(false),
325 Teuchos::RCP<const PureBasis> basis
326 = p.get< Teuchos::RCP<BasisIRLayout> >(
"Basis")->getBasis();
329 TEUCHOS_TEST_FOR_EXCEPTION(!basis->supportsCurl(),std::logic_error,
330 "DOFCurl: Basis of type \"" << basis->name() <<
"\" does not support CURL");
331 TEUCHOS_TEST_FOR_EXCEPTION(!basis->requiresOrientations(),std::logic_error,
332 "DOFCurl: Basis of type \"" << basis->name() <<
"\" in DOF Curl should require orientations. So we are throwing.");
337 dof_curl_scalar = PHX::MDField<ScalarT,Cell,Point>(p.get<std::string>(
"Curl Name"),
338 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_scalar );
342 dof_curl_vector = PHX::MDField<ScalarT,Cell,Point,Dim>(p.get<std::string>(
"Curl Name"),
343 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_vector );
347 { TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"DOFCurl only works for 2D and 3D basis functions"); }
357template<
typename EvalT,
typename TRAITS>
359DOFCurl(
const PHX::FieldTag & input,
360 const PHX::FieldTag & output,
364 : use_descriptors_(true)
383 { TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"DOFCurl only works for 2D and 3D basis functions"); }
393template<
typename EvalT,
typename TRAITS>
399 if(basis_dimension==3)
400 this->utils.setFieldData(dof_curl_vector,fm);
402 this->utils.setFieldData(dof_curl_scalar,fm);
404 if(not use_descriptors_)
409template<
typename EvalT,
typename TRAITS>
414 : *this->wda(workset).bases[basis_index];
416 if(basis_dimension==3) {
417 EvaluateCurlWithSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,3> functor(
dof_value,dof_curl_vector,basisValues.
curl_basis_vector);
418 Kokkos::parallel_for(workset.num_cells,functor);
421 EvaluateCurlWithSens_Scalar<ScalarT,typename BasisValues2<double>::Array_CellBasisIP> functor(
dof_value,dof_curl_scalar,basisValues.
curl_basis_scalar);
422 Kokkos::parallel_for(workset.num_cells,functor);
433template<
typename TRAITS>
435DOFCurl(
const Teuchos::ParameterList & p) :
436 use_descriptors_(false),
441 Teuchos::RCP<const PureBasis> basis
442 = p.get< Teuchos::RCP<BasisIRLayout> >(
"Basis")->getBasis();
446 if(p.isType<Teuchos::RCP<
const std::vector<int> > >(
"Jacobian Offsets Vector")) {
447 offsets = *p.get<Teuchos::RCP<const std::vector<int> > >(
"Jacobian Offsets Vector");
450 PHX::View<int*> offsets_array_nc(
"offsets",
offsets.size());
451 auto offsets_array_nc_h = Kokkos::create_mirror_view(offsets_array_nc);
452 for(std::size_t i=0;i<
offsets.size();i++)
453 offsets_array_nc_h(i) =
offsets[i];
454 Kokkos::deep_copy(offsets_array_nc, offsets_array_nc_h);
455 offsets_array = offsets_array_nc;
457 accelerate_jacobian =
true;
460 accelerate_jacobian =
false;
463 TEUCHOS_TEST_FOR_EXCEPTION(!basis->supportsCurl(),std::logic_error,
464 "DOFCurl: Basis of type \"" << basis->name() <<
"\" does not support CURL");
465 TEUCHOS_TEST_FOR_EXCEPTION(!basis->requiresOrientations(),std::logic_error,
466 "DOFCurl: Basis of type \"" << basis->name() <<
"\" in DOF Curl should require orientations. So we are throwing.");
471 dof_curl_scalar = PHX::MDField<ScalarT,Cell,Point>(p.get<std::string>(
"Curl Name"),
472 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_scalar );
476 dof_curl_vector = PHX::MDField<ScalarT,Cell,Point,Dim>(p.get<std::string>(
"Curl Name"),
477 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_vector );
481 { TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"DOFCurl only works for 2D and 3D basis functions"); }
491template<
typename TRAITS>
493DOFCurl(
const PHX::FieldTag & input,
494 const PHX::FieldTag & output,
498 : use_descriptors_(true)
507 accelerate_jacobian =
false;
519 { TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"DOFCurl only works for 2D and 3D basis functions"); }
529template<
typename TRAITS>
535 if(basis_dimension==3)
536 this->utils.setFieldData(dof_curl_vector,fm);
538 this->utils.setFieldData(dof_curl_scalar,fm);
540 if(not use_descriptors_)
544template<
typename TRAITS>
549 : *this->wda(workset).bases[basis_index];
551 if(!accelerate_jacobian) {
552 if(basis_dimension==3) {
555 EvaluateCurlWithSens_Vector<ScalarT,Array,3> functor(
dof_value,dof_curl_vector,curl_basis_vector);
556 Kokkos::parallel_for(workset.num_cells,functor);
561 EvaluateCurlWithSens_Scalar<ScalarT,Array> functor(
dof_value,dof_curl_scalar,curl_basis_scalar);
562 Kokkos::parallel_for(workset.num_cells,functor);
569 if(basis_dimension==3) {
572 EvaluateCurlFastSens_Vector<ScalarT,Array,3> functor(
dof_value,dof_curl_vector,offsets_array,curl_basis_vector);
573 Kokkos::parallel_for(workset.num_cells,functor);
578 EvaluateCurlFastSens_Scalar<ScalarT,Array> functor(
dof_value,dof_curl_scalar,offsets_array,curl_basis_scalar);
579 Kokkos::parallel_for(workset.num_cells,functor);
PHX::MDField< ScalarT, Cell, Point, Dim > dof_curl
PHX::View< const int * > offsets
PHX::MDField< const ScalarT, Cell, Point > dof_value
const std::string & getType() const
Get type of basis.
PHX::MDField< const Scalar, Cell, BASIS, IP > ConstArray_CellBasisIP
Array_CellBasisIPDim curl_basis_vector
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.
Array_CellBasisIP curl_basis_scalar
ConstArray_CellBasisIP getBasisValues(const bool weighted, const bool cache=true, const bool force=false) const
Get the basis values evaluated at mesh points.
PHX::MDField< const Scalar, Cell, BASIS, IP, Dim > ConstArray_CellBasisIPDim
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.
DOFCurl(const Teuchos::ParameterList &p)
PHX::MDField< ScalarT, Cell, Point > dof_curl_scalar
void evaluateFields(typename TRAITS::EvalData d)
PHX::MDField< const ScalarT, Cell, Point > dof_value
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &fm)
panzer::BasisDescriptor bd_
PHX::MDField< ScalarT, Cell, Point, Dim > dof_curl_vector
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.