43 #ifndef PANZER_DOF_DIV_IMPL_HPP 44 #define PANZER_DOF_DIV_IMPL_HPP 49 #include "Intrepid2_FunctionSpaceTools.hpp" 56 template<
typename ScalarT,
typename ArrayT>
57 void evaluateDiv_withSens(
int numCells,
58 PHX::MDField<ScalarT,Cell,IP> & dof_div,
59 PHX::MDField<const ScalarT,Cell,Point> &
dof_value,
60 const ArrayT & div_basis)
68 for (
int cell=0; cell<numCells; cell++) {
73 dof_div(cell,pt) =
dof_value(cell, 0) * div_basis(cell, 0, pt);
75 dof_div(cell,pt) +=
dof_value(cell, bf) * div_basis(cell, bf, pt);
88 template<
typename EvalT,
typename TRAITS>
91 use_descriptors_(false),
101 "DOFDiv: Basis of type \"" << basis->name() <<
"\" does not support DIV");
103 "DOFDiv: Basis of type \"" << basis->name() <<
"\" in DOF Div should require orientations. So we are throwing.");
106 dof_div = PHX::MDField<ScalarT,Cell,IP>(p.
get<std::string>(
"Div Name"),
110 this->addEvaluatedField(
dof_div);
113 std::string n =
"DOFDiv: " +
dof_div.fieldTag().name() +
" ("+PHX::typeAsString<EvalT>()+
")";
118 template<
typename EvalT,
typename TRAITS>
121 const PHX::FieldTag & output,
124 : use_descriptors_(true)
135 this->addEvaluatedField(
dof_div);
138 std::string n =
"DOFDiv: " +
dof_div.fieldTag().name() +
" ("+PHX::typeAsString<EvalT>()+
")";
143 template<
typename EvalT,
typename TRAITS>
149 this->utils.setFieldData(dof_div,fm);
151 if(not use_descriptors_)
156 template<
typename EvalT,
typename TRAITS>
161 : *this->wda(workset).bases[basis_index];
163 evaluateDiv_withSens<ScalarT>(workset.num_cells,dof_div,
dof_value,basisValues.
div_basis);
173 template<
typename TRAITS>
176 use_descriptors_(false),
188 accelerate_jacobian =
true;
191 accelerate_jacobian =
false;
192 accelerate_jacobian =
false;
196 "DOFDiv: Basis of type \"" << basis->name() <<
"\" does not support DIV");
198 "DOFDiv: Basis of type \"" << basis->name() <<
"\" in DOF Div should require orientations. So we are throwing.");
201 dof_div = PHX::MDField<ScalarT,Cell,IP>(p.
get<std::string>(
"Div Name"),
205 this->addEvaluatedField(
dof_div);
208 std::string n =
"DOFDiv: " +
dof_div.fieldTag().name() +
" ("+PHX::typeAsString<panzer::Traits::Jacobian>()+
")";
213 template<
typename TRAITS>
216 const PHX::FieldTag & output,
219 : use_descriptors_(true)
229 accelerate_jacobian =
false;
232 this->addEvaluatedField(
dof_div);
235 std::string n =
"DOFDiv: " +
dof_div.fieldTag().name() +
" ("+PHX::typeAsString<panzer::Traits::Jacobian>()+
")";
240 template<
typename TRAITS>
246 this->utils.setFieldData(
dof_div,fm);
252 template<
typename TRAITS>
259 if(!accelerate_jacobian) {
panzer::BasisDescriptor bd_
panzer::IntegrationDescriptor id_
DOFDiv(const Teuchos::ParameterList &p)
T & get(const std::string &name, T def_value)
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void postRegistrationSetup(typename TRAITS::SetupData d, PHX::FieldManager< TRAITS > &fm)
Array_CellBasisIP div_basis
Teuchos::RCP< PHX::DataLayout > dl_scalar
Data layout for scalar fields.
PHX::MDField< const ScalarT, Cell, Point > dof_value
bool isType(const std::string &name) const
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)
WorksetDetailsAccessor wda
#define TEUCHOS_ASSERT(assertion_test)
PHX::MDField< const ScalarT, Cell, Point > dof_value
const std::string & getType() const
Get type of basis.
Kokkos::View< const int *, PHX::Device > offsets
PHX::MDField< ScalarT, Cell, IP > dof_div