dune-fem  2.4.1-rc
pyramidpoints.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_PYRAMIDPOINTS_HH
2 #define DUNE_FEM_PYRAMIDPOINTS_HH
3 
4 #include <dune/common/fvector.hh>
5 
6 namespace Dune
7 {
8 
9  namespace Fem
10  {
11 
16  {
17  public:
18  enum { numQuads = 2 };
19  enum { MAXP = 8 };
20  enum { highest_order = 2 };
21 
23  inline static const PyramidPoints& instance() {
24  static PyramidPoints pyramidPoints;
25  return pyramidPoints;
26  }
27 
29  const FieldVector<double, 3>& point(int m, int i) const
30  {
31  return G[m][i];
32  }
33 
35  double weight (int m, int i) const
36  {
37  return W[m][i];
38  }
39 
41  int order (int m) const
42  {
43  return O[m];
44  }
45 
47  int numPoints(int m) const
48  {
49  return N[m];
50  }
51 
52  private:
55  {
56  int m = 0;
57  O[m] = 0;
58  N[m] = 0;
59 
60  // polynom degree 2 ???
61  m = 1;
62  G[m][0][0] =0.58541020;
63  G[m][0][1] =0.72819660;
64  G[m][0][2] =0.13819660;
65 
66  G[m][1][0] =0.13819660;
67  G[m][1][1] =0.72819660;
68  G[m][1][2] =0.13819660;
69 
70  G[m][2][0] =0.13819660;
71  G[m][2][1] =0.27630920;
72  G[m][2][2] =0.58541020;
73 
74  G[m][3][0] =0.13819660;
75  G[m][3][1] =0.27630920;
76  G[m][3][2] =0.13819660;
77 
78  G[m][4][0] =0.72819660;
79  G[m][4][1] =0.13819660;
80  G[m][4][2] =0.13819660;
81 
82  G[m][5][0] =0.72819660;
83  G[m][5][1] =0.58541020;
84  G[m][5][2] =0.13819660;
85 
86  G[m][6][0] =0.27630920;
87  G[m][6][1] =0.13819660;
88  G[m][6][2] =0.58541020;
89 
90  G[m][7][0] =0.27630920;
91  G[m][7][1] =0.13819660;
92  G[m][7][2] =0.13819660;
93 
94  W[m][0] = 0.041666667;
95  W[m][1] = 0.041666667;
96  W[m][2] = 0.041666667;
97  W[m][3] = 0.041666667;
98  W[m][4] = 0.041666667;
99  W[m][5] = 0.041666667;
100  W[m][6] = 0.041666667;
101  W[m][7] = 0.041666667;
102 
103  O[m] = 2;// verify ????????
104  N[m] = 8;
105  }
106 
107  private:
108  FieldVector<double, 3> G[numQuads][MAXP];
109  double W[numQuads][MAXP]; // weights associated with points
110  int O[numQuads]; // order of the rule
111  int N[numQuads]; // number of points in quadrature rule
112 
113  };
114 
115  } // namespace Fem
116 
117 } // namespace Dune
118 
119 #endif // #ifndef DUNE_FEM_PYRAMIDPOINTS_HH
double weight(int m, int i) const
Access to the ith weight of quadrature rule m.
Definition: pyramidpoints.hh:35
Definition: pyramidpoints.hh:20
Definition: pyramidpoints.hh:19
int numPoints(int m) const
Number of points in the quadrature rule m.
Definition: pyramidpoints.hh:47
Definition: coordinate.hh:4
Definition: pyramidpoints.hh:18
const FieldVector< double, 3 > & point(int m, int i) const
Access to the ith point of quadrature rule m.
Definition: pyramidpoints.hh:29
int order(int m) const
Actual order of quadrature rule m.
Definition: pyramidpoints.hh:41
static const PyramidPoints & instance()
Access to the singleton object.
Definition: pyramidpoints.hh:23
Definition: pyramidpoints.hh:15