42 #include "Intrepid2_config.h" 43 #if defined( HAVE_INTREPID2_KOKKOS_DYNRANKVIEW ) 46 #include "Teuchos_Assert.hpp" 47 #include "Intrepid2_CellTools.hpp" 48 #include "Shards_CellTopology.hpp" 56 const shards::CellTopology ct =
intrepidBasis_->getBaseCellTopology();
57 return ct.getSubcellCount(dim);
60 const std::vector<int> &
63 subcellIndices_.clear();
70 for (
int i=0;i<ndof;++i)
71 subcellIndices_.push_back(tag(dim, cellIndex, i));
73 return subcellIndices_;
87 const shards::CellTopology ct =
intrepidBasis_->getBaseCellTopology();
89 std::set<std::pair<unsigned,unsigned> > closure;
93 std::set<std::pair<unsigned,unsigned> >::const_iterator itr;
94 for (itr=closure.begin();itr!=closure.end();++itr) {
96 const std::vector<int> & subcellIndices =
getSubcellIndices(itr->first,itr->second);
99 indices.insert(indices.end(),subcellIndices.begin(),subcellIndices.end());
120 getSubcellNodes(
const shards::CellTopology & cellTopo,
unsigned dim,
unsigned subCell,
121 std::vector<unsigned> & nodes)
124 nodes.push_back(subCell);
129 unsigned subCellNodeCount = cellTopo.getNodeCount(dim,subCell);
130 for(
unsigned node=0;node<subCellNodeCount;++node)
131 nodes.push_back(cellTopo.getNodeMap(dim,subCell,node));
134 std::sort(nodes.begin(),nodes.end());
140 const std::vector<unsigned> & nodes,
141 std::set<std::pair<unsigned,unsigned> > & subCells)
143 unsigned subCellCount = cellTopo.getSubcellCount(dim);
144 for(
unsigned subCellOrd=0;subCellOrd<subCellCount;++subCellOrd) {
146 std::vector<unsigned> subCellNodes;
150 bool isSubset = std::includes( nodes.begin(), nodes.end(),
151 subCellNodes.begin(), subCellNodes.end());
153 subCells.insert(std::make_pair(dim,subCellOrd));
167 std::set<std::pair<unsigned,unsigned> > & closure)
171 closure.insert(std::make_pair(0,subCell));
174 closure.insert(std::make_pair(0,cellTopo.getNodeMap(dim,subCell,0)));
175 closure.insert(std::make_pair(0,cellTopo.getNodeMap(dim,subCell,1)));
176 closure.insert(std::make_pair(1,subCell));
180 unsigned cnt = (shards::CellTopology(cellTopo.getCellTopologyData(dim,subCell))).getSubcellCount(dim-1);
181 for(
unsigned i=0;i<cnt;i++) {
182 int edge = mapCellFaceEdge(cellTopo.getCellTopologyData(),subCell,i);
185 closure.insert(std::make_pair(2,subCell));
224 Kokkos::DynRankView<double,PHX::Device> & coords)
const 228 int numCells = cellVertices.dimension(0);
231 Kokkos::DynRankView<double,PHX::Device> localCoords;
235 coords = Kokkos::DynRankView<double,PHX::Device>(
"coords",numCells,localCoords.dimension(0),
getDimension());
238 Intrepid2::CellTools<PHX::Device> cellTools;
239 cellTools.mapToPhysicalFrame(coords,localCoords,cellVertices,
intrepidBasis_->getBaseCellTopology());
virtual int getSubcellCount(int dim) const
virtual const std::vector< int > & getSubcellIndices(int dim, int cellIndex) const
bool supportsInterpolatoryCoordinates() const
Does this field pattern support interpolatory coordinates?
virtual int getDimension() const
static void buildSubcellClosure(const shards::CellTopology &cellTopo, unsigned dim, unsigned subCell, std::set< std::pair< unsigned, unsigned > > &closure)
Teuchos::RCP< Intrepid2::Basis< double, Kokkos::DynRankView< double, PHX::Device > > > intrepidBasis_
virtual void getSubcellClosureIndices(int dim, int cellIndex, std::vector< int > &indices) const
static void getSubcellNodes(const shards::CellTopology &cellTopo, unsigned dim, unsigned subCell, std::vector< unsigned > &nodes)
void getInterpolatoryCoordinates(Kokkos::DynRankView< double, PHX::Device > &coords) const
virtual shards::CellTopology getCellTopology() const
static void findContainedSubcells(const shards::CellTopology &cellTopo, unsigned dim, const std::vector< unsigned > &nodes, std::set< std::pair< unsigned, unsigned > > &subCells)
#define TEUCHOS_ASSERT(assertion_test)