Panzer Version of the Day
Loading...
Searching...
No Matches
Panzer_IntrepidFieldPattern.hpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Panzer: A partial differential equation assembly
5// engine for strongly coupled complex multiphysics systems
6// Copyright (2011) Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact Roger P. Pawlowski (rppawlo@sandia.gov) and
39// Eric C. Cyr (eccyr@sandia.gov)
40// ***********************************************************************
41// @HEADER
42
43#ifndef __Panzer_IntrepidFieldPattern_hpp__
44#define __Panzer_IntrepidFieldPattern_hpp__
45
47
48// Trilinos includes
49#include "Kokkos_Core.hpp"
50#include "Kokkos_DynRankView.hpp"
51#include "Intrepid2_Basis.hpp"
52
53#include "Phalanx_KokkosDeviceTypes.hpp"
54#include "Teuchos_RCP.hpp"
55#include <set>
56
57namespace panzer {
58
63 public:
64 Intrepid2FieldPattern(const Teuchos::RCP< Intrepid2::Basis<PHX::Device,double,double> > &intrepidBasis);
65
66 virtual int getSubcellCount(int dim) const;
67 virtual const std::vector<int> & getSubcellIndices(int dim, int cellIndex) const;
68 virtual int getDimension() const;
69 virtual shards::CellTopology getCellTopology() const;
70
71 virtual void getSubcellClosureIndices(int dim,int cellIndex,std::vector<int> & indices) const;
72
73 // static functions for examining shards objects
74
86 static void buildSubcellClosure(const shards::CellTopology & cellTopo,unsigned dim,unsigned subCell,
87 std::set<std::pair<unsigned,unsigned> > & closure);
88
100 static void findContainedSubcells(const shards::CellTopology & cellTopo,unsigned dim,
101 const std::vector<unsigned> & nodes,
102 std::set<std::pair<unsigned,unsigned> > & subCells);
103
111 static void getSubcellNodes(const shards::CellTopology & cellTopo,unsigned dim,unsigned subCell,
112 std::vector<unsigned> & nodes);
113
122
128 void getInterpolatoryCoordinates(Kokkos::DynRankView<double,PHX::Device> & coords) const;
129
135 void getInterpolatoryCoordinates(const Kokkos::DynRankView<double,PHX::Device> & cellVertices,
136 Kokkos::DynRankView<double,PHX::Device> & coords) const;
137
139 Teuchos::RCP< Intrepid2::Basis<PHX::Device,double,double> > getIntrepidBasis() const;
140
141 protected:
142 Teuchos::RCP< Intrepid2::Basis<PHX::Device,double,double> > intrepidBasis_;
143
144 //mutable std::vector<int> subcellIndices_;
145 mutable std::vector<std::vector<std::vector<int> > > subcellIndicies_;
146 std::vector<int> empty_;
147 };
148
149}
150
151#endif
Teuchos::RCP< Intrepid2::Basis< PHX::Device, double, double > > 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
Teuchos::RCP< Intrepid2::Basis< PHX::Device, double, double > > getIntrepidBasis() const
Returns the underlying Intrepid2::Basis object.
bool supportsInterpolatoryCoordinates() const
Does this field pattern support interpolatory coordinates?
std::vector< std::vector< std::vector< int > > > subcellIndicies_
static void findContainedSubcells(const shards::CellTopology &cellTopo, unsigned dim, const std::vector< unsigned > &nodes, std::set< std::pair< unsigned, unsigned > > &subCells)
virtual const std::vector< int > & getSubcellIndices(int dim, int cellIndex) const
static void buildSubcellClosure(const shards::CellTopology &cellTopo, unsigned dim, unsigned subCell, std::set< std::pair< unsigned, unsigned > > &closure)
virtual int getSubcellCount(int dim) const