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,
59 PHX::MDField<ScalarT,Cell,IP> & dof_div,
60 PHX::MDField<ScalarT,Cell,Point> &
dof_value,
61 const ArrayT & div_basis)
69 for (
int cell=0; cell<numCells; cell++) {
74 dof_div(cell,pt) =
dof_value(cell, 0) * div_basis(cell, 0, pt);
76 dof_div(cell,pt) +=
dof_value(cell, bf) * div_basis(cell, bf, pt);
89 template<
typename EvalT,
typename TRAITS>
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.");
107 dof_div = PHX::MDField<ScalarT,Cell,IP>(p.
get<std::string>(
"Div Name"),
111 this->addEvaluatedField(
dof_div);
114 std::string n =
"DOFDiv: " +
dof_div.fieldTag().name() +
" ("+PHX::typeAsString<EvalT>()+
")";
119 template<
typename EvalT,
typename TRAITS>
125 this->utils.setFieldData(dof_div,fm);
131 template<
typename EvalT,
typename TRAITS>
135 evaluateDiv_withSens<ScalarT>(workset.num_cells,basis_dimension,dof_div,
dof_value,this->wda(workset).bases[
basis_index]->div_basis);
145 template<
typename TRAITS>
159 accelerate_jacobian =
true;
162 accelerate_jacobian =
false;
163 accelerate_jacobian =
false;
167 "DOFDiv: Basis of type \"" <<
basis->
name() <<
"\" does not support DIV");
169 "DOFDiv: Basis of type \"" <<
basis->
name() <<
"\" in DOF Div should require orientations. So we are throwing.");
173 dof_div = PHX::MDField<ScalarT,Cell,IP>(p.
get<std::string>(
"Div Name"),
177 this->addEvaluatedField(
dof_div);
180 std::string n =
"DOFDiv: " +
dof_div.fieldTag().name() +
" ("+PHX::typeAsString<panzer::Traits::Jacobian>()+
")";
185 template<
typename TRAITS>
191 this->utils.setFieldData(
dof_div,fm);
196 template<
typename TRAITS>
200 if(!accelerate_jacobian) {
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.
std::string name() const
A unique key that is the combination of the basis type and basis order.
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)
bool requiresOrientations() const
Teuchos::RCP< PHX::DataLayout > dl_scalar
Data layout for scalar fields.
PHX::MDField< ScalarT, Cell, Point > dof_value
bool isType(const std::string &name) const
int dimension() const
Returns the dimension of the basis from the topology.
PHX::MDField< ScalarT, Cell, Point > dof_value
Interpolates basis DOF values to IP DOF Gradient values.
void evaluateFields(typename TRAITS::EvalData d)
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
#define TEUCHOS_ASSERT(assertion_test)
Kokkos::View< const int *, PHX::Device > offsets
PHX::MDField< ScalarT, Cell, IP > dof_div