43 #ifndef PANZER_DIRICHLET_RESIDUAL_EDGEBASIS_IMPL_HPP 44 #define PANZER_DIRICHLET_RESIDUAL_EDGEBASIS_IMPL_HPP 50 #include "Intrepid2_CellTools.hpp" 53 #include "Kokkos_ViewFactory.hpp" 60 std::string residual_name = p.get<std::string>(
"Residual Name");
62 basis = p.get<Teuchos::RCP<const panzer::PureBasis> >(
"Basis");
63 pointRule = p.get<Teuchos::RCP<const panzer::PointRule> >(
"Point Rule");
65 std::string field_name = p.get<std::string>(
"DOF Name");
66 std::string dof_name = field_name+
"_"+
pointRule->getName();
67 std::string value_name = p.get<std::string>(
"Value Name");
69 Teuchos::RCP<PHX::DataLayout> basis_layout =
basis->functional;
70 Teuchos::RCP<PHX::DataLayout> vector_layout_dof =
pointRule->dl_vector;
71 Teuchos::RCP<PHX::DataLayout> vector_layout_vector =
basis->functional_grad;
74 TEUCHOS_ASSERT(
basis->isVectorBasis());
75 TEUCHOS_ASSERT(basis_layout->dimension(0)==vector_layout_dof->dimension(0));
76 TEUCHOS_ASSERT(basis_layout->dimension(1)==vector_layout_dof->dimension(1));
77 TEUCHOS_ASSERT(Teuchos::as<unsigned>(
basis->dimension())==vector_layout_dof->dimension(2));
78 TEUCHOS_ASSERT(vector_layout_vector->dimension(0)==vector_layout_dof->dimension(0));
79 TEUCHOS_ASSERT(vector_layout_vector->dimension(1)==vector_layout_dof->dimension(1));
80 TEUCHOS_ASSERT(vector_layout_vector->dimension(2)==vector_layout_dof->dimension(2));
82 residual = PHX::MDField<ScalarT,Cell,BASIS>(residual_name, basis_layout);
83 dof = PHX::MDField<const ScalarT,Cell,Point,Dim>(dof_name, vector_layout_dof);
84 value = PHX::MDField<const ScalarT,Cell,Point,Dim>(value_name, vector_layout_vector);
87 std::string orientationFieldName =
basis->name() +
" Orientation";
88 dof_orientation = PHX::MDField<const ScalarT,Cell,BASIS>(orientationFieldName,
102 this->addDependentField(
dof);
104 this->addDependentField(
value);
106 std::string n =
"Dirichlet Residual Edge Basis Evaluator";
113 this->utils.setFieldData(
residual,fm);
114 this->utils.setFieldData(
dof,fm);
116 this->utils.setFieldData(
value,fm);
125 if(workset.num_cells<=0)
130 if(workset.subcell_dim==1) {
131 Intrepid2::CellTools<ScalarT>::getPhysicalEdgeTangents(
edgeTan,
133 this->wda(workset).subcell_index,
134 *
basis->getCellTopology());
136 for(index_t c=0;c<workset.num_cells;c++) {
137 for(
int b=0;b<
dof.extent_int(1);b++) {
138 for(
int d=0;d<
dof.extent_int(2);d++)
143 else if(workset.subcell_dim==2) {
146 const shards::CellTopology & parentCell = *
basis->getCellTopology();
147 int cellDim = parentCell.getDimension();
148 int numEdges =
dof.extent_int(1);
152 for(
int i=0;i<numEdges;i++) {
153 Kokkos::DynRankView<double,PHX::Device> refEdgeTan_local(
"refEdgeTan_local",cellDim);
154 Intrepid2::CellTools<double>::getReferenceEdgeTangent(refEdgeTan_local, i, parentCell);
156 for(
int d=0;d<cellDim;d++)
161 for(index_t c=0;c<workset.num_cells;c++) {
162 for(
int pt = 0; pt < numEdges; pt++) {
165 for(
int i = 0; i < cellDim; i++) {
167 for(
int j = 0; j < cellDim; j++){
174 for(index_t c=0;c<workset.num_cells;c++) {
175 for(
int b=0;b<
dof.extent_int(1);b++) {
176 for(
int d=0;d<
dof.extent_int(2);d++)
184 TEUCHOS_ASSERT(
false);
190 for(index_t c=0;c<workset.num_cells;c++) {
191 for(
int b=0;b<
dof.extent_int(1);b++) {
PHX::MDField< ScalarT > residual
Evaluates a Dirichlet BC residual corresponding to a field value.
PointValues2< ScalarT, PHX::MDField > pointValues
Interpolates basis DOF values to IP DOF values.
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
PHX::MDField< const ScalarT > dof
Teuchos::RCP< const panzer::PointRule > pointRule
PHX_EVALUATOR_CTOR(BasisValues_Evaluator, p)
Kokkos::DynRankView< ScalarT, PHX::Device > refEdgeTan
PHX_EVALUATE_FIELDS(BasisValues_Evaluator, workset)
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.
PHX_POST_REGISTRATION_SETUP(BasisValues_Evaluator, sd, fm)