dune-fem 2.8.0
Loading...
Searching...
No Matches
gausspoints.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_GAUSSPOINTS_HH
2#define DUNE_FEM_GAUSSPOINTS_HH
3
4#include <cassert>
5#include <vector>
6
7#include <dune/common/visibility.hh>
8
9namespace Dune
10{
11
12 namespace Fem
13 {
14
26 {
27 public:
29 enum { MAXP=10 };
30
32 enum { highestOrder=19 };
33
34 private:
35 std::vector< std::vector< double > > G; //[MAXP+1][MAXP]; // positions of Gauss points
36 std::vector< std::vector< double > > W; //[MAXP+1][MAXP]; // weights associated with points
37 std::vector< int > O; //[MAXP+1]; // order of the rule
38
39 public:
42 GaussPts ();
43
51 double point ( int m, int i ) const
52 {
53 assert(m > 0 && i < m);
54 return G[m][i];
55 }
56
64 double weight ( int m, int i ) const
65 {
66 assert(m > 0 && i < m);
67 return W[m][i];
68 }
69
76 int order ( int m ) const
77 {
78 return O[m];
79 }
80
90 int power ( int y, int d ) const
91 {
92 int m = 1;
93 for( int i = 0; i < d; ++i )
94 m *= y;
95 return m;
96 }
97 };
98
99 } // namespace Fem
100
101} // namespace Dune
102
104#endif // #ifndef DUNE_FEM_GAUSSPOINTS_HH
Definition: bindguard.hh:11
one-dimensional Gauss points and their weights
Definition: gausspoints.hh:26
double point(int m, int i) const
obtain the i-th point of the m-th quadratre
Definition: gausspoints.hh:51
@ highestOrder
Definition: gausspoints.hh:32
int order(int m) const
obtain the order of the m-th quadratre
Definition: gausspoints.hh:76
@ MAXP
Definition: gausspoints.hh:29
int power(int y, int d) const
a simple power method
Definition: gausspoints.hh:90
double weight(int m, int i) const
obtain the i-th weight of the m-th quadratre
Definition: gausspoints.hh:64
GaussPts()
constructor initializing the Gauss points for all orders
Definition: gausspoints_implementation.hh:12