60 const ArrayPoint& cellVertices,
62 : degree_(degree), cubDimension_(2), cellTopology_(cellTopology), cellVertices_(cellVertices){
65 ">>> ERROR (CubaturePolygon): No direct cubature rule implemented for the desired polynomial degree.");
68 std::vector<Scalar> centroid(2,0);
71 for (
int i=0;i<numNodes;i++){
79 centroid[0] /= (6*area);
80 centroid[1] /= (6*area);
83 CubatureDirectTriDefault<Scalar,ArrayPoint,ArrayWeight> cubatureTri(
degree_);
84 int numCubPointsPerTri = cubatureTri.getNumPoints();
85 int cubDim = cubatureTri.getDimension();
87 FieldContainer<Scalar> cubatureTriPoints(numCubPointsPerTri,cubDim);
89 FieldContainer<Scalar> cubatureTriWeights(numCubPointsPerTri);
90 cubatureTri.getCubature(cubatureTriPoints,cubatureTriWeights);
94 FieldContainer<Scalar> cubatureCellPoints(numCells,numCubPointsPerTri,cubDim);
95 for (
int k=0;k<numCells;k++){
96 for (
int i=0;i<numCubPointsPerTri;i++){
97 for (
int j=0;j<cubDim;j++){
98 cubatureCellPoints(k,i,j) = cubatureTriPoints(i,j);
105 shards::CellTopology triangleTopology(shards::getCellTopologyData<shards::Triangle<3> >());
106 int totalCubPoints = numCubPointsPerTri*
cellTopology_.getEdgeCount();
111 FieldContainer<Scalar> physicalPoints(numCells,numCubPointsPerTri,cubDim);
112 FieldContainer<Scalar> trianglePoints(numCells,3,cubDim);
114 for (
int i=0;i<numCells;i++){
115 for (
int j=0;j<cubDim;j++){
118 trianglePoints(i,2,j) = centroid[j];
125 FieldContainer<Scalar> jacobians(numCells,numCubPointsPerTri,cubDim,cubDim);
126 FieldContainer<Scalar> detJacobians(numCells, numCubPointsPerTri);
131 for (
int i=0;i<numCells;i++){
132 for (
int j=0;j<numCubPointsPerTri;j++){
133 for (
int k=0;k<cubDim;k++){
144 ArrayWeight& cubWeights)
const{
145 int numCubPoints = numPoints_;
146 int cellDim = cubDimension_;
148 TEUCHOS_TEST_FOR_EXCEPTION ( ( cubPoints.size() < numCubPoints*cellDim || cubWeights.size() < numCubPoints ),
150 ">>> ERROR (CubaturePolygon): Insufficient space allocated for cubature points or weights.");
152 for (
int pointId = 0; pointId < numCubPoints; pointId++){
153 for (
int dim = 0; dim < cellDim; dim++){
154 cubPoints(pointId,dim) = cubaturePoints_(pointId,dim);
156 cubWeights(pointId) = cubatureWeights_(pointId);