46 #include "Teuchos_ArrayRCP.hpp" 47 #include "Teuchos_Assert.hpp" 48 #include "Phalanx_DataLayout_MDALayout.hpp" 49 #include "Intrepid2_DefaultCubatureFactory.hpp" 50 #include "Intrepid2_CubatureControlVolume.hpp" 51 #include "Intrepid2_CubatureControlVolumeSide.hpp" 52 #include "Intrepid2_CubatureControlVolumeBoundary.hpp" 66 setup(in_cubature_degree,cell_data);
76 if(in_cv_type ==
"volume"){
82 }
else if(in_cv_type ==
"side"){
88 }
else if(in_cv_type ==
"boundary"){
95 TEUCHOS_ASSERT(
false);
102 const Teuchos::RCP<const shards::CellTopology> & cell_topology,
125 TEUCHOS_ASSERT(num_faces!=-1);
134 TEUCHOS_ASSERT(
false);
141 cubature_degree = in_cubature_degree;
145 std::stringstream ss;
146 ss <<
"CubaturePoints (Degree=" << cubature_degree;
151 if(cell_data.
isSide() && spatialDimension==1) {
159 Teuchos::RCP<shards::CellTopology> sideTopo = getSideTopology(cell_data);
161 Intrepid2::DefaultCubatureFactory cubature_factory;
162 Teuchos::RCP<Intrepid2::Cubature<PHX::Device::execution_space,double,double>> intrepid_cubature;
165 if (Teuchos::is_null(sideTopo)) {
167 intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(topo, cubature_degree);
171 intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(*sideTopo, cubature_degree);
180 const int cell_dim = cell_topology->getDimension();
181 const int subcell_dim = cell_dim-1;
182 const int num_faces_per_cell = cell_topology->getSubcellCount(subcell_dim);
186 std::string point_rule_name;
188 std::stringstream ss;
189 ss <<
"CubaturePoints (Degree=" << getOrder() <<
",surface)";
190 point_rule_name = ss.str();
195 const int num_points_per_cell = num_faces_per_cell;
196 const int num_points_per_face = 1;
197 PointRule::setup(point_rule_name, num_cells, num_points_per_cell, num_faces, num_points_per_face, cell_topology);
198 _point_offsets.resize(3,0);
199 _point_offsets[0] = 0;
200 _point_offsets[1] = num_points_per_face;
201 _point_offsets[2] = _point_offsets[1]+num_points_per_face;
205 Intrepid2::DefaultCubatureFactory cubature_factory;
207 _point_offsets.resize(num_faces_per_cell+1,0);
208 int test_face_size = -1;
209 for(
int subcell_index=0; subcell_index<num_faces_per_cell; ++subcell_index){
210 Teuchos::RCP<shards::CellTopology> face_topology = Teuchos::rcp(
new shards::CellTopology(cell_topology->getCellTopologyData(subcell_dim,subcell_index)));
211 const auto & intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(*face_topology, getOrder());
212 const int num_face_points = intrepid_cubature->getNumPoints();
213 _point_offsets[subcell_index+1] = _point_offsets[subcell_index] + num_face_points;
216 if(test_face_size==-1){
217 test_face_size = num_face_points;
219 TEUCHOS_ASSERT(num_face_points == test_face_size);
223 const int num_points_per_cell = _point_offsets.back();
224 const int num_points_per_face = _point_offsets[1];
226 PointRule::setup(point_rule_name, num_cells, num_points_per_cell, num_faces, num_points_per_face, cell_topology);
234 cv_type = in_cv_type;
235 if (cv_type ==
"volume") {
236 cubature_degree = 75;
238 if (cv_type ==
"side") {
239 cubature_degree = 85;
241 if (cv_type ==
"boundary") {
242 cubature_degree = 95;
247 std::stringstream ss;
248 ss <<
"CubaturePoints ControlVol (Index=" << cubature_degree;
252 Teuchos::RCP<Intrepid2::Cubature<PHX::Device::execution_space,double,double> > intrepid_cubature;
254 int tmp_num_points = 0;
255 if (cv_type ==
"volume") {
257 intrepid_cubature = Teuchos::rcp(
new Intrepid2::CubatureControlVolume<PHX::Device::execution_space,double,double>(topo));
258 tmp_num_points = intrepid_cubature->getNumPoints();
260 else if (cv_type ==
"side") {
262 intrepid_cubature = Teuchos::rcp(
new Intrepid2::CubatureControlVolumeSide<PHX::Device::execution_space,double,double>(topo));
263 tmp_num_points = intrepid_cubature->getNumPoints();
265 else if (cv_type ==
"boundary") {
267 intrepid_cubature = Teuchos::rcp(
new Intrepid2::CubatureControlVolumeBoundary<PHX::Device::execution_space,double,double>(topo,cell_data.
side()));
268 tmp_num_points = intrepid_cubature->getNumPoints();
275 {
return cubature_degree; }
281 TEUCHOS_ASSERT(getType() == SURFACE);
282 return _point_offsets[subcell_index];
288 os <<
"IntegrationRule ( " 289 <<
"Name = " << getName()
290 <<
", Degree = " << cubature_degree
291 <<
", Dimension = " << spatial_dimension
292 <<
", Workset Size = " << workset_size
293 <<
", Num Points = " << num_points
294 <<
", Side = " << side
301 Teuchos::RCP< Intrepid2::Cubature<PHX::Device::execution_space,double,double> > intrepid_cubature;
302 Intrepid2::DefaultCubatureFactory cubature_factory;
305 intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(*(topology),cubature_degree);
307 intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(*(side_topology),cubature_degree);
309 int num_ip = intrepid_cubature->getNumPoints();
310 Kokkos::DynRankView<double,PHX::Device> cub_weights(
"cub_weights",num_ip);
314 cub_points = Kokkos::DynRankView<double,PHX::Device>(
"cub_points", num_ip, topology->getDimension());
315 intrepid_cubature->getCubature(cub_points, cub_weights);
318 cub_points = Kokkos::DynRankView<double,PHX::Device>(
"cub_points", num_ip, side_topology->getDimension());
319 intrepid_cubature->getCubature(cub_points, cub_weights);
Control volume side integral.
const int & getSide() const
Get side associated with integration - this is for backward compatibility.
const int & getType() const
Get type of integrator.
void setup(const std::string &ptName, int np, const panzer::CellData &cell_data)
Integral over a specific side of cells (side must be set)
No integral specified - default state.
void referenceCoordinates(Kokkos::DynRankView< double, PHX::Device > &container)
Construct an array containing the reference coordinates.
int baseCellDimension() const
Dimension of the base cell. NOT the dimension of the local side, even if the side() method returns tr...
Teuchos::RCP< const shards::CellTopology > getCellTopology() const
Get CellTopology for the base cell.
void setup(int cubature_degree, const panzer::CellData &cell_data)
void setup(const int cubature_order, const int integration_type, const int side=-1)
Setup function.
Data for determining cell topology and dimensionality.
int getPointOffset(const int subcell_index) const
Returns the integration point offset for a given subcell_index (i.e. local face index) ...
void setup_cv(const panzer::CellData &cell_data, std::string cv_type)
void setup_surface(const Teuchos::RCP< const shards::CellTopology > &cell_topology, const int num_cells, const int num_faces)
Setup a surface integration.
Integral over all sides of cells (closed surface integral)
virtual void print(std::ostream &os)
print information about the integration rule
int order() const
Returns the order of integration (cubature degree in intrepid lingo)
IntegrationRule(int cubature_degree, const panzer::CellData &cell_data)
if side = -1 then we use the cell volume integration rule.
const int & getOrder() const
Get order of integrator.