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

« back to all changes in this revision

Viewing changes to nfem/plib/field_element.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
[ 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:
26
26
#include "rheolef/point.h"
27
27
#include "rheolef/field.h"
28
28
#include "rheolef/quadrature.h"
29
 
#include "rheolef/basis_on_lattice.h"
 
29
#include "rheolef/basis_on_pointset.h"
30
30
#include "rheolef/piola.h"
 
31
#include "rheolef/band.h"
31
32
#include <boost/functional.hpp>
32
33
 
33
34
namespace rheolef { 
43
44
// allocators:
44
45
 
45
46
    field_element();
46
 
    ~field_element();
47
47
 
48
48
// accessors
49
49
 
50
50
    template <class Function>
51
 
    void operator() (const geo_element& K, Function f, std::vector<T>& uk) const;
 
51
    void operator() (const geo_element& K, Function f, std::vector<T>& lk) const;
52
52
 
53
53
    // specializations for field as special Function classes (avoid localization octree search)
54
 
    void operator() (const geo_element& K, const field_basic<T,M>& fh, std::vector<T>& uk) const;
 
54
    void operator() (const geo_element& K, const field_basic<T,M>& fh, std::vector<T>& lk) const;
55
55
 
56
 
    void eval (const geo_element& K, const std::vector<T>& fq, std::vector<T>& uk) const;
57
 
    void eval (const geo_element& K, const std::vector<point_basic<T> >& fq, std::vector<T>& uk) const;
 
56
    void eval (const geo_element& K, const std::vector<T>&               fq, std::vector<T>& lk) const;
 
57
    void eval (const geo_element& K, const std::vector<point_basic<T> >& fq, std::vector<T>& lk) const;
58
58
 
59
59
// modifiers (const, since data are mutable: avoid copies):
60
60
 
62
62
        const quadrature_option_type& qopt
63
63
         = quadrature_option_type(quadrature_option_type::gauss)) const;
64
64
 
 
65
    void set_band (const band_basic<T,M>& gh) const { _band = gh; _is_on_band = true; }
 
66
 
65
67
// accessors (for building elementary fields):
66
68
 
67
 
    const space_basic<T,M>&   get_space() const;
68
 
    const basis&              get_basis() const;
69
 
    const numbering<T,M>&     get_numbering() const;
70
 
    const basis&              get_piola_basis () const;
71
 
 
72
 
    size_type                 dimension() const;
 
69
    const space_basic<T,M>&   get_space() const { return _X; }
 
70
    const numbering<T,M>&     get_numbering() const { return _X.get_numbering(); }
 
71
    const basis_basic<T>&     get_basis() const { return get_numbering().get_basis(); }
 
72
    bool                      is_on_band() const { return _is_on_band; }
 
73
    const band_basic<T,M>&    get_band () const { return _band; }
 
74
    const geo_basic<T,M>&     get_geo() const { return _is_on_band ? _band.level_set() : get_space().get_geo(); }
 
75
    const basis_basic<T>&     get_piola_basis () const { return get_geo().get_piola_basis(); }
 
76
    size_type                 dimension() const { return get_geo().dimension(); }
73
77
 
74
78
// data:
75
79
protected:
76
80
    mutable space_basic<T,M>       _X;  // origin dof numbering
77
 
    mutable quadrature             _quad;
78
 
    mutable basis_on_quadrature    _basis_on_quad;
79
 
    mutable basis_on_quadrature    _piola_on_quad;
 
81
    mutable quadrature<T>          _quad;
 
82
    mutable basis_on_pointset<T>   _basis_on_quad;
 
83
    mutable basis_on_pointset<T>   _piola_on_quad;
80
84
    mutable bool                   _initialized;
 
85
    mutable bool                   _is_on_band;
 
86
    mutable band_basic<T,M>        _band;
81
87
};
82
88
// -----------------------------------------------------------------------
83
89
// inlined
84
90
// -----------------------------------------------------------------------
85
 
 
86
91
template<class T, class M>
87
92
inline
88
93
field_element<T,M>::field_element()
89
94
 : _X(),
90
95
   _quad(),
91
96
   _basis_on_quad(),
92
 
   _initialized(false)
93
 
{
94
 
}
95
 
template<class T, class M>
96
 
inline
97
 
field_element<T,M>::~field_element()
98
 
{
99
 
}
100
 
template<class T, class M>
101
 
inline
102
 
const space_basic<T,M>&
103
 
field_element<T,M>::get_space() const
104
 
{
105
 
  return _X;
106
 
}
107
 
template<class T, class M>
108
 
inline
109
 
const numbering<T,M>& 
110
 
field_element<T,M>::get_numbering() const
111
 
{
112
 
  return get_space().get_numbering();
113
 
}
114
 
template<class T, class M>
115
 
inline
116
 
const basis&
117
 
field_element<T,M>::get_basis() const
118
 
{
119
 
  return get_numbering().get_basis();
120
 
}
121
 
template<class T, class M>
122
 
inline
123
 
const basis&
124
 
field_element<T,M>::get_piola_basis () const
125
 
{
126
 
  return get_space().get_geo().get_piola_basis();
127
 
}
128
 
template<class T, class M>
129
 
inline
130
 
typename field_element<T,M>::size_type
131
 
field_element<T,M>::dimension() const
132
 
{
133
 
  return get_space().get_geo().dimension();
134
 
}
 
97
   _initialized(false),
 
98
   _is_on_band(false),
 
99
   _band()
 
100
{
 
101
}
 
102
// evaluate a template class-function f at quadrature points
 
103
// then call a compiled function for the integrations over element K
135
104
template<class T, class M>
136
105
template<class Function>
137
 
inline
138
106
void 
139
107
field_element<T,M>::operator() (
140
108
    const geo_element& K,
141
109
    Function           f,
142
 
    std::vector<T>&    uk) const
 
110
    std::vector<T>&    lk) const
143
111
{
144
112
  typedef typename boost::unary_traits<Function>::result_type result_type;
145
 
  std::vector<result_type> fq (_quad.size(K));
 
113
  const geo_basic<T,M>& omega = get_geo();
 
114
  reference_element hat_K (K.variant());
 
115
  std::vector<result_type> fq (_quad.size(hat_K));
146
116
  std::vector<size_type> dis_inod;
147
 
  get_space().get_geo().dis_inod (K, dis_inod);
148
 
  for (size_type q = 0, nq = _piola_on_quad.size(K); q < nq; q++) {
149
 
    point_basic<T> xq = piola_transformation (get_space().get_geo(), _piola_on_quad, K, dis_inod, q);
 
117
  omega.dis_inod (K, dis_inod);
 
118
  for (size_type q = 0, nq = _piola_on_quad.size(hat_K); q < nq; q++) {
 
119
    point_basic<T> xq = piola_transformation (omega, _piola_on_quad, K, dis_inod, q);
150
120
    fq[q] = f(xq);
151
121
  }
152
 
  eval (K, fq, uk);
 
122
  eval (K, fq, lk);
153
123
}
154
124
 
155
125
}// namespace rheolef