43 #ifndef PANZER_GATHER_TANGENTS_IMPL_HPP 44 #define PANZER_GATHER_TANGENTS_IMPL_HPP 46 #include "Teuchos_Assert.hpp" 47 #include "Phalanx_DataLayout.hpp" 51 #include "Kokkos_ViewFactory.hpp" 53 #include "Teuchos_FancyOStream.hpp" 55 template<
typename EvalT,
typename Traits>
60 dof_name = (p.
get< std::string >(
"DOF Name"));
76 std::string orientationFieldName =
basis->
name() +
" Orientation";
77 dof_orientation = PHX::MDField<const ScalarT,Cell,NODE>(orientationFieldName, basis_layout);
89 gatherFieldTangents = PHX::MDField<ScalarT,Cell,NODE,Dim>(dof_name+
"_Tangents",vector_layout_vector);
90 this->addEvaluatedField(gatherFieldTangents);
92 this->setName(
"Gather Tangents");
96 template<
typename EvalT,
typename Traits>
102 this->utils.setFieldData(gatherFieldTangents,fm);
106 edgeTan =
Kokkos::createDynRankView(gatherFieldTangents.get_static_view(),
"edgeTan",gatherFieldTangents.dimension(0),gatherFieldTangents.dimension(1),gatherFieldTangents.dimension(2));
110 template<
typename EvalT,
typename Traits>
119 int cellDim = parentCell.getDimension();
120 int numEdges = gatherFieldTangents.dimension(1);
124 for(
int i=0;i<numEdges;i++) {
125 Kokkos::DynRankView<double,PHX::Device> refEdgeTan_local(
"refEdgeTan_local",cellDim);
126 Intrepid2::CellTools<double>::getReferenceEdgeTangent(refEdgeTan_local, i, parentCell);
128 for(
int d=0;d<cellDim;d++)
133 for(index_t c=0;c<workset.
num_cells;c++) {
134 for(
int pt = 0; pt < numEdges; pt++) {
137 for(
int i = 0; i < cellDim; i++) {
139 for(
int j = 0; j < cellDim; j++){
147 for(index_t c=0;c<workset.
num_cells;c++) {
148 for(
int b=0;b<gatherFieldTangents.extent_int(1);b++) {
149 for(
int d=0;d<gatherFieldTangents.extent_int(2);d++) {
std::string name() const
A unique key that is the combination of the basis type and basis order.
T & get(const std::string &name, T def_value)
PointValues2< ScalarT, PHX::MDField > pointValues
Interpolates basis DOF values to IP DOF values.
bool isVectorBasis() const
Kokkos::DynRankView< ScalarT, PHX::Device > edgeTan
Kokkos::DynRankView< typename InputArray::value_type, PHX::Device > createDynRankView(const InputArray &a, const std::string &name, const DimensionPack... dims)
Wrapper to simplify Panzer use of Sacado ViewFactory.
PHX::MDField< const ScalarT, Cell, BASIS > dof_orientation
Teuchos::RCP< const panzer::PointRule > pointRule
Teuchos::RCP< PHX::DataLayout > functional_grad
<Cell,Basis,Dim>
bool isType(const std::string &name) const
Kokkos::DynRankView< ScalarT, PHX::Device > refEdgeTan
Teuchos::RCP< const shards::CellTopology > getCellTopology() const
Teuchos::RCP< const panzer::PureBasis > basis
Interpolates basis DOF values to IP DOF values.
Array< Scalar, Cell, IP, Dim, Dim, void, void, void, void > jac
void setupArrays(const Teuchos::RCP< const panzer::PointRule > &pr, const ArrayFactory &af)
Sizes/allocates memory for arrays.
const std::string & getName() const
#define TEUCHOS_ASSERT(assertion_test)
Teuchos::RCP< PHX::DataLayout > functional
<Cell,Basis> or <Cell,Basis>
void evaluateFields(typename Traits::EvalData d)
void postRegistrationSetup(typename Traits::SetupData d, PHX::FieldManager< Traits > &vm)