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);
75 if(in_cv_type ==
"volume"){
77 "IntegrationRule::IntegrationRule : Control Volume 'volume' type requested, but CellData is setup for sides.");
79 }
else if(in_cv_type ==
"side"){
81 }
else if(in_cv_type ==
"boundary"){
83 "IntegrationRule::IntegrationRule : Control Volume 'boundary' type requested, but CellData is not setup for sides.");
132 cubature_degree = in_cubature_degree;
136 std::stringstream ss;
137 ss <<
"CubaturePoints (Degree=" << cubature_degree;
142 if(cell_data.
isSide() && spatialDimension==1) {
152 Intrepid2::DefaultCubatureFactory cubature_factory;
158 intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(topo, cubature_degree);
162 intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(*sideTopo, cubature_degree);
171 const int cell_dim = cell_topology->getDimension();
172 const int subcell_dim = cell_dim-1;
173 const int num_faces_per_cell = cell_topology->getSubcellCount(subcell_dim);
177 std::string point_rule_name;
179 std::stringstream ss;
180 ss <<
"CubaturePoints (Degree=" << getOrder() <<
",surface)";
181 point_rule_name = ss.str();
186 const int num_points_per_cell = num_faces_per_cell;
187 const int num_points_per_face = 1;
188 PointRule::setup(point_rule_name, num_cells, num_points_per_cell, num_faces, num_points_per_face, cell_topology);
189 _point_offsets.resize(3,0);
190 _point_offsets[0] = 0;
191 _point_offsets[1] = num_points_per_face;
192 _point_offsets[2] = _point_offsets[1]+num_points_per_face;
196 Intrepid2::DefaultCubatureFactory cubature_factory;
198 _point_offsets.resize(num_faces_per_cell+1,0);
199 int test_face_size = -1;
200 for(
int subcell_index=0; subcell_index<num_faces_per_cell; ++subcell_index){
202 const auto & intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(*face_topology, getOrder());
203 const int num_face_points = intrepid_cubature->getNumPoints();
204 _point_offsets[subcell_index+1] = _point_offsets[subcell_index] + num_face_points;
207 if(test_face_size==-1){
208 test_face_size = num_face_points;
214 const int num_points_per_cell = _point_offsets.back();
215 const int num_points_per_face = _point_offsets[1];
217 PointRule::setup(point_rule_name, num_cells, num_points_per_cell, num_faces, num_points_per_face, cell_topology);
225 cv_type = in_cv_type;
226 if (cv_type ==
"volume") {
227 cubature_degree = 75;
229 if (cv_type ==
"side") {
230 cubature_degree = 85;
232 if (cv_type ==
"boundary") {
233 cubature_degree = 95;
238 std::stringstream ss;
239 ss <<
"CubaturePoints ControlVol (Index=" << cubature_degree;
245 int tmp_num_points = 0;
246 if (cv_type ==
"volume") {
248 intrepid_cubature =
Teuchos::rcp(
new Intrepid2::CubatureControlVolume<PHX::Device::execution_space,double,double>(topo));
249 tmp_num_points = intrepid_cubature->getNumPoints();
251 else if (cv_type ==
"side") {
253 intrepid_cubature =
Teuchos::rcp(
new Intrepid2::CubatureControlVolumeSide<PHX::Device::execution_space,double,double>(topo));
254 tmp_num_points = intrepid_cubature->getNumPoints();
256 else if (cv_type ==
"boundary") {
258 intrepid_cubature =
Teuchos::rcp(
new Intrepid2::CubatureControlVolumeBoundary<PHX::Device::execution_space,double,double>(topo,cell_data.
side()));
259 tmp_num_points = intrepid_cubature->getNumPoints();
266 {
return cubature_degree; }
273 return _point_offsets[subcell_index];
279 os <<
"IntegrationRule ( "
280 <<
"Name = " << getName()
281 <<
", Degree = " << cubature_degree
282 <<
", Dimension = " << spatial_dimension
283 <<
", Workset Size = " << workset_size
284 <<
", Num Points = " << num_points
285 <<
", Side = " << side
293 Intrepid2::DefaultCubatureFactory cubature_factory;
296 intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(*(topology),cubature_degree);
298 intrepid_cubature = cubature_factory.create<PHX::Device::execution_space,double,
double>(*(side_topology),cubature_degree);
300 int num_ip = intrepid_cubature->getNumPoints();
301 Kokkos::DynRankView<double,PHX::Device> cub_weights(
"cub_weights",num_ip);
305 cub_points = Kokkos::DynRankView<double,PHX::Device>(
"cub_points", num_ip, topology->getDimension());
306 intrepid_cubature->getCubature(cub_points, cub_weights);
309 cub_points = Kokkos::DynRankView<double,PHX::Device>(
"cub_points", num_ip, side_topology->getDimension());
310 intrepid_cubature->getCubature(cub_points, cub_weights);
Control volume side integral.
int baseCellDimension() const
Dimension of the base cell. NOT the dimension of the local side, even if the side() method returns tr...
void setup(const std::string &ptName, int np, const panzer::CellData &cell_data)
bool is_null(const std::shared_ptr< T > &p)
Integral over a specific side of cells (side must be set)
Teuchos::RCP< const shards::CellTopology > getCellTopology() const
Get CellTopology for the base cell.
No integral specified - default state.
const int & getType() const
Get type of integrator.
void referenceCoordinates(Kokkos::DynRankView< double, PHX::Device > &container)
Construct an array containing the reference coordinates.
void setup(int cubature_degree, const panzer::CellData &cell_data)
int getPointOffset(const int subcell_index) const
Returns the integration point offset for a given subcell_index (i.e. local face index) ...
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void setup(const int cubature_order, const int integration_type, const int side=-1)
Setup function.
#define TEUCHOS_TEST_FOR_EXCEPT_MSG(throw_exception_test, msg)
Data for determining cell topology and dimensionality.
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.
const int & getSide() const
Get side associated with integration - this is for backward compatibility.
Integral over all sides of cells (closed surface integral)
virtual void print(std::ostream &os)
print information about the integration rule
IntegrationRule(int cubature_degree, const panzer::CellData &cell_data)
if side = -1 then we use the cell volume integration rule.
int order() const
Returns the order of integration (cubature degree in intrepid lingo)
#define TEUCHOS_ASSERT(assertion_test)
const int & getOrder() const
Get order of integrator.