43 #ifndef PANZER_DOF_IMPL_HPP 44 #define PANZER_DOF_IMPL_HPP 53 #include "Intrepid2_FunctionSpaceTools.hpp" 66 template<
typename EvalT,
typename TRAITS>
70 const std::string fieldName = p.
get<std::string>(
"Name");
73 is_vector_basis = basis->isVectorBasis();
75 std::string evalName = fieldName+
"_"+pointRule->getName();
76 if(p.
isType<
bool>(
"Use DOF Name")) {
77 if(p.
get<
bool>(
"Use DOF Name"))
81 dof_basis = PHX::MDField<const ScalarT,Cell,Point>(fieldName, basis->functional);
83 this->addDependentField(dof_basis);
88 basisValues->setupArrays(layout,
false);
92 if(basis->isScalarBasis()) {
93 dof_ip_scalar = PHX::MDField<ScalarT,Cell,Point>(
95 pointRule->dl_scalar);
96 this->addEvaluatedField(dof_ip_scalar);
97 constBasisRefScalar_ = basisValues->basis_ref_scalar;
98 constBasisScalar_ = basisValues->basis_scalar;
99 this->addDependentField(constBasisRefScalar_);
100 this->addDependentField(constBasisScalar_);
102 else if(basis->isVectorBasis()) {
103 dof_ip_vector = PHX::MDField<ScalarT,Cell,Point,Dim>(
105 pointRule->dl_vector);
106 this->addEvaluatedField(dof_ip_vector);
107 constBasisRefVector_ = basisValues->basis_ref_vector;
108 constBasisVector_ = basisValues->basis_vector;
109 this->addDependentField(constBasisRefVector_);
110 this->addDependentField(constBasisVector_);
115 std::string n =
"DOF_PointValues: " + dof_basis.fieldTag().name();
120 template<
typename EvalT,
typename TRAITS>
125 if(!is_vector_basis) {
126 this->utils.setFieldData(basisValues->basis_ref_scalar,fm);
127 this->utils.setFieldData(basisValues->basis_scalar,fm);
130 this->utils.setFieldData(basisValues->basis_ref_vector,fm);
131 this->utils.setFieldData(basisValues->basis_vector,fm);
136 template<
typename EvalT,
typename TRAITS>
140 if(is_vector_basis) {
141 int spaceDim = basisValues->basis_vector.extent(3);
144 Kokkos::parallel_for(workset.num_cells,functor);
148 Kokkos::parallel_for(workset.num_cells,functor);
153 Kokkos::parallel_for(workset.num_cells,functor);
164 template<
typename TRAITS>
168 const std::string fieldName = p.
get<std::string>(
"Name");
177 offsets_array = Kokkos::View<int*,PHX::Device>(
"offsets",
offsets.size());
178 for(std::size_t i=0;i<
offsets.size();i++)
181 accelerate_jacobian =
true;
184 accelerate_jacobian =
false;
186 std::string evalName = fieldName+
"_"+pointRule->getName();
187 if(p.
isType<
bool>(
"Use DOF Name")) {
188 if(p.
get<
bool>(
"Use DOF Name"))
189 evalName = fieldName;
192 dof_basis = PHX::MDField<const ScalarT,Cell,Point>(fieldName,
basis->functional);
203 if(
basis->isScalarBasis()) {
206 pointRule->dl_scalar);
213 else if(
basis->isVectorBasis()) {
216 pointRule->dl_vector);
226 std::string n =
"DOF_PointValues: " +
dof_basis.fieldTag().name() +
" Jacobian";
231 template<
typename TRAITS>
237 this->utils.setFieldData(
basisValues->basis_ref_scalar,fm);
238 this->utils.setFieldData(
basisValues->basis_scalar,fm);
241 this->utils.setFieldData(
basisValues->basis_ref_vector,fm);
242 this->utils.setFieldData(
basisValues->basis_vector,fm);
247 template<
typename TRAITS>
252 if(accelerate_jacobian) {
253 int spaceDim =
basisValues->basis_vector.extent(3);
256 Kokkos::parallel_for(workset.num_cells,functor);
260 Kokkos::parallel_for(workset.num_cells,functor);
264 int spaceDim =
basisValues->basis_vector.extent(3);
267 Kokkos::parallel_for(workset.num_cells,functor);
271 Kokkos::parallel_for(workset.num_cells,functor);
276 if(accelerate_jacobian) {
278 Kokkos::parallel_for(workset.num_cells,functor);
282 Kokkos::parallel_for(workset.num_cells,functor);
Teuchos::RCP< BasisValues2< ScalarT > > basisValues
T & get(const std::string &name, T def_value)
void evaluateFields(typename TRAITS::EvalData d)
PHX::MDField< const ScalarT, BASIS, IP, void, void, void, void, void, void > constBasisRefScalar_
PHX::MDField< ScalarT, Cell, Point, Dim > dof_ip_vector
PHX::MDField< const ScalarT, BASIS, IP, Dim, void, void, void, void, void > constBasisRefVector_
PHX::MDField< const ScalarT, Cell, Point > dof_basis
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
Teuchos::RCP< const PureBasis > basis
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &fm)
bool isType(const std::string &name) const
PHX::MDField< const ScalarT, Cell, BASIS, IP, void, void, void, void, void > constBasisScalar_
PHX::MDField< const ScalarT, Cell, BASIS, IP, Dim, void, void, void, void > constBasisVector_
PHX::MDField< ScalarT, Cell, Point > dof_ip_scalar
#define TEUCHOS_ASSERT(assertion_test)
DOF_PointValues(const Teuchos::ParameterList &p)
Kokkos::View< const int *, PHX::Device > offsets