← Back to branch summary

~ubuntu-branches/ubuntu/raring/rheolef/raring-proposed

« back to all changes in this revision

Viewing changes to nfem/geo_element/quadrature.h

  • Committer: Package Import Robot
  • Author(s): Pierre Saramito, Pierre Saramito, Sylvestre Ledru
  • Date: 2012-05-14 14:02:09 UTC
  • mfrom: (1.1.6)
  • Revision ID: package-import@ubuntu.com-20120514140209-dzbdlidkotyflf9e
Tags: 6.1-1
http://bugs.debian.org/671996
[ Pierre Saramito ]
* New upstream release 6.1 (minor changes):
  - support arbitrarily polynomial order Pk
  - source code supports g++-4.7 (closes: #671996)

[ Sylvestre Ledru ]
* update of the watch file

Show diffs side-by-side

added added

removed removed

Lines of Context:
nfem/geo_element/quadrature.h
21
21
/// 
22
22
/// =========================================================================
23
23
 
 
24
#include "rheolef/reference_element.h"
 
25
#include "rheolef/point.h"
 
26
 
 
27
namespace rheolef { 
 
28
 
24
29
/*Class:quadrature
25
30
NAME: @code{quadrature} - quadrature formulae on the reference lement
26
31
@cindex  quadrature formulae
58
63
End:
59
64
*/
60
65
 
61
 
#include "rheolef/reference_element.h"
62
 
#include "rheolef/point.h"
63
 
 
64
 
namespace rheolef { 
65
 
struct weighted_point {
66
 
    weighted_point ();
67
 
    weighted_point (const point& x, Float w);
68
 
    point x;
69
 
    Float w;
70
 
};
71
 
 
72
66
class quadrature_option_type {
73
67
public:
74
 
    typedef enum {
 
68
// typedefs:
 
69
 
 
70
  typedef size_t size_type;
 
71
 
 
72
  typedef enum {
75
73
        gauss           = 0,
76
74
        gauss_lobatto   = 1,
77
75
        gauss_radau     = 2,
78
76
        middle_edge     = 3,
79
77
        superconvergent = 4,
80
78
        max_family      = 5
81
 
    } family_type; // update also family_name[] in quatrature.cc
82
 
 
83
 
    typedef size_t size_type;
84
 
    quadrature_option_type(
 
79
  } family_type; // update also family_name[] in quatrature.cc
 
80
 
 
81
// allocators:
 
82
 
 
83
  quadrature_option_type(
85
84
        family_type ft = quadrature_option_type::gauss,
86
85
        size_type k = std::numeric_limits<size_type>::max());
87
 
    size_t         get_order() const;
88
 
    family_type    get_family() const;
89
 
    std::string    get_family_name() const;
90
 
    void set_order (size_t r);
91
 
    void set_family (family_type ft);
 
86
 
 
87
  quadrature_option_type (const quadrature_option_type& qopt);
 
88
  quadrature_option_type& operator= (const quadrature_option_type& qopt);
 
89
 
 
90
// accessors & modifiers:
 
91
 
 
92
  size_t         get_order() const;
 
93
  family_type    get_family() const;
 
94
  std::string    get_family_name() const;
 
95
  void set_order (size_t r);
 
96
  void set_family (family_type ft);
 
97
// data:
92
98
protected:
93
 
    family_type   _family;
94
 
    size_t        _order;
95
 
};
96
 
class quadrature_on_geo : public std::vector<weighted_point> {
 
99
  family_type   _family;
 
100
  size_t        _order;
 
101
};
 
102
template<class T>
 
103
struct weighted_point {
 
104
  weighted_point () : x(), w(0) {}
 
105
  weighted_point (const point_basic<T>& x1, const T& w1) : x(x1), w(w1) {}
 
106
  point_basic<T> x;
 
107
  T              w;
 
108
};
 
109
template<class T>
 
110
class quadrature_on_geo : public std::vector<weighted_point<T> > {
97
111
public:
98
112
 
99
113
// typedefs:
100
114
 
101
 
    typedef std::vector<weighted_point>:: size_type size_type;
102
 
    typedef point                                   x;
 
115
    typedef std::vector<weighted_point<T> > base;
 
116
    typedef typename base::size_type        size_type;
 
117
    typedef point_basic<T>                  x;
103
118
 
104
119
// alocators/deallocators:
105
120
 
107
122
 
108
123
// modifier:
109
124
 
110
 
    void initialize (reference_element hat_K, quadrature_option_type opt);
111
 
    void init_point (quadrature_option_type opt);
112
 
    void init_edge (quadrature_option_type opt);
113
 
    void init_triangle (quadrature_option_type opt);
114
 
    void init_square (quadrature_option_type opt);
 
125
    void initialize       (reference_element hat_K, quadrature_option_type opt);
 
126
    void init_point       (quadrature_option_type opt);
 
127
    void init_edge        (quadrature_option_type opt);
 
128
    void init_triangle    (quadrature_option_type opt);
 
129
    void init_square      (quadrature_option_type opt);
115
130
    void init_tetrahedron (quadrature_option_type opt);
116
 
    void init_prism (quadrature_option_type opt);
117
 
    void init_hexahedron (quadrature_option_type opt);
118
 
 
119
 
    void wx (const point& x, const Float& w) {
120
 
            push_back (weighted_point(x,w)); }
121
 
 
122
 
    static size_type n_node_gauss (size_t r);
123
 
 
124
 
    friend std::ostream& operator<< (std::ostream&, const quadrature_on_geo&);
 
131
    void init_prism       (quadrature_option_type opt);
 
132
    void init_hexahedron  (quadrature_option_type opt);
 
133
 
 
134
    void wx (const point_basic<T>& x, const T& w) {
 
135
            base::push_back (weighted_point<T>(x,w)); }
 
136
 
 
137
    static size_type n_node_gauss (size_type r);
 
138
    
 
139
    template<class U>
 
140
    friend std::ostream& operator<< (std::ostream&, const quadrature_on_geo<U>&);
125
141
private:
126
142
    quadrature_on_geo (const quadrature_on_geo&);
127
143
    quadrature_on_geo operator= (const quadrature_on_geo&);
128
144
};
129
145
//<quadrature:
 
146
template<class T>
130
147
class quadrature {
131
148
public:
132
149
 
133
150
// typedefs:
134
151
 
135
 
    typedef quadrature_on_geo::size_type size_type;
136
 
    typedef quadrature_option_type::family_type family_type;
137
 
    typedef std::vector<weighted_point>::const_iterator const_iterator;
 
152
    typedef typename quadrature_on_geo<T>::size_type                 size_type;
 
153
    typedef quadrature_option_type::family_type                      family_type;
 
154
    typedef typename std::vector<weighted_point<T> >::const_iterator const_iterator;
138
155
 
139
156
// allocators:
140
157
 
142
159
 
143
160
// modifiers:
144
161
 
145
 
    void set_order (size_type order);
 
162
    void set_order  (size_type order);
146
163
    void set_family (family_type ft);
147
164
 
148
165
// accessors:
153
170
    size_type      size  (reference_element hat_K) const;
154
171
    const_iterator begin (reference_element hat_K) const;
155
172
    const_iterator end   (reference_element hat_K) const;
156
 
    friend std::ostream& operator<< (std::ostream&, const quadrature&);
 
173
    template<class U>
 
174
    friend std::ostream& operator<< (std::ostream&, const quadrature<U>&);
157
175
protected:
158
 
    quadrature_option_type     _options;
159
 
    mutable quadrature_on_geo  _quad        [reference_element::max_variant];
160
 
    mutable std::vector<bool>  _initialized;
 
176
    quadrature_option_type        _options;
 
177
    mutable quadrature_on_geo<T>  _quad [reference_element::max_variant];
 
178
    mutable std::vector<bool>     _initialized;
161
179
    void _initialize (reference_element hat_K) const;
162
180
private:
163
 
    quadrature (const quadrature&);
164
 
    quadrature operator= (const quadrature&);
 
181
    quadrature (const quadrature<T>&);
 
182
    quadrature operator= (const quadrature<T>&);
165
183
};
166
184
//>quadrature:
167
185
// ------------------------------------------------------------
174
192
{
175
193
}
176
194
inline
 
195
quadrature_option_type::quadrature_option_type (const quadrature_option_type& qopt)
 
196
 : _family(qopt._family),
 
197
   _order(qopt._order)
 
198
{
 
199
}
 
200
inline
 
201
quadrature_option_type&
 
202
quadrature_option_type::operator= (const quadrature_option_type& qopt)
 
203
{
 
204
  _family = qopt._family;
 
205
  _order  = qopt._order;
 
206
  return *this;
 
207
}
 
208
inline
177
209
quadrature_option_type::size_type
178
210
quadrature_option_type::get_order () const
179
211
{
197
229
{
198
230
    _family = ft;
199
231
}
200
 
inline
201
 
weighted_point::weighted_point ()
202
 
 : x(), w(0) {}
203
 
 
204
 
inline
205
 
weighted_point::weighted_point (const point& x1, Float w1)
206
 
 : x(x1), w(w1) {}
207
 
 
208
 
inline 
209
 
quadrature_on_geo::quadrature_on_geo ()
210
 
 : std::vector<weighted_point>()
211
 
{
212
 
}
213
 
inline 
214
 
quadrature::size_type
215
 
quadrature::get_order () const
 
232
template<class T>
 
233
inline 
 
234
quadrature_on_geo<T>::quadrature_on_geo ()
 
235
 : std::vector<weighted_point<T> >()
 
236
{
 
237
}
 
238
template<class T>
 
239
inline
 
240
typename quadrature_on_geo<T>::size_type
 
241
quadrature_on_geo<T>::n_node_gauss (size_type r)
 
242
{
 
243
    // when using n nodes : gauss quadrature formulae order is r=2*n-1
 
244
    if (r == 0) return 1;
 
245
    size_type n = (r % 2 == 0) ? r/2+1 : (r+1)/2;
 
246
    return std::max(size_t(1), n);
 
247
}
 
248
template<class T>
 
249
inline 
 
250
typename quadrature<T>::size_type
 
251
quadrature<T>::get_order () const
216
252
{
217
253
    return _options.get_order();
218
254
}
 
255
template<class T>
219
256
inline 
220
 
quadrature::family_type
221
 
quadrature::get_family () const
 
257
typename quadrature<T>::family_type
 
258
quadrature<T>::get_family () const
222
259
{
223
260
    return _options.get_family();
224
261
}
 
262
template<class T>
225
263
inline 
226
264
std::string
227
 
quadrature::get_family_name () const
 
265
quadrature<T>::get_family_name () const
228
266
{
229
267
    return _options.get_family_name();
230
268
}
 
269
template<class T>
231
270
inline 
232
 
quadrature::quadrature (quadrature_option_type opt)
 
271
quadrature<T>::quadrature (quadrature_option_type opt)
233
272
 : _options(opt),
234
273
   _quad(),
235
274
   _initialized(reference_element::max_variant, false)
236
275
{
237
276
}
238
277
 
 
278
template<class T>
239
279
inline
240
280
void
241
 
quadrature::set_order (size_type r)
 
281
quadrature<T>::set_order (size_type r)
242
282
{
243
283
  if (get_order() != r) {
244
284
    // do not re-initialize nodes-weights if unchanged
246
286
    std::fill (_initialized.begin(), _initialized.end(), false);
247
287
  }
248
288
}
 
289
template<class T>
249
290
inline
250
291
void
251
 
quadrature::set_family (family_type ft)
 
292
quadrature<T>::set_family (family_type ft)
252
293
{
253
294
  if (get_family() != ft) {
254
295
    // do not re-initialize nodes-weights if unchanged
256
297
    std::fill (_initialized.begin(), _initialized.end(), false);
257
298
  }
258
299
}
259
 
inline
260
 
quadrature_on_geo::size_type
261
 
quadrature_on_geo::n_node_gauss (size_type r)
262
 
{
263
 
    // when using n nodes : gauss quadrature formulae order is r=2*n-1
264
 
    if (r == 0) return 1;
265
 
    size_type n = (r % 2 == 0) ? r/2+1 : (r+1)/2;
266
 
    return std::max(size_t(1), n);
267
 
}
 
300
 
268
301
}// namespace rheolef
269
302
#endif // _RHEO_QUADRATURE_H