43 #ifndef PANZER_DOF_CURL_IMPL_HPP 44 #define PANZER_DOF_CURL_IMPL_HPP 49 #include "Intrepid2_FunctionSpaceTools.hpp" 50 #include "Phalanx_KokkosDeviceTypes.hpp" 57 template <
typename ScalarT,
typename Array,
int spaceDim>
58 class 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
81 for (
int d=0; d<spaceDim; d++) {
92 template <
typename ScalarT,
typename ArrayT>
93 void 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++) {
119 template <
typename ScalarT,
typename Array>
120 class 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 152 template <
typename ScalarT,
typename ArrayT>
153 void 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++) {
176 template <
typename ScalarT,
typename Array,
int spaceDim>
177 class 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 Kokkos::View<const int*,PHX::Device> in_offsets,
198 KOKKOS_INLINE_FUNCTION
199 void operator()(
const unsigned int cell)
const 202 for (
int d=0; d<spaceDim; d++) {
215 template <
typename ScalarT,
typename ArrayT>
216 void 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++) {
245 template <
typename ScalarT,
typename Array>
246 class EvaluateCurlFastSens_Scalar {
247 PHX::MDField<const ScalarT,Cell,Point>
dof_value;
248 PHX::MDField<ScalarT,Cell,Point>
dof_curl;
249 Kokkos::View<const int*,PHX::Device>
offsets;
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 Kokkos::View<const int*,PHX::Device> in_offsets,
267 KOKKOS_INLINE_FUNCTION
268 void operator()(
const unsigned int cell)
const 282 template <
typename ScalarT,
typename ArrayT>
283 void 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++) {
317 template<
typename EvalT,
typename TRAITS>
324 Teuchos::RCP<const PureBasis>
basis 325 = p.get< Teuchos::RCP<BasisIRLayout> >(
"Basis")->getBasis();
328 TEUCHOS_TEST_FOR_EXCEPTION(!
basis->supportsCurl(),std::logic_error,
329 "DOFCurl: Basis of type \"" <<
basis->name() <<
"\" does not support CURL");
330 TEUCHOS_TEST_FOR_EXCEPTION(!
basis->requiresOrientations(),std::logic_error,
331 "DOFCurl: Basis of type \"" <<
basis->name() <<
"\" in DOF Curl should require orientations. So we are throwing.");
336 dof_curl_scalar = PHX::MDField<ScalarT,Cell,Point>(p.get<std::string>(
"Curl Name"),
337 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_scalar );
341 dof_curl_vector = PHX::MDField<ScalarT,Cell,Point,Dim>(p.get<std::string>(
"Curl Name"),
342 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_vector );
346 { TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"DOFCurl only works for 2D and 3D basis functions"); }
356 template<
typename EvalT,
typename TRAITS>
362 if(basis_dimension==3)
363 this->utils.setFieldData(dof_curl_vector,fm);
365 this->utils.setFieldData(dof_curl_scalar,fm);
371 template<
typename EvalT,
typename TRAITS>
377 if(basis_dimension==3) {
378 EvaluateCurlWithSens_Vector<ScalarT,typename BasisValues2<double>::Array_CellBasisIPDim,3> functor(
dof_value,dof_curl_vector,
basisValues.curl_basis_vector);
379 Kokkos::parallel_for(workset.num_cells,functor);
382 EvaluateCurlWithSens_Scalar<ScalarT,typename BasisValues2<double>::Array_CellBasisIP> functor(
dof_value,dof_curl_scalar,
basisValues.curl_basis_scalar);
383 Kokkos::parallel_for(workset.num_cells,functor);
394 template<
typename TRAITS>
401 Teuchos::RCP<const PureBasis>
basis 402 = p.get< Teuchos::RCP<BasisIRLayout> >(
"Basis")->getBasis();
406 if(p.isType<Teuchos::RCP<
const std::vector<int> > >(
"Jacobian Offsets Vector")) {
407 offsets = *p.get<Teuchos::RCP<const std::vector<int> > >(
"Jacobian Offsets Vector");
410 Kokkos::View<int*,PHX::Device> offsets_array_nc
411 = Kokkos::View<int*,PHX::Device>(
"offsets",
offsets.size());
412 for(std::size_t i=0;i<
offsets.size();i++)
413 offsets_array_nc(i) =
offsets[i];
414 offsets_array = offsets_array_nc;
416 accelerate_jacobian =
true;
419 accelerate_jacobian =
false;
422 TEUCHOS_TEST_FOR_EXCEPTION(!
basis->supportsCurl(),std::logic_error,
423 "DOFCurl: Basis of type \"" <<
basis->name() <<
"\" does not support CURL");
424 TEUCHOS_TEST_FOR_EXCEPTION(!
basis->requiresOrientations(),std::logic_error,
425 "DOFCurl: Basis of type \"" <<
basis->name() <<
"\" in DOF Curl should require orientations. So we are throwing.");
430 dof_curl_scalar = PHX::MDField<ScalarT,Cell,Point>(p.get<std::string>(
"Curl Name"),
431 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_scalar );
435 dof_curl_vector = PHX::MDField<ScalarT,Cell,Point,Dim>(p.get<std::string>(
"Curl Name"),
436 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_vector );
440 { TEUCHOS_TEST_FOR_EXCEPTION(
true,std::logic_error,
"DOFCurl only works for 2D and 3D basis functions"); }
450 template<
typename TRAITS>
464 template<
typename TRAITS>
470 if(!accelerate_jacobian) {
473 Kokkos::parallel_for(workset.num_cells,functor);
477 Kokkos::parallel_for(workset.num_cells,functor);
486 Kokkos::parallel_for(workset.num_cells,functor);
490 Kokkos::parallel_for(workset.num_cells,functor);
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.
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &fm)
void operator()(const FieldMultTag< NUM_FIELD_MULT > &, const std::size_t &cell) const
PHX::MDField< ScalarT, Cell, Point > dof_curl_scalar
PHX::MDField< ScalarT, Cell, Point, Dim > dof_curl_vector
Teuchos::RCP< BasisValues2< ScalarT > > basisValues
DOFCurl(const Teuchos::ParameterList &p)
PHX::MDField< ScalarT, Cell, Point > dof_value
Interpolates basis DOF values to IP DOF Gradient values.
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
WorksetDetailsAccessor wda
void evaluateFields(typename TRAITS::EvalData d)
PHX::MDField< ScalarT, Cell, Point, Dim > dof_curl
Kokkos::View< const int *, PHX::Device > offsets
PHX::MDField< const ScalarT, Cell, Point > dof_value