Rheolef  7.2
an efficient C++ finite element environment
 
Loading...
Searching...
No Matches
geo_element.h
Go to the documentation of this file.
1#ifndef _RHEO_GEO_ELEMENT_H
2#define _RHEO_GEO_ELEMENT_H
3//
4// This file is part of Rheolef.
5//
6// Copyright (C) 2000-2009 Pierre Saramito <Pierre.Saramito@imag.fr>
7//
8// Rheolef is free software; you can redistribute it and/or modify
9// it under the terms of the GNU General Public License as published by
10// the Free Software Foundation; either version 2 of the License, or
11// (at your option) any later version.
12//
13// Rheolef is distributed in the hope that it will be useful,
14// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16// GNU General Public License for more details.
17//
18// You should have received a copy of the GNU General Public License
19// along with Rheolef; if not, write to the Free Software
20// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
21//
22// =========================================================================
23// AUTHOR: Pierre.Saramito@imag.fr
24// DATE: 10 oct 2001, initial: 10 jan 1998
25
26namespace rheolef {
58/* not sure that this example still compiles...
59Example
60--------
61
62 geo_element_auto<> K;
63 K.set_name('t') ;
64 cout << "n_vertices: " << K.size() << endl
65 << "n_edges : " << K.n_edges() << endl
66 << "dimension : " << K.dimension() << endl << endl;
67 for(geo_element::size_type i = 0; i < K.size(); i++)
68 K[i] = i*10 ;
69 for(geo_element::size_type i = 0; i < K.n_edges(); i++)
70 K.set_edge(i, i*10+5) ;
71 cout << "vertices: local -> global" << endl;
72 for (geo_element::size_type vloc = 0; vloc < K.size(); vloc++)
73 cout << vloc << "-> " << K[vloc] << endl;
74 cout << endl
75 << "edges: local -> global" << endl;
76 for (geo_element::size_type eloc = 0; eloc < K.n_edges(); eloc++) {
77 geo_element::size_type vloc1 = subgeo_local_vertex(1, eloc, 0);
78 geo_element::size_type vloc2 = subgeo_local_vertex(1, eloc, 1);
79 cout << eloc << "-> " << K.edge(eloc) << endl
80 << "local_vertex_from_edge(" << eloc
81 << ") -> (" << vloc1 << ", " << vloc2 << ")" << endl;
82 }
83*/
84} // namespace rheolef
85
86#include "rheolef/reference_element.h"
87#include "rheolef/geo_element_indirect.h"
88#include "rheolef/heap_allocator.h"
89#include "rheolef/reference_element_face_transformation.h"
90
91#include <boost/serialization/serialization.hpp>
92#include <boost/serialization/base_object.hpp>
93
94namespace rheolef {
95
96// --------------------------------------------------------------------
97// geo_element: abstract class
98// --------------------------------------------------------------------
99template <class A> class geo_element_auto;
100
101// [verbatim_geo_element]
103public:
104
105// typedefs:
106
107 enum {
108 _variant_offset = 0, // i.e. type, as triangle(t) or tetra(T), etc
109 _order_offset = 1, // i.e. k, when Pk curved element
110 _dis_ie_offset = 2, // internal numbering, depend upon partitionand nproc
111 _ios_dis_ie_offset = 3, // i/o numbering, independent of parition and nproc
112 _master_offset = 4, // (d-1)-side has one or two master d-element that contains it
113 _last_offset = 6 // here starts node indexes, face indexes, etc
114 };
115 // Implementation note: _master_offset reserve 2 size_type but is used only for sides,
116 // i.e. tri or quad in 3d mesh, edge in 2d mesh, or point in 1d
117 // => waste a lot of place
118 // it would be better with a polymorphic class
119 // and the geo class would define an array of smart_pointers on this class
120 // so, we could define edge with or without the 2 size_type for master elements
121 // then, hack_array will become obsolete (good thing)
122 // and reference_element could be also polymorphic, avoiding large swich (variant)
123 // in all internal loops. This change of implementation will be considered
124 // in the future.
128 typedef const size_type* const_iterator;
130
133
135 typedef geo_element_indirect::shift_type shift_type; // for 0..3 face shift
142
143// affectation:
144
146 {
147 reset (K.variant(), K.order()); // resize auto, nothing for hack
148 std::copy (K._data_begin(), K._data_begin() + _data_size(), _data_begin());
149 reset (K.variant(), K.order()); // reset order=1 for hack, resize nothing for auto
150 return *this;
151 }
152 virtual ~geo_element() {}
154
155// implicit conversion:
156
157 operator reference_element () const { return reference_element(variant()); }
158
159// accessors & modifiers:
160
162 size_type order() const { return *(_data_begin() + _order_offset); }
163 size_type dis_ie() const { return *(_data_begin() + _dis_ie_offset); }
165 size_type master (bool i) const { return *(_data_begin() + _master_offset + i); }
166
169 char name() const { return reference_element::name (variant()); }
171
174 void set_master (bool i, size_type dis_ie) const {
175 const_iterator p = _data_begin() + _master_offset + i; // mutable member fct
176 *(const_cast<iterator>(p)) = dis_ie;
177 }
178
181 iterator end() { return begin() + size(); }
182 const_iterator end() const { return begin() + size(); }
183 size_type& operator[] (size_type loc_inod) { return *(begin() + loc_inod); }
184 size_type operator[] (size_type loc_inod) const { return *(begin() + loc_inod); }
185 size_type& node (size_type loc_inod) { return operator[] (loc_inod); }
186 size_type node (size_type loc_inod) const { return operator[] (loc_inod); }
187
188 iterator begin(size_type node_subgeo_dim) { return begin() + first_inod (node_subgeo_dim); }
189 const_iterator begin(size_type node_subgeo_dim) const { return begin() + first_inod (node_subgeo_dim); }
190 iterator end (size_type node_subgeo_dim) { return begin() + last_inod (node_subgeo_dim); }
191 const_iterator end (size_type node_subgeo_dim) const { return begin() + last_inod (node_subgeo_dim); }
192
195 return *(reinterpret_cast<const geo_element_indirect*>(p));
196 }
199 return *(reinterpret_cast<const geo_element_indirect*>(p));
200 }
203 return *(reinterpret_cast<geo_element_indirect*>(p));
204 }
207 return *(reinterpret_cast<geo_element_indirect*>(p));
208 }
209 size_type edge (size_type i) const { return (dimension() <= 1) ? dis_ie() : edge_indirect(i).index(); }
210 size_type face (size_type i) const { return (dimension() <= 2) ? dis_ie() : face_indirect(i).index(); }
211
212 size_type n_subgeo (size_type subgeo_dim) const {
213 return reference_element::n_subgeo (variant(), subgeo_dim); }
215 return (subgeo_dim == 0) ? operator[](i) : (subgeo_dim == 1) ? edge(i) : (subgeo_dim == 2) ? face(i) : dis_ie(); }
216
217 size_type subgeo_n_node (size_type subgeo_dim, size_type loc_isid) const {
218 return reference_element::subgeo_n_node (variant(), order(), subgeo_dim, loc_isid); }
219 size_type subgeo_local_node (size_type subgeo_dim, size_type loc_isid, size_type loc_jsidnod) const {
220 return reference_element::subgeo_local_node (variant(), order(), subgeo_dim, loc_isid, loc_jsidnod); }
221 size_type subgeo_size (size_type subgeo_dim, size_type loc_isid) const {
222 return reference_element::subgeo_n_node (variant(), 1, subgeo_dim, loc_isid); }
223 size_type subgeo_local_vertex(size_type subgeo_dim, size_type i_subgeo, size_type i_subgeo_vertex) const {
224 return reference_element::subgeo_local_node (variant(), 1, subgeo_dim, i_subgeo, i_subgeo_vertex); }
225 size_type first_inod (size_type subgeo_dim) const {
226 return reference_element::first_inod (variant(), order(), subgeo_dim); }
227 size_type last_inod (size_type subgeo_dim) const {
228 return reference_element::last_inod (variant(), order(), subgeo_dim); }
229
230 size_type n_edge () const { return n_subgeo (1); }
231 size_type n_face () const { return n_subgeo (2); }
232
233// orientation accessors:
234
235 // search S in all sides of K
237 const geo_element& S,
238 size_type& loc_isid,
239 size_type& shift) const;
240
242 const geo_element& S,
243 side_information_type& sid) const;
244
246
247 // compare two sides: S and *this
249 orientation_type& orient, shift_type& shift) const;
252 size_type dis_iv0, size_type dis_iv1, size_type dis_iv2,
253 orientation_type& orient,
254 shift_type& shift) const;
256 size_type dis_iv0, size_type dis_iv1, size_type dis_iv2, size_type dis_iv3,
257 orientation_type& orient,
258 shift_type& shift) const;
259
260// i/o;
261
262 void put (std::ostream& is) const;
263 void get (std::istream& os);
264// [verbatim_geo_element]
265
266// internals:
267// static: fix orientation & shift helpers, for 2d edges & 3d faces:
268
270 const geo_element& K,
271 size_type loc_iedg,
272 size_type loc_iedg_j,
274
276 orientation_type orient,
278 size_type loc_iedg_j);
279
280 static void loc_tri_inod2lattice (
281 size_type loc_tri_inod,
283 point_basic<size_type>& ij_lattice);
284
285 static void loc_qua_inod2lattice (
286 size_type loc_qua_inod,
288 point_basic<size_type>& ij_lattice);
289
291 const geo_element& K,
292 size_type loc_itri,
293 size_type loc_itri_j,
295
297 orientation_type orient,
298 shift_type shift,
300 size_type loc_itri_j);
301
303 const geo_element& K,
304 size_type loc_iqua,
305 size_type loc_iqua_j,
307
309 orientation_type orient,
310 shift_type shift,
312 size_type loc_iqua_j);
313
314 // main call: fix on any element type by switch:
315 static size_type fix_indirect (
316 const geo_element& K,
317 size_type subgeo_variant,
318 size_type loc_ige,
319 size_type loc_comp_idof_on_subgeo,
321
322//protected:
323
328 static size_type _data_size (const parameter_type& p) { return _data_size (p.variant,p.order); }
329
330 size_type _data_size() const { return _data_size (variant(),order()); }
331
332 virtual iterator _data_begin() = 0;
333 virtual iterator _data_end() = 0;
334 virtual const_iterator _data_begin() const = 0;
335 virtual const_iterator _data_end() const = 0;
336
337 template<class Archive>
338 void serialize (Archive& ar, const unsigned int version) {
339 }
340#ifdef TO_CLEAN
341 template<class A = std::allocator<std::vector<int>::size_type> > class geo_element_auto;
342#endif // TO_CLEAN
343// [verbatim_geo_element_cont]
344};
345// [verbatim_geo_element_cont]
346
347#ifdef TO_CLEAN
348// specialization, from disarray.h (because of a g++ bug when T=float128)
349template<>
350struct _disarray_put_element_type<geo_element> {
351 std::ostream& operator() (std::ostream& os, const geo_element& K) const { K.put(os); return os; }
352};
353template<>
354struct _disarray_get_element_type<geo_element> {
355 std::istream& operator() (std::istream& is, geo_element& K) const { K.get(is); return is; }
356};
357#endif // TO_CLEAN
358inline
359std::istream&
360operator>> (std::istream& is, geo_element& K)
361{
362 K.get (is);
363 return is;
364}
365inline
366std::ostream&
367operator<< (std::ostream& os, const geo_element& K)
368{
369 K.put (os);
370 return os;
371}
372// --------------------------------------------------------------------
373// geo_element_auto: generic dynamically allocated class
374// --------------------------------------------------------------------
375template<class A = std::allocator<std::vector<int>::size_type> >
377public:
378
379// typedefs:
380
389
390// allocators:
391
392 explicit geo_element_auto (const A& alloc = A())
393 : _data (_last_offset, std::numeric_limits<size_type>::max(), alloc)
394 {
396 _data [_order_offset] = 0;
397 }
398 explicit geo_element_auto (variant_type variant, size_type order = 1, const A& alloc = A())
399 : _data (_data_size(variant,order), std::numeric_limits<size_type>::max(), alloc)
400 {
403 }
404 explicit geo_element_auto (parameter_type p, const A& alloc = A())
405 : _data (_data_size(p), std::numeric_limits<size_type>::max(), alloc)
406 {
407 _data [_variant_offset] = p.variant;
408 _data [_order_offset] = p.order;
409 }
411 : _data (K._data_size(), size_type(0), A()) // cree un nouvel allocateur
412 { std::copy (K._data_begin(), K._data_end(), _data.begin()); }
414 : _data (K._data.size(), size_type(0), K._data.get_allocator()) // re-utilise l'allocateur precedent
415 { std::copy (K._data.begin(), K._data.end(), _data.begin()); }
416 template <class A2>
418 : _data (K._data.size(), size_type(0), A()) // cree un nouvel allocateur
419 { std::copy (K._data.begin(), K._data.end(), _data.begin()); }
421 {
422 _data.resize(K._data_size());
423 std::copy (K._data_begin(), K._data_end(), _data.begin());
424 return *this;
425 }
427 _data.resize (_data_size(variant,order), std::numeric_limits<size_type>::max());
430 }
431 void reset (const parameter_type& param) { reset (param.variant, param.order); }
432
433// internals:
434
435 template<class Archive>
436 void serialize (Archive& ar, const unsigned int version) {
437 ar & boost::serialization::base_object<geo_element>(*this);
438 ar & _data;
439 }
440
441//protected:
442
443 iterator _data_begin() { return _data.begin().operator->(); }
444 const_iterator _data_begin() const { return _data.begin().operator->(); }
445 iterator _data_end() { return _data.end().operator->(); }
446 const_iterator _data_end() const { return _data.end().operator->(); }
447
448 template <class A2> friend class geo_element_auto;
449
450// data:
451
452 std::vector<size_type,A> _data;
453};
454// -------------------------------------------------------------------
455// base raw class
456// --------------------------------------------------------------------
458public:
459
460// constants & typedefs:
461
462 enum {
463 _vtable_size = 1 /* = sizeof(geo_element_X_hack)/sizeof(size_type) */
464 };
465
474
475// allocators:
476
478 template <class A>
480 {
481 check_macro (K.variant() == variant(), "incompatible conversion");
482#ifdef TO_CLEAN
483 check_macro (K.order() == order(), "incompatible conversion");
484#endif // TO_CLEAN
485 std::copy (K._data.begin(), K._data.begin() + _data_size(), _data_begin());
486 }
487
488// accesors & modifiers
489
490 void reset (variant_type variant1, size_type order1) {
491 check_macro (variant1 == variant(), "cannot change variant from "<<variant()<<" to "<<variant1<<" in a raw element");
492#ifdef TO_CLEAN
493 check_macro (order1 == 1, "cannot change order "<<order1<< " > 1 in a raw element");
494#endif // TO_CLEAN
496 }
497// internals:
498protected:
499
501
502 iterator _data_begin() { return reinterpret_cast<iterator> (this) + _vtable_size; }
503 const_iterator _data_begin() const { return reinterpret_cast<const_iterator>(this) + _vtable_size; }
506
507 size_type _get_data (size_type i) const { return *(_data_begin() + i); }
509 void _set_data (size_type i, size_type value) { _get_data_ref(i) = value; }
510
512 check_macro (order == 1, "cannot set order "<<order<< " > 1 in a raw element");
515 for (size_type i = _order_offset+1, n = _data_size (variant,order); i < n; i++) {
516 _set_data (i, std::numeric_limits<size_type>::max());
517 }
518 }
519 void _set_parameter (const parameter_type& p) { _reset (p.variant, p.order); }
520
521 template <class T, class A> friend class hack_array_seq_rep;
522};
523// -------------------------------------------------------------------
524// permuted io helper
525// --------------------------------------------------------------------
528 geo_element_permuted_put (const std::vector<size_type>& perm1) : perm(perm1) {}
529 std::ostream& operator() (std::ostream& os, const geo_element& K) {
530 static const bool do_verbose = true;
531 if (do_verbose || K.size() > 2 || K.order() > 1) { os << K.name() << "\t"; }
532 if (K.order() > 1) { os << "p" << K.order() << " "; }
533 for (geo_element::size_type iloc = 0; iloc < K.n_node(); iloc++) {
534 os << perm [K[iloc]];
535 if (iloc < K.n_node() - 1) os << " ";
536 }
537 return os;
538 }
539 const std::vector<size_type>& perm;
540};
541
542}// namespace rheolef
543#endif // _RHEO_GEO_ELEMENT_H
field::size_type size_type
Definition branch.cc:430
see the edge page for the full documentation
geo_element_auto(const geo_element_auto< A2 > &K)
geo_element_auto(const A &alloc=A())
geo_element_auto(variant_type variant, size_type order=1, const A &alloc=A())
geo_element::automatic_type automatic_type
geo_element::parameter_type parameter_type
geo_element::iterator iterator
std::vector< size_type, A > _data
geo_element::const_iterator const_iterator
geo_element_auto(parameter_type p, const A &alloc=A())
void reset(const parameter_type &param)
void serialize(Archive &ar, const unsigned int version)
reference_element::size_type size_type
const geo_element_auto< A > & operator=(const geo_element &K)
geo_element_auto(const geo_element_auto< A > &K)
const_iterator _data_end() const
void reset(variant_type variant, size_type order)
geo_element_auto(const geo_element &K)
reference_element::variant_type variant_type
const_iterator _data_begin() const
geo_element_hack(const geo_element_auto< A > &K)
void reset(variant_type variant1, size_type order1)
void _reset(variant_type variant, size_type order)
geo_element::automatic_type automatic_type
geo_element::parameter_type parameter_type
geo_element::iterator iterator
geo_element::const_iterator const_iterator
size_type & _get_data_ref(size_type i)
const_iterator _data_end() const
size_type _get_data(size_type i) const
void _set_parameter(const parameter_type &p)
static size_type _size_of(const parameter_type &p)
void _set_data(size_type i, size_type value)
geo_element::variant_type variant_type
geo_element::size_type size_type
const_iterator _data_begin() const
see the geo_element page for the full documentation
geo_element_auto< heap_allocator< size_type > > automatic_type
size_type subgeo_n_node(size_type subgeo_dim, size_type loc_isid) const
geo_element_indirect::orientation_type orientation_type
iterator end(size_type node_subgeo_dim)
const_iterator end(size_type node_subgeo_dim) const
bool get_orientation_and_shift(const geo_element &S, orientation_type &orient, shift_type &shift) const
return orientation and shift between *this element and S
virtual iterator _data_begin()=0
size_type & operator[](size_type loc_inod)
const_iterator begin() const
size_type edge(size_type i) const
size_type face(size_type i) const
size_type n_node() const
size_type subgeo_local_vertex(size_type subgeo_dim, size_type i_subgeo, size_type i_subgeo_vertex) const
size_type first_inod(size_type subgeo_dim) const
static size_type _data_size(variant_type variant, size_type order)
static size_type fix_edge_indirect(const geo_element &K, size_type loc_iedg, size_type loc_iedg_j, size_type order)
const_iterator begin(size_type node_subgeo_dim) const
size_type size() const
iterator begin(size_type node_subgeo_dim)
static size_type fix_triangle_indirect(const geo_element &K, size_type loc_itri, size_type loc_itri_j, size_type order)
void serialize(Archive &ar, const unsigned int version)
reference_element::size_type size_type
size_type n_edge() const
void set_ios_dis_ie(size_type ios_dis_ie)
void put(std::ostream &is) const
size_type node(size_type loc_inod) const
size_type dimension() const
size_type ios_dis_ie() const
static size_type fix_edge_indirect(orientation_type orient, size_type order, size_type loc_iedg_j)
static size_type _face_offset(variant_type variant, size_type order)
void get(std::istream &os)
const size_type * const_iterator
size_type subgeo_size(size_type subgeo_dim, size_type loc_isid) const
static size_type _data_size(const parameter_type &p)
const geo_element_indirect & edge_indirect(size_type i) const
orientation_type get_side_orientation(const geo_element &S) const
geo_element_indirect & edge_indirect(size_type i)
static size_type fix_indirect(const geo_element &K, size_type subgeo_variant, size_type loc_ige, size_type loc_comp_idof_on_subgeo, size_type order)
size_type master(bool i) const
size_type & node(size_type loc_inod)
variant_type variant() const
geo_element_indirect & face_indirect(size_type i)
static void loc_qua_inod2lattice(size_type loc_qua_inod, size_type order, point_basic< size_type > &ij_lattice)
virtual iterator _data_end()=0
geo_element & operator=(const geo_element &K)
virtual const_iterator _data_end() const =0
const_iterator end() const
orientation_type get_edge_orientation(size_type dis_iv0, size_type dis_iv1) const
static size_type _node_offset(variant_type variant, size_type order)
virtual const_iterator _data_begin() const =0
static void loc_tri_inod2lattice(size_type loc_tri_inod, size_type order, point_basic< size_type > &ij_lattice)
size_type _data_size() const
geo_element_indirect::shift_type shift_type
orientation_type get_side_informations(const geo_element &S, size_type &loc_isid, size_type &shift) const
size_type n_subgeo(size_type subgeo_dim) const
geo_element generic_type
size_type dis_ie() const
void set_master(bool i, size_type dis_ie) const
static size_type fix_quadrangle_indirect(const geo_element &K, size_type loc_iqua, size_type loc_iqua_j, size_type order)
size_type n_face() const
const geo_element_indirect & face_indirect(size_type i) const
virtual void reset(variant_type variant, size_type order)=0
reference_element::variant_type variant_type
size_type subgeo_local_node(size_type subgeo_dim, size_type loc_isid, size_type loc_jsidnod) const
size_type last_inod(size_type subgeo_dim) const
size_type order() const
size_type subgeo_dis_index(size_type subgeo_dim, size_type i) const
void set_dis_ie(size_type dis_ie)
static size_type _edge_offset(variant_type variant, size_type order)
see the reference_element page for the full documentation
static size_type subgeo_local_node(variant_type variant, size_type order, size_type subgeo_dim, size_type loc_isid, size_type loc_jsidnod)
static size_type first_inod(variant_type variant, size_type order, size_type subgeo_dim)
static const variant_type max_variant
static size_type n_sub_edge(variant_type variant)
static size_type n_sub_face(variant_type variant)
static size_type last_inod(variant_type variant, size_type order, size_type subgeo_dim)
static size_type subgeo_n_node(variant_type variant, size_type order, size_type subgeo_dim, size_type loc_isid)
std::vector< int >::size_type size_type
static size_type n_node(variant_type variant, size_type order)
size_type n_subgeo(size_type subgeo_dim) const
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
const size_t face[n_face][4]
This file is part of Rheolef.
std::ostream & operator<<(std::ostream &os, const catchmark &m)
Definition catchmark.h:99
std::istream & operator>>(std::istream &is, const catchmark &m)
Definition catchmark.h:88
STL namespace.
Definition sphere.icc:25
std::istream & operator()(std::istream &is, T &x) const
Definition disarray.h:206
std::ostream & operator()(std::ostream &os, const T &x) const
Definition disarray.h:197
parameter_type(variant_type v=reference_element::max_variant, size_type o=0)
geo_element_permuted_put(const std::vector< size_type > &perm1)
std::ostream & operator()(std::ostream &os, const geo_element &K)
const std::vector< size_type > & perm
geo_element::size_type size_type