51template <
class Scalar,
class ArrayPo
int,
class ArrayWeight>
53 unsigned numCubs = cubatures.size();
54 TEUCHOS_TEST_FOR_EXCEPTION( (numCubs < 1),
56 ">>> ERROR (CubatureTensor): Input cubature array must be of size 1 or larger.");
58 cubatures_ = cubatures;
60 unsigned numDegrees = 0;
61 for (
unsigned i=0; i<numCubs; i++) {
63 cubatures[i]->getAccuracy(tmp);
64 numDegrees += tmp.size();
67 degree_.assign(numDegrees, 0);
70 for (
unsigned i=0; i<numCubs; i++) {
72 cubatures[i]->getAccuracy(tmp);
73 for (
unsigned j=0; j<tmp.size(); j++) {
74 degree_[count] = tmp[j];
77 dimension_ += cubatures[i]->getDimension();
83template <
class Scalar,
class ArrayPo
int,
class ArrayWeight>
87 cubatures_[0] = cubature1;
88 cubatures_[1] = cubature2;
91 for (
unsigned i=0; i<2; i++){
93 cubatures_[i]->getAccuracy(d); degree_[i] = d[0];
96 dimension_ = cubatures_[0]->getDimension() + cubatures_[1]->getDimension();
101template <
class Scalar,
class ArrayPo
int,
class ArrayWeight>
105 cubatures_.resize(3);
106 cubatures_[0] = cubature1;
107 cubatures_[1] = cubature2;
108 cubatures_[2] = cubature3;
110 degree_.assign(3, 0);
111 for (
unsigned i=0; i<3; i++){
113 cubatures_[i]->getAccuracy(d); degree_[i] = d[0];
116 dimension_ = cubatures_[0]->getDimension() + cubatures_[1]->getDimension() + cubatures_[2]->getDimension();
121template <
class Scalar,
class ArrayPo
int,
class ArrayWeight>
123 cubatures_.resize(n);
124 for (
int i=0; i<n; i++) {
125 cubatures_[i] = cubature;
129 cubatures_[0]->getAccuracy(d);
130 degree_.assign(n,d[0]);
132 dimension_ = cubatures_[0]->getDimension()*n;
137template <
class Scalar,
class ArrayPo
int,
class ArrayWeight>
139 ArrayWeight & cubWeights)
const {
140 int numCubPoints = getNumPoints();
141 int cubDim = getDimension();
143 TEUCHOS_TEST_FOR_EXCEPTION( ( ( (
int)cubPoints.size() < numCubPoints*cubDim ) || ( (
int)cubWeights.size() < numCubPoints ) ),
145 ">>> ERROR (CubatureTensor): Insufficient space allocated for cubature points or weights.");
147 unsigned numCubs = cubatures_.size();
148 std::vector<unsigned> numLocPoints(numCubs);
149 std::vector<unsigned> locDim(numCubs);
150 std::vector< FieldContainer<Scalar> > points(numCubs);
151 std::vector< FieldContainer<Scalar> > weights(numCubs);
154 for (
unsigned i=0; i<numCubs; i++) {
156 numLocPoints[i] = cubatures_[i]->getNumPoints();
157 locDim[i] = cubatures_[i]->getDimension();
158 points[i].resize(numLocPoints[i], locDim[i]);
159 weights[i].resize(numLocPoints[i]);
162 cubatures_[i]->getCubature(cubPoints, cubWeights);
163 for (
unsigned pt=0; pt<numLocPoints[i]; pt++) {
164 for (
unsigned d=0; d<locDim[i]; d++) {
165 points[i](pt,d) = cubPoints(pt,d);
166 weights[i](pt) = cubWeights(pt);
173 for (
int i=0; i<numCubPoints; i++) {
174 cubWeights(i) = (Scalar)1.0;
178 int globDimCounter = 0;
180 for (
unsigned i=0; i<numCubs; i++) {
182 for (
int j=0; j<numCubPoints; j++) {
184 int itmp = (j % (numCubPoints/shift))*shift + (j / (numCubPoints/shift));
185 for (
unsigned k=0; k<locDim[i]; k++) {
186 cubPoints(itmp , globDimCounter+k) = points[i](j % numLocPoints[i], k);
188 cubWeights( itmp ) *= weights[i](j % numLocPoints[i]);
191 shift *= numLocPoints[i];
192 globDimCounter += locDim[i];
197template<
class Scalar,
class ArrayPo
int,
class ArrayWeight>
199 ArrayWeight& cubWeights,
200 ArrayPoint& cellCoords)
const
202 TEUCHOS_TEST_FOR_EXCEPTION( (
true), std::logic_error,
203 ">>> ERROR (CubatureTensor): Cubature defined in reference space calling method for physical space cubature.");
207template <
class Scalar,
class ArrayPo
int,
class ArrayWeight>
209 unsigned numCubs = cubatures_.size();
210 int numCubPoints = 1;
211 for (
unsigned i=0; i<numCubs; i++) {
212 numCubPoints *= cubatures_[i]->getNumPoints();
218template <
class Scalar,
class ArrayPo
int,
class ArrayWeight>
225template <
class Scalar,
class ArrayPo
int,
class ArrayWeight>
Defines direct cubature (integration) rules in Intrepid.
virtual int getNumPoints() const
Returns the number of cubature points.
virtual void getAccuracy(std::vector< int > °ree) const
Returns max. degree of polynomials that are integrated exactly. The return vector has the size of the...
virtual int getDimension() const
Returns dimension of integration domain.
virtual void getCubature(ArrayPoint &cubPoints, ArrayWeight &cubWeights) const
Returns cubature points and weights (return arrays must be pre-sized/pre-allocated).
CubatureTensor(std::vector< Teuchos::RCP< Cubature< Scalar, ArrayPoint, ArrayWeight > > > cubatures)
Constructor.
Defines the base class for cubature (integration) rules in Intrepid.