43 #ifndef PANZER_DOF_GRADIENT_IMPL_HPP 44 #define PANZER_DOF_GRADIENT_IMPL_HPP 49 #include "Intrepid2_FunctionSpaceTools.hpp" 55 template <
typename ScalarT,
typename ArrayT>
56 struct evaluateGrad_withSens {
57 evaluateGrad_withSens (PHX::MDField<ScalarT> dof_grad,
58 PHX::MDField<const ScalarT,Cell,Point>
dof_value,
59 const ArrayT grad_basis) :
62 KOKKOS_INLINE_FUNCTION
63 void operator() (
const size_t &cell)
const 70 for (
int d=0; d<spaceDim; d++) {
88 template<
typename EvalT,
typename Traits>
91 const Teuchos::ParameterList& p) :
92 use_descriptors_(false),
95 dof_gradient( p.get<
std::string>(
"Gradient Name"),
99 Teuchos::RCP<const PureBasis> basis
100 = p.get< Teuchos::RCP<BasisIRLayout> >(
"Basis")->getBasis();
103 TEUCHOS_TEST_FOR_EXCEPTION(!basis->supportsGrad(),std::logic_error,
104 "DOFGradient: Basis of type \"" << basis->name() <<
"\" does not support GRAD");
109 std::string n =
"DOFGradient: " +
dof_gradient.fieldTag().name() +
" ("+PHX::typeAsString<EvalT>()+
")";
114 template<
typename EvalT,
typename TRAITS>
117 const PHX::FieldTag & output,
120 : use_descriptors_(true)
124 , dof_gradient(output)
127 TEUCHOS_TEST_FOR_EXCEPTION(
bd_.
getType()==
"HGrad",std::logic_error,
128 "DOFGradient: Basis of type \"" <<
bd_.
getType() <<
"\" does not support GRAD");
133 std::string n =
"DOFGradient: " +
dof_gradient.fieldTag().name() +
" ("+PHX::typeAsString<EvalT>()+
")";
138 template<
typename EvalT,
typename Traits>
146 this->utils.setFieldData(dof_gradient,fm);
148 if(not use_descriptors_)
153 template<
typename EvalT,
typename Traits>
163 : *this->wda(workset).bases[basis_index];
165 typedef decltype(basisValues.
grad_basis) ArrayT;
167 evaluateGrad_withSens<ScalarT, ArrayT> eval(dof_gradient,
dof_value,basisValues.
grad_basis);
168 Kokkos::parallel_for(workset.
num_cells, eval);
DOFGradient(const Teuchos::ParameterList &p)
PHX::MDField< ScalarT > dof_gradient
PHX::MDField< ScalarT > dof_grad_
panzer::BasisDescriptor bd_
PHX::MDField< const ScalarT, Cell, Point > dof_value
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.
void evaluateFields(typename TRAITS::EvalData d)
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &fm)
PHX::MDField< const ScalarT, Cell, Point > dof_value_
Array_CellBasisIPDim grad_basis
PHX::MDField< const ScalarT, Cell, Point > dof_value
const std::string & getType() const
Get type of basis.
Teuchos::RCP< const std::vector< panzer::Workset > > worksets_