Intrepid
Intrepid_CubatureSparse.hpp
Go to the documentation of this file.
1 // @HEADER
2 // ************************************************************************
3 //
4 // Intrepid Package
5 // Copyright (2007) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Pavel Bochev (pbboche@sandia.gov)
38 // Denis Ridzal (dridzal@sandia.gov), or
39 // Kara Peterson (kjpeter@sandia.gov)
40 //
41 // ************************************************************************
42 // @HEADER
43 
49 #ifndef INTREPID_CUBATURE_SPARSE_HPP
50 #define INTREPID_CUBATURE_SPARSE_HPP
51 
52 #include "Intrepid_ConfigDefs.hpp"
53 #include "Intrepid_Cubature.hpp"
55 #include "Intrepid_CubatureSparseHelper.hpp"
56 #include "Teuchos_Assert.hpp"
57 
58 
63 #define INTREPID_CUBATURE_SPARSE2D_GAUSS_MAX 59
64 
69 #define INTREPID_CUBATURE_SPARSE3D_GAUSS_MAX 57
70 
71 
72 namespace Intrepid{
73 
74 template<class Scalar, int dimension_, class ArrayPoint = FieldContainer<Scalar>, class ArrayWeight = ArrayPoint>
75 class CubatureSparse : public Intrepid::Cubature<Scalar,ArrayPoint,ArrayWeight> {
76  private:
77 
78  int level_;
79 
80  int numPoints_;
81 
82  const int degree_;
83 
84 
85  public:
86 
87  ~CubatureSparse() {}
88 
89 
90  CubatureSparse(const int degree);
91 
98  virtual void getCubature(ArrayPoint & cubPoints,
99  ArrayWeight & cubWeights) const;
100 
108  virtual void getCubature(ArrayPoint& cubPoints,
109  ArrayWeight& cubWeights,
110  ArrayPoint& cellCoords) const;
111 
114  virtual int getNumPoints() const;
115 
118  virtual int getDimension() const;
119 
123  virtual void getAccuracy(std::vector<int> & accuracy) const;
124 
125 }; // end class CubatureSparse
126 
127 // helper functions
128 
129 template<class Scalar, int DIM>
130 void iterateThroughDimensions(int level,
131  int dims_left,
132  SGNodes<Scalar,DIM> & cubPointsND,
133  Teuchos::Array<Scalar> & partial_node,
134  Scalar partial_weight);
135 
136 inline int factorial(int num)
137 {
138  int answer = 1;
139  if(num >= 1)
140  {
141  while(num > 0)
142  {
143  answer = answer*num;
144  num--;
145  }
146  }
147  else if(num == 0)
148  answer = 1;
149  else
150  answer = -1;
151 
152  return answer;
153 }
154 
155 inline double combination(int top, int bot)
156 {
157  double answer = factorial(top)/(factorial(bot) * factorial(top-bot));
158  return answer;
159 }
160 
161 inline int iterateThroughDimensionsForNumCalc(int dims_left,
162  int level,
163  int levels_left,
164  int level_so_far,
165  Teuchos::Array<int> & nodes,
166  int product,
167  bool no_uni_quad)
168 {
169  int numNodes = 0;
170  for(int j = 1; j <= levels_left; j++)
171  {
172  bool temp_bool = no_uni_quad;
173  int temp_knots = nodes[j-1]*product;
174  int temp_lsf = level_so_far + j;
175 
176  if(j==1)
177  temp_bool = false;
178 
179  if(dims_left == 1)
180  {
181  if(temp_lsf < level && temp_bool == true)
182  numNodes += 0;
183  else
184  {
185  numNodes += temp_knots;
186  }
187  }
188  else
189  {
190  numNodes += iterateThroughDimensionsForNumCalc(dims_left-1,level, levels_left-j+1, temp_lsf, nodes, temp_knots, temp_bool);
191  }
192  }
193  return numNodes;
194 }
195 
196 inline int calculateNumPoints(int dim, int level)
197 {
198  //int* uninum = new int[level];
199  Teuchos::Array<int> uninum(level);
200  uninum[0] = 1;
201  for(int i = 1; i <= level-1; i++)
202  {
203  uninum[i] = 2*i;
204  }
205 
206  int numOfNodes = iterateThroughDimensionsForNumCalc(dim, level, level, 0, uninum, 1, true);
207  return numOfNodes;
208 }
209 
210 
211 } // end namespace Intrepid
212 
213 
214 // include templated definitions
216 
217 #endif
Header file for the Intrepid::CubatureDirectLineGauss class.
virtual void getAccuracy(std::vector< int > &accuracy) const
Returns algebraic accuracy (e.g. max. degree of polynomial that is integrated exactly).
virtual int getNumPoints() const
Returns the number of cubature points.
virtual int getDimension() const
Returns dimension of the integration domain.
Definition file for the Intrepid::CubatureSparse class.
virtual void getCubature(ArrayPoint &cubPoints, ArrayWeight &cubWeights) const
Returns cubature points and weights (return arrays must be pre-sized/pre-allocated).
Defines the base class for cubature (integration) rules in Intrepid.
Header file for the Intrepid::Cubature class.