43#ifndef PANZER_DOF_IMPL_HPP
44#define PANZER_DOF_IMPL_HPP
54#include "Intrepid2_FunctionSpaceTools.hpp"
67template<
typename EvalT,
typename TRAITS>
69DOF(
const Teuchos::ParameterList & p) :
70 use_descriptors_(false),
71 dof_basis( p.get<
std::string>(
"Name"),
75 Teuchos::RCP<const PureBasis> basis
76 = p.get< Teuchos::RCP<BasisIRLayout> >(
"Basis")->getBasis();
80 if(basis->isScalarBasis()) {
82 p.get<std::string>(
"Name"),
83 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_scalar);
86 else if(basis->isVectorBasis()) {
88 p.get<std::string>(
"Name"),
89 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_vector);
93 { TEUCHOS_ASSERT(
false); }
97 std::string n =
"DOF: " +
dof_basis.fieldTag().name() +
" ("+PHX::print<EvalT>()+
")";
102template<
typename EvalT,
typename TRAITS>
104DOF(
const PHX::FieldTag & input,
105 const PHX::FieldTag & output,
108 : use_descriptors_(true)
130 std::string n =
"DOF: " +
dof_basis.fieldTag().name() +
" ("+PHX::print<EvalT>()+
")";
135template<
typename EvalT,
typename TRAITS>
140 this->utils.setFieldData(dof_basis,fm);
142 this->utils.setFieldData(dof_ip_vector,fm);
144 this->utils.setFieldData(dof_ip_scalar,fm);
147 if(not use_descriptors_)
152template<
typename EvalT,
typename TRAITS>
157 : *this->wda(workset).bases[basis_index];
162 if(is_vector_basis) {
165 const int spaceDim = array.extent(3);
168 Kokkos::parallel_for(this->getName(),policy,functor);
172 Kokkos::parallel_for(this->getName(),policy,functor);
180 Kokkos::parallel_for(workset.num_cells,functor);
191template<
typename TRAITS>
193DOF(
const Teuchos::ParameterList & p) :
194 use_descriptors_(false),
195 dof_basis( p.get<
std::string>(
"Name"),
199 Teuchos::RCP<const PureBasis> basis
200 = p.get< Teuchos::RCP<BasisIRLayout> >(
"Basis")->getBasis();
203 if(p.isType<Teuchos::RCP<
const std::vector<int> > >(
"Jacobian Offsets Vector")) {
204 const std::vector<int> &
offsets = *p.get<Teuchos::RCP<const std::vector<int> > >(
"Jacobian Offsets Vector");
207 offsets_array = PHX::View<int*>(
"offsets",
offsets.size());
208 auto offsets_array_h = Kokkos::create_mirror_view(offsets_array);
209 for(std::size_t i=0;i<
offsets.size();i++)
210 offsets_array_h(i) =
offsets[i];
211 Kokkos::deep_copy(offsets_array, offsets_array_h);
213 accelerate_jacobian_enabled =
true;
216 sensitivities_name =
true;
217 if (p.isType<std::string>(
"Sensitivities Name"))
218 sensitivities_name = p.get<std::string>(
"Sensitivities Name");
221 accelerate_jacobian_enabled =
false;
224 if(basis->isScalarBasis()) {
226 p.get<std::string>(
"Name"),
227 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_scalar);
230 else if(basis->isVectorBasis()) {
232 p.get<std::string>(
"Name"),
233 p.get< Teuchos::RCP<panzer::IntegrationRule> >(
"IR")->dl_vector);
237 { TEUCHOS_ASSERT(
false); }
241 std::string n =
"DOF: " +
dof_basis.fieldTag().name()
242 + ( accelerate_jacobian_enabled ?
" accel_jac " :
"slow_jac" )
243 +
" ("+PHX::print<panzer::Traits::Jacobian>()+
")";
248template<
typename TRAITS>
250DOF(
const PHX::FieldTag & input,
251 const PHX::FieldTag & output,
254 : use_descriptors_(true)
262 accelerate_jacobian_enabled = false;
278 std::string n =
"DOF: " +
dof_basis.fieldTag().name() +
" slow_jac(descriptor) ("+PHX::print<typename TRAITS::Jacobian>()+
")";
283template<
typename TRAITS>
288 this->utils.setFieldData(dof_basis,fm);
290 this->utils.setFieldData(dof_ip_vector,fm);
292 this->utils.setFieldData(dof_ip_scalar,fm);
295 if(not use_descriptors_)
300template<
typename TRAITS>
301void DOF<typename TRAITS::Jacobian, TRAITS>::
302preEvaluate(
typename TRAITS::PreEvalData d)
306 accelerate_jacobian =
false;
307 if(accelerate_jacobian_enabled && d.first_sensitivities_name==sensitivities_name) {
308 accelerate_jacobian =
true;
313template<
typename TRAITS>
318 : *this->wda(workset).bases[basis_index];
320 if(is_vector_basis) {
321 if(accelerate_jacobian) {
324 const int spaceDim = array.extent(3);
326 dof_functors::EvaluateDOFFastSens_Vector<ScalarT,Array,3> functor(dof_basis,dof_ip_vector,offsets_array,array);
327 Kokkos::parallel_for(workset.num_cells,functor);
330 dof_functors::EvaluateDOFFastSens_Vector<ScalarT,Array,2> functor(dof_basis,dof_ip_vector,offsets_array,array);
331 Kokkos::parallel_for(workset.num_cells,functor);
339 const int spaceDim = array.extent(3);
341 dof_functors::EvaluateDOFWithSens_Vector<ScalarT,Array,3> functor(dof_basis.get_static_view(),dof_ip_vector.get_static_view(),array,use_shared_memory);
342 Kokkos::parallel_for(this->getName(),policy,functor);
345 dof_functors::EvaluateDOFWithSens_Vector<ScalarT,Array,2> functor(dof_basis.get_static_view(),dof_ip_vector.get_static_view(),array,use_shared_memory);
346 Kokkos::parallel_for(this->getName(),policy,functor);
353 if(accelerate_jacobian) {
354 dof_functors::EvaluateDOFFastSens_Scalar<ScalarT,Array> functor(dof_basis,dof_ip_scalar,offsets_array,array);
355 Kokkos::parallel_for(workset.num_cells,functor);
358 dof_functors::EvaluateDOFWithSens_Scalar<ScalarT,Array> functor(dof_basis,dof_ip_scalar,array);
359 Kokkos::parallel_for(workset.num_cells,functor);
PHX::View< const int * > offsets
const std::string & getType() const
Get type of basis.
PHX::MDField< const Scalar, Cell, BASIS, IP > ConstArray_CellBasisIP
Array_CellBasisIP 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
Array_CellBasisIPDim basis_vector
ConstArray_CellBasisIPDim getVectorBasisValues(const bool weighted, const bool cache=true, const bool force=false) const
Get the vector basis values evaluated at mesh points.
PHX::MDField< const ScalarT, Cell, Point > dof_basis
PHX::MDField< ScalarT, Cell, Point, Dim > dof_ip_vector
DOF(const Teuchos::ParameterList &p)
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &fm)
void evaluateFields(typename TRAITS::EvalData d)
PHX::MDField< ScalarT, Cell, Point > dof_ip_scalar
static HP & inst()
Private ctor.
bool useSharedMemory() const
Kokkos::TeamPolicy< TeamPolicyProperties... > teamPolicy(const int &league_size)
Returns a TeamPolicy for hierarchic parallelism.
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.