rheolef  6.3
geo.h
Go to the documentation of this file.
1 #ifndef _RHEOLEF_GEO_H
2 #define _RHEOLEF_GEO_H
3 
4 #include "rheolef/geo_element.h"
5 #include "rheolef/array.h"
6 #include "rheolef/hack_array.h"
7 #include "rheolef/geo_size.h"
8 #include "rheolef/point.h"
9 #include "rheolef/domain_indirect.h"
10 #include "rheolef/geo_header.h"
11 #include "rheolef/numbering.h"
12 #include "rheolef/space_constant.h"
13 #include "rheolef/geo_locate.h"
14 #include "rheolef/geo_trace_ray_boundary.h"
15 
16 namespace rheolef {
17 
18 template <class T>
19 struct _point_get {
20  typedef typename point_basic<T>::size_type size_type;
21  _point_get (size_type d1) : d(d1) {}
22  size_type d;
23  std::istream& operator() (std::istream& is, point_basic<T>& x) { return x.get (is, d); }
24 };
25 template <class T>
26 struct _point_put {
27  typedef typename point_basic<T>::size_type size_type;
28  _point_put (size_type d1) : d(d1) {}
29  size_type d;
30  std::ostream& operator() (std::ostream& os, const point_basic<T>& x) { return x.put (os, d); }
31 };
32 
33 template <class T, class M> class geo_basic;
34 template <class T, class M> class geo_domain_indirect_rep;
35 
36 template <class T> class geo_abstract_base_rep;
37 
38 template<class T, class Ref, class Ptr, class IteratorByVariant>
39 struct geo_iterator {
40  typedef geo_iterator<T,Ref,Ptr,IteratorByVariant> _self;
41  typedef geo_iterator<T,T&,T*,typename hack_array<T>::iterator> _nonconst_iterator;
42 
43 
44  // see std::deque<T>::iterator : TODO: complete to a full random iterator
45  typedef std::random_access_iterator_tag iterator_category;
46  typedef T value_type;
47  typedef Ptr pointer;
48  typedef Ref reference;
49  typedef typename T::size_type size_type;
50  typedef ptrdiff_t difference_type;
51 
52 
53  template<class Geo>
54  geo_iterator (size_type dim, size_type variant, IteratorByVariant iter, Geo& omega);
55  geo_iterator (const _nonconst_iterator& y); // conversion from iter to const_iter
56 
57 
58  reference operator* () const { return *_iter_by_var; }
59  pointer operator->() const { return _iter_by_var.operator->(); }
60 
61  _self& operator++ () {
62  ++_iter_by_var;
63  _reset_to_next_or_last();
64  return *this;
65  }
66  void _reset_to_next_or_last () {
67  while (_variant < _variant_max && _iter_by_var == _last_by_var [_variant]) {
68  ++_variant; // then variant <= variant_max
69  if (_variant < _variant_max) {
70  _iter_by_var = _first_by_var [_variant];
71  }
72  }
73  }
74  _self operator++ (int) { _self tmp = *this; operator++(); return tmp; }
75 
76  bool operator== (const _self& y) const { return _variant == y._variant && _iter_by_var == y._iter_by_var; }
77  bool operator!= (const _self& y) const { return ! operator== (y); }
78 
79 
80  size_type _variant;
81  size_type _variant_max;
82  IteratorByVariant _iter_by_var;
83  IteratorByVariant _first_by_var [reference_element::max_variant];
84  IteratorByVariant _last_by_var [reference_element::max_variant];
85 };
86 template <class T>
87 class geo_abstract_base_rep {
88 public:
89 
90 
91  typedef enum {
92  geo = 0,
93  geo_domain = 1,
94  geo_domain_indirect = 2,
95  max_variant = 3
96  } geo_rep_variant_type;
97 
98  typedef geo_element_hack::size_type size_type;
99  typedef point_basic<T> node_type;
100  typedef reference_element::variant_type variant_type;
101  typedef geo_element& reference;
102  typedef const geo_element& const_reference;
103  typedef typename hack_array<geo_element_hack>::iterator
104  iterator_by_variant;
105  typedef typename hack_array<geo_element_hack>::const_iterator
106  const_iterator_by_variant;
107  typedef geo_iterator<geo_element, geo_element&, geo_element*, iterator_by_variant>
108  iterator;
109  typedef geo_iterator<geo_element, const geo_element&, const geo_element*, const_iterator_by_variant>
110  const_iterator;
111  typedef space_constant::coordinate_type coordinate_type;
112 
113 
114  geo_abstract_base_rep () {}
115  virtual ~geo_abstract_base_rep () {}
116 
117 
118  virtual size_type variant() const = 0;
119  virtual std::string name() const = 0;
120  virtual std::string familyname() const = 0;
121  virtual size_type dimension() const = 0;
122  virtual size_type serial_number() const = 0;
123  virtual size_type map_dimension() const = 0;
124  virtual coordinate_type coordinate_system() const = 0;
125  virtual const basis_basic<T>& get_piola_basis() const = 0;
126  virtual const node_type& xmin() const = 0;
127  virtual const node_type& xmax() const = 0;
128  virtual const geo_size& sizes() const = 0;
129  virtual const geo_size& ios_sizes() const = 0;
130 
131  virtual const distributor& geo_element_ownership (size_type dim) const = 0;
132  virtual const_reference get_geo_element (size_type dim, size_type ige) const = 0;
133  virtual size_type neighbour (size_type ie, size_type loc_isid) const = 0;
134  virtual void neighbour_guard() const = 0;
135 
136  virtual const_iterator_by_variant begin_by_variant (variant_type variant) const = 0;
137  virtual const_iterator_by_variant end_by_variant (variant_type variant) const = 0;
138 
139  const_iterator begin (size_type dim) const;
140  const_iterator end (size_type dim) const;
141 
142  virtual size_type n_node() const = 0;
143  virtual const node_type& node (size_type inod) const = 0;
144  virtual const node_type& dis_node (size_type dis_inod) const = 0;
145  virtual node_type piola (const geo_element& K, const node_type& hat_x) const = 0;
146  virtual void dis_inod (const geo_element& K, std::vector<size_type>& dis_inod) const = 0;
147  virtual size_type dis_inod2dis_iv (size_type dis_inod) const = 0;
148 
149  virtual size_type n_domain_indirect () const = 0;
150  virtual bool have_domain_indirect (const std::string& name) const = 0;
151  virtual void reset_order (size_type order) = 0;
152 
153  virtual size_type seq_locate (
154  const point_basic<T>& x,
155  size_type dis_ie_guest = std::numeric_limits<size_type>::max()) const = 0;
156  virtual size_type dis_locate (
157  const point_basic<T>& x,
158  size_type dis_ie_guest = std::numeric_limits<size_type>::max()) const = 0;
159  virtual size_type seq_trace_move (
160  const point_basic<T>& x,
161  const point_basic<T>& v,
162  point_basic<T>& y) const = 0;
163  virtual size_type dis_trace_move (
164  const point_basic<T>& x,
165  const point_basic<T>& v,
166  point_basic<T>& y) const = 0;
167 
168 // virtual i/o:
169 
170  virtual odiststream& put (odiststream& ops) const = 0;
171  virtual bool check(bool verbose) const = 0;
172 
173 
174  bool operator== (const geo_abstract_base_rep<T>& omega2) const {
175  return name() == omega2.name(); }
176 };
177 template<class T, class Ref, class Ptr, class IteratorByVariant>
178 template<class Geo>
179 geo_iterator<T,Ref,Ptr,IteratorByVariant>::geo_iterator (
180  size_type dim,
181  size_type variant,
182  IteratorByVariant iter,
183  Geo& omega)
184  : _variant (variant),
185  _variant_max (reference_element::last_variant_by_dimension(dim)),
186  _iter_by_var (iter),
187  _first_by_var(),
188  _last_by_var()
189 {
190  for (size_type variant = reference_element::first_variant_by_dimension(dim);
191  variant < reference_element::last_variant_by_dimension(dim); variant++) {
192  _first_by_var [variant] = omega.begin_by_variant (variant);
193  _last_by_var [variant] = omega. end_by_variant (variant);
194  }
195  _reset_to_next_or_last();
196 }
197 template<class T, class Ref, class Ptr, class IteratorByVariant>
198 geo_iterator<T,Ref,Ptr,IteratorByVariant>::geo_iterator (const _nonconst_iterator& y)
199  : _variant (y._variant),
200  _variant_max (y._variant_max),
201  _iter_by_var (y._iter_by_var),
202  _first_by_var(),
203  _last_by_var()
204 {
205  std::copy (y._first_by_var, y._first_by_var + reference_element::max_variant, _first_by_var);
206  std::copy (y. _last_by_var, y. _last_by_var + reference_element::max_variant, _last_by_var);
207 }
208 template <class T, class M>
209 class geo_abstract_rep {};
210 
211 template <class T>
212 class geo_abstract_rep<T,sequential> : public geo_abstract_base_rep<T> {
213 public:
214 
215 
216  typedef geo_abstract_base_rep<T> base;
217  typedef typename base::size_type size_type;
218  typedef typename base::node_type node_type;
219  typedef typename base::variant_type variant_type;
220  typedef typename base::iterator iterator;
221  typedef typename base::const_iterator const_iterator;
222  typedef typename base::iterator_by_variant iterator_by_variant;
223  typedef typename base::const_iterator_by_variant const_iterator_by_variant;
224  typedef typename base::reference reference;
225  typedef typename base::const_reference const_reference;
226 
227 
228  geo_abstract_rep () {}
229  virtual geo_abstract_rep<T,sequential>* clone() const = 0;
230  virtual ~geo_abstract_rep () {}
231 
232 
233  virtual const domain_indirect_basic<sequential>& get_domain_indirect (size_type i) const = 0;
234  virtual const domain_indirect_basic<sequential>& get_domain_indirect (const std::string& name) const = 0;
235  virtual void insert_domain_indirect (const domain_indirect_basic<sequential>& dom) const = 0;
236 
237  virtual const array<node_type,sequential>& get_nodes() const = 0;
238  virtual void set_nodes (const array<node_type,sequential>&) = 0;
239 
240  virtual void locate (
241  const array<point_basic<T>, sequential>& x,
242  array<size_type, sequential>& dis_ie,
243  bool do_check = false) const = 0;
244  virtual void trace_ray_boundary (
245  const array<point_basic<T>,sequential>& x,
246  const array<point_basic<T>,sequential>& v,
247  array<size_type, sequential>& dis_ie,
248  array<point_basic<T>,sequential>& y,
249  bool do_check = false) const = 0;
250  virtual void trace_move (
251  const array<point_basic<T>,sequential>& x,
252  const array<point_basic<T>,sequential>& v,
253  array<size_type, sequential>& dis_ie,
254  array<point_basic<T>,sequential>& y) const = 0;
255 
256 };
257 #ifdef _RHEOLEF_HAVE_MPI
258 template <class T>
259 class geo_abstract_rep<T,distributed> : public geo_abstract_base_rep<T> {
260 public:
261 
262 
263  typedef geo_abstract_base_rep<T> base;
264  typedef typename base::size_type size_type;
265  typedef typename base::node_type node_type;
266  typedef typename base::const_reference const_reference;
267  typedef std::map <size_type, node_type, std::less<size_type>,
268  heap_allocator<std::pair<size_type,node_type> > > node_map_type;
269 
270 
271  geo_abstract_rep () {}
272  virtual geo_abstract_rep<T,distributed>* clone() const = 0;
273  virtual ~geo_abstract_rep () {}
274 
275 
276  virtual distributor geo_element_ios_ownership (size_type dim) const = 0;
277  virtual const_reference dis_get_geo_element (size_type dim, size_type dis_ige) const = 0;
278  virtual size_type ige2ios_dis_ige (size_type dim, size_type ige) const = 0;
279  virtual size_type dis_ige2ios_dis_ige (size_type dim, size_type dis_ige) const = 0;
280  virtual size_type ios_ige2dis_ige (size_type dim, size_type ios_ige) const = 0;
281 
282  virtual const domain_indirect_basic<distributed>& get_domain_indirect (size_type i) const = 0;
283  virtual const domain_indirect_basic<distributed>& get_domain_indirect (const std::string& name) const = 0;
284  virtual void insert_domain_indirect (const domain_indirect_basic<distributed>& dom) const = 0;
285 
286  virtual const array<node_type,distributed>& get_nodes() const = 0;
287  virtual void set_nodes (const array<node_type,distributed>&) = 0;
288 
289  virtual void locate (
290  const array<point_basic<T>,distributed>& x,
291  array<size_type, distributed>& dis_ie,
292  bool do_check = true) const = 0;
293  virtual void trace_ray_boundary (
294  const array<point_basic<T>,distributed>& x,
295  const array<point_basic<T>,distributed>& v,
296  array<size_type, distributed>& dis_ie,
297  array<point_basic<T>,distributed>& y,
298  bool do_check = false) const = 0;
299  virtual void trace_move (
300  const array<point_basic<T>,distributed>& x,
301  const array<point_basic<T>,distributed>& v,
302  array<size_type, distributed>& dis_ie,
303  array<point_basic<T>,distributed>& y) const = 0;
304 
305 
306  virtual void set_ios_permutation (
307  boost::array<size_type,reference_element::max_variant>& loc_ndof_by_variant,
308  array<size_type,distributed>& idof2ios_dis_idof) const = 0;
309 };
310 #endif // _RHEOLEF_HAVE_MPI
311 // NOTE: since geo_rep<seq> contains sequential arrays for vertices and elts,
312 // the geo_rep<mpi> cannot derive from geo_rep<seq>. The solution is to
313 // derive both geo_rep<seq> and geo_rep<mpi> classes from a generic base class
314 template <class T, class M>
315 class geo_base_rep : public geo_abstract_rep<T,M> {
316 public:
317 
318  typedef geo_abstract_rep<T,M> base;
319  typedef typename base::size_type size_type;
320  typedef typename base::node_type node_type;
321  typedef typename base::variant_type variant_type;
322  typedef typename base::iterator iterator;
323  typedef typename base::const_iterator const_iterator;
324  typedef typename base::iterator_by_variant iterator_by_variant;
325  typedef typename base::const_iterator_by_variant const_iterator_by_variant;
326  typedef typename base::reference reference;
327  typedef typename base::const_reference const_reference;
328  typedef typename base::coordinate_type coordinate_type;
329 
330 
331  geo_base_rep ();
332  geo_base_rep (const geo_base_rep<T,M>&);
333 
334  void build_from_list (
335  const geo_basic<T,M>& lambda,
336  const array<point_basic<T>,M>& node_list,
337  const boost::array<array<geo_element_auto<heap_allocator<size_type> >,M>,
338  reference_element::max_variant>& elt_list);
339 
340  ~geo_base_rep ();
341 
342 
343  size_type variant() const { return geo_abstract_base_rep<T>::geo; }
344  std::string familyname() const { return _name; }
345  std::string name() const;
346  size_type serial_number() const { return _serial_number; }
347  size_type dimension() const { return _dimension; }
348  size_type map_dimension() const { return _map_dimension; }
349  coordinate_type coordinate_system() const { return _sys_coord; }
350  void set_coordinate_system (coordinate_type sys_coord) { _sys_coord = sys_coord; }
351  void set_name (std::string name) { _name = name; }
352  void set_dimension (size_type dim) { _dimension = dim; }
353  void set_serial_number(size_type i) { _serial_number = i; }
354  const basis_basic<T>& get_piola_basis() const { return _numbering.get_basis(); }
355  const node_type& xmin() const { return _xmin; }
356  const node_type& xmax() const { return _xmax; }
357  const geo_size& sizes() const { return _gs; }
358  const geo_size& ios_sizes() const { return _gs; }
359 
360  const distributor& geo_element_ownership (size_type dim) const { return _gs.ownership_by_dimension[dim]; }
361  const_reference get_geo_element (size_type dim, size_type ige) const;
362  const_reference dis_get_geo_element (size_type dim, size_type dis_ige) const;
363 
364  const_iterator_by_variant begin_by_variant (variant_type variant) const;
365  iterator_by_variant begin_by_variant (variant_type variant);
366  const_iterator_by_variant end_by_variant (variant_type variant) const;
367  iterator_by_variant end_by_variant (variant_type variant);
368 
369  const node_type& node (size_type inod) const { return _node [inod]; }
370  const node_type& dis_node (size_type dis_inod) const { return _node.dis_at (dis_inod); }
371  const node_type& node (const geo_element& K, size_type loc_inod) const;
372  void dis_inod (const geo_element& K, std::vector<size_type>& dis_inod) const
373  { _numbering.dis_idof (_gs, K, dis_inod); }
374  node_type piola (const geo_element& K, const node_type& hat_x) const;
375  const array<node_type,M>& get_nodes() const { return _node; }
376 
377  void set_nodes (const array<node_type,M>& x) { _node = x; _node.reset_dis_indexes(); compute_bbox(); }
378  size_type n_domain_indirect () const { return _domains.size(); }
379  bool have_domain_indirect (const std::string& name) const;
380  const domain_indirect_basic<M>& get_domain_indirect (size_type i) const { return _domains[i]; }
381  const domain_indirect_basic<M>& get_domain_indirect (const std::string& name) const;
382  void insert_domain_indirect (const domain_indirect_basic<M>& dom) const;
383 
384  size_type seq_locate (
385  const point_basic<T>& x,
386  size_type dis_ie_guest = std::numeric_limits<size_type>::max()) const;
387  size_type dis_locate (
388  const point_basic<T>& x,
389  size_type dis_ie_guest = std::numeric_limits<size_type>::max()) const;
390  size_type seq_trace_move (
391  const point_basic<T>& x,
392  const point_basic<T>& v,
393  point_basic<T>& y) const;
394  size_type dis_trace_move (
395  const point_basic<T>& x,
396  const point_basic<T>& v,
397  point_basic<T>& y) const;
398 
399  size_type neighbour (size_type ie, size_type loc_isid) const;
400  void neighbour_guard() const;
401 
402 
403  reference get_geo_element (size_type dim, size_type ige);
404  iterator begin (size_type dim);
405  iterator end (size_type dim);
406 
407 
408  size_type size(size_type dim) const;
409  size_type dis_size(size_type dim) const;
410  const distributor& ownership() const { return geo_element_ownership (map_dimension()); }
411  const distributor& vertex_ownership() const { return geo_element_ownership (0); }
412  const communicator& comm() const { return ownership().comm(); }
413  size_type order() const { return get_piola_basis().degree(); }
414 
415  size_type n_node() const { return _node. size(); }
416  size_type dis_n_node() const { return _node.dis_size(); }
417  size_type n_vertex() const { return size (0); }
418  size_type size() const { return size (map_dimension()); }
419  size_type dis_n_vertex() const { return dis_size (0); }
420  size_type dis_size() const { return dis_size (map_dimension()); }
421  size_type dis_n_edge() const { return dis_size (1); }
422  size_type dis_n_face() const { return dis_size (2); }
423  size_type dis_inod2dis_iv (size_type dis_inod) const;
424  size_type dis_iv2dis_inod (size_type dis_iv) const;
425 
426  const_reference operator[] (size_type ie) const { return get_geo_element (map_dimension(), ie); }
427  reference operator[] (size_type ie) { return get_geo_element (map_dimension(), ie); }
428 
429  const_iterator begin (size_type dim) const { return base::begin (dim); }
430  const_iterator end (size_type dim) const { return base::end (dim); }
431 
432  const_iterator begin() const { return base::begin (map_dimension()); }
433  const_iterator end() const { return base::end (map_dimension()); }
434 
435  const_iterator begin_edge() const { return base::begin (1); }
436  const_iterator end_edge() const { return base::end (1); }
437  const_iterator begin_face() const { return base::begin (2); }
438  const_iterator end_face() const { return base::end (2); }
439 
440 protected:
441  void compute_bbox();
442  void init_neighbour() const;
443  std::string _name;
444  size_type _version;
445  size_type _serial_number;
446  hack_array<geo_element_hack,M> _geo_element [reference_element::max_variant];
447  size_type _map_dimension;
448  geo_size _gs; // counters by geo_element dimension: 0,1,2,3
449  mutable std::vector<domain_indirect_basic<M> > _domains;
450  bool _have_connectivity; // e.g.list of edges in a 2d triangular mesh
451  mutable bool _have_neighbour; // inter-element connectivity
452  array<node_type, M> _node;
453  size_type _dimension;
454  coordinate_type _sys_coord;
455  node_type _xmin; // bounding box
456  node_type _xmax;
457  numbering<T,M> _numbering;
458  geo_locate<T,M> _locator;
459  geo_trace_ray_boundary<T,M> _tracer_ray_boundary;
460 };
461 template <class T, class M>
462 inline
463 geo_base_rep<T,M>::geo_base_rep ()
464  : geo_abstract_rep<T,M>(),
465  _name("*nogeo*"),
466  _version(0),
467  _serial_number(0),
468  _geo_element(),
469  _map_dimension(0),
470  _gs(),
471  _domains(),
472  _have_connectivity(false),
473  _have_neighbour(false),
474  _node(),
475  _dimension(0),
476  _sys_coord(space_constant::cartesian),
477  _xmin(),
478  _xmax(),
479  _numbering("P1"),
480  _locator(),
481  _tracer_ray_boundary()
482 {
483 }
484 template <class T, class M>
485 inline
486 geo_base_rep<T,M>::geo_base_rep (const geo_base_rep<T,M>& o)
487  : geo_abstract_rep<T,M>(o),
488  _name (o._name),
489  _version (o._version),
490  _serial_number (o._serial_number),
491  _geo_element (o._geo_element),
492  _map_dimension (o._map_dimension),
493  _gs (o._gs),
494  _domains (o._domains),
495  _have_connectivity (o._have_connectivity),
496  _have_neighbour (o._have_neighbour),
497  _node (o._node),
498  _dimension (o._dimension),
499  _sys_coord (o._sys_coord),
500  _xmin (o._xmin),
501  _xmax (o._xmax),
502  _numbering (o._numbering),
503  _locator (o._locator),
504  _tracer_ray_boundary(o._tracer_ray_boundary)
505 {
506  trace_macro ("physical copy of geo_base_rep: dim="<<o.dimension()<<", dis_size="<<o.dis_size());
507 }
508 template <class T, class M>
509 inline
510 void
511 geo_base_rep<T,M>::neighbour_guard() const
512 {
513  if (_have_neighbour) return;
514  _have_neighbour = true;
515  init_neighbour();
516 }
517 template <class T, class M>
518 typename geo_base_rep<T,M>::iterator_by_variant
519 geo_base_rep<T,M>::begin_by_variant (variant_type variant)
520 {
521  return _geo_element [variant].begin();
522 }
523 template <class T, class M>
524 typename geo_base_rep<T,M>::const_iterator_by_variant
525 geo_base_rep<T,M>::begin_by_variant (variant_type variant) const
526 {
527  return _geo_element [variant].begin();
528 }
529 template <class T, class M>
530 typename geo_base_rep<T,M>::iterator_by_variant
531 geo_base_rep<T,M>::end_by_variant (variant_type variant)
532 {
533  return _geo_element [variant].end();
534 }
535 template <class T, class M>
536 typename geo_base_rep<T,M>::const_iterator_by_variant
537 geo_base_rep<T,M>::end_by_variant (variant_type variant) const
538 {
539  return _geo_element [variant].end();
540 }
541 template <class T, class M>
542 inline
543 typename geo_base_rep<T,M>::iterator
544 geo_base_rep<T,M>::begin (size_type dim)
545 {
546  variant_type variant = reference_element::first_variant_by_dimension(dim);
547  iterator_by_variant iter = begin_by_variant (variant);
548  iterator res = iterator (dim, variant, iter, *this);
549  return res;
550 }
551 template <class T, class M>
552 inline
553 typename geo_base_rep<T,M>::iterator
554 geo_base_rep<T,M>::end (size_type dim)
555 {
556  variant_type variant = reference_element::last_variant_by_dimension(dim);
557  iterator_by_variant iter = end_by_variant (variant - 1);
558  iterator res = iterator (dim, variant, iter, *this);
559  return res;
560 }
561 template <class T>
562 inline
563 typename geo_abstract_base_rep<T>::const_iterator
564 geo_abstract_base_rep<T>::begin (size_type dim) const
565 {
566  variant_type variant = reference_element::first_variant_by_dimension(dim);
567  const_iterator_by_variant iter = begin_by_variant (variant);
568  return const_iterator (dim, variant, iter, *this);
569 }
570 template <class T>
571 inline
572 typename geo_abstract_base_rep<T>::const_iterator
573 geo_abstract_base_rep<T>::end (size_type dim) const
574 {
575  variant_type variant = reference_element::last_variant_by_dimension(dim);
576  const_iterator_by_variant iter = end_by_variant (variant - 1);
577  return const_iterator (dim, variant, iter, *this);
578 }
579 template <class T, class M> class geo_rep {};
580 
581 template <class T>
582 class geo_rep<T,sequential> : public geo_base_rep<T,sequential> {
583 public:
584 
585  typedef geo_base_rep<T,sequential> base;
586  typedef typename base::size_type size_type;
587  typedef typename base::node_type node_type;
588  typedef typename base::variant_type variant_type;
589  typedef typename base::reference reference;
590  typedef typename base::const_reference const_reference;
591  typedef typename base::iterator iterator;
592  typedef typename base::const_iterator const_iterator;
593  typedef typename base::iterator_by_variant iterator_by_variant;
594  typedef typename base::const_iterator_by_variant const_iterator_by_variant;
595  typedef typename base::coordinate_type coordinate_type;
596 
597 
598  geo_rep();
599  geo_rep (const geo_rep<T,sequential>&);
600  geo_abstract_rep<T,sequential>* clone() const;
601 
602  geo_rep (
603  const geo_basic<T,sequential>& lambda,
604  const array<point_basic<T>,sequential>& node_list,
605  const boost::array<array<geo_element_auto<heap_allocator<size_type> >,sequential>,
606  reference_element::max_variant>& elt_list);
607 
608  void build_from_domain (
609  const domain_indirect_rep<sequential>& indirect,
610  const geo_abstract_rep<T,sequential>& omega,
611  std::map<size_type,size_type>& bgd_ie2dom_ie);
612 
613 
614  void locate (
615  const array<point_basic<T>, sequential>& x,
616  array<size_type, sequential>& dis_ie,
617  bool do_check = false) const;
618  void trace_ray_boundary (
619  const array<point_basic<T>,sequential>& x,
620  const array<point_basic<T>,sequential>& v,
621  array<size_type, sequential>& dis_ie,
622  array<point_basic<T>,sequential>& y,
623  bool do_check = false) const;
624  void trace_move (
625  const array<point_basic<T>,sequential>& x,
626  const array<point_basic<T>,sequential>& v,
627  array<size_type, sequential>& dis_ie,
628  array<point_basic<T>,sequential>& y) const;
629 
630 
631  size_type map_dimension () const { return base::map_dimension(); }
632  const distributor& geo_element_ownership(size_type dim) const { return base::geo_element_ownership (dim); }
633  const_reference get_geo_element (size_type dim, size_type ige) const { return base::get_geo_element (dim, ige); }
634  reference get_geo_element (size_type dim, size_type ige) { return base::get_geo_element (dim, ige); }
635  iterator begin (size_type dim) { return base::begin(dim); }
636  iterator end (size_type dim) { return base::end (dim); }
637  const_iterator begin (size_type dim) const { return base::begin(dim); }
638  const_iterator end (size_type dim) const { return base::end (dim); }
639 
640  const domain_indirect_basic<sequential>& get_domain_indirect (size_type idom) const { return base::get_domain_indirect (idom); }
641  const domain_indirect_basic<sequential>& get_domain_indirect (const std::string& name) const { return base::get_domain_indirect (name); }
642 
643 
644  const distributor& vertex_ownership() const { return geo_element_ownership(0); }
645  const_reference operator[] (size_type ie) const { return get_geo_element (map_dimension(), ie); }
646  reference operator[] (size_type ie) { return get_geo_element (map_dimension(), ie); }
647 #ifdef TO_CLEAN
648  size_type dis_inod2dis_iv (size_type dis_inod) const { return dis_inod; }
649  size_type dis_iv2dis_inod (size_type dis_iv) const { return dis_iv; }
650 #endif // TO_CLEAN
651 
652 // i/o:
653 
654  idiststream& get (idiststream&);
655  odiststream& put_geo (odiststream&) const;
656  odiststream& put (odiststream& ops) const { return put_geo(ops); }
657 
658  void dump (std::string name) const;
659  void load (std::string name, const communicator& = communicator());
660  bool check(bool verbose) const;
661 
662 
663  void reset_order (size_type order);
664  void build_from_data (
665  const geo_header& hdr,
666  const array<node_type, sequential>& node,
667  boost::array<array<geo_element_auto<heap_allocator<size_type> >,sequential, heap_allocator<size_type> >, reference_element::max_variant>&
668  tmp_geo_element,
669  bool do_upgrade);
670 
671 protected:
672  idiststream& get_standard (idiststream&, const geo_header&);
673  idiststream& get_upgrade (idiststream&, const geo_header&);
674  void build_connectivity (
675  boost::array<array<geo_element_auto<heap_allocator<size_type> >,sequential, heap_allocator<size_type> >, reference_element::max_variant>& tmp_geo_element);
676 
677  void build_connectivity_sides (
678  size_type side_dim,
679  boost::array<array<geo_element_auto<heap_allocator<size_type> >,sequential, heap_allocator<size_type> >, reference_element::max_variant>& tmp_geo_element);
680  void set_element_side_index (size_type side_dim);
681  void domain_set_side_part1 (
682  const domain_indirect_rep<sequential>& indirect,
683  const geo_abstract_rep<T,sequential>& bgd_omega,
684  size_type sid_dim,
685  array<size_type,sequential>& bgd_isid2dom_dis_isid,
686  array<size_type,sequential>& dom_isid2bgd_isid,
687  array<size_type,sequential>& dom_isid2dom_ios_dis_isid,
688  size_type size_by_variant [reference_element::max_variant]);
689  void domain_set_side_part2 (
690  const domain_indirect_rep<sequential>& indirect,
691  const geo_abstract_rep<T,sequential>& bgd_omega,
692  array<size_type,sequential>& bgd_iv2dom_dis_iv,
693  size_type sid_dim,
694  array<size_type,sequential>& bgd_isid2dom_dis_isid,
695  array<size_type,sequential>& dom_isid2bgd_isid,
696  array<size_type,sequential>& dom_isid2dom_ios_dis_isid,
697  size_type size_by_variant [reference_element::max_variant]);
698  void build_external_entities () {} // for distributed compat
699 };
700 template <class T>
701 inline
702 geo_rep<T,sequential>::geo_rep ()
703  : geo_base_rep<T,sequential>()
704 {
705 }
706 template <class T>
707 inline
708 geo_rep<T,sequential>::geo_rep (const geo_rep<T,sequential>& omega)
709  : geo_base_rep<T,sequential>(omega)
710 {
711 }
712 template <class T>
713 inline
714 geo_abstract_rep<T,sequential>*
715 geo_rep<T,sequential>::clone() const
716 {
717  typedef geo_rep<T,sequential> rep;
718  return new_macro (rep(*this));
719 }
720 
721 #ifdef _RHEOLEF_HAVE_MPI
722 template <class T>
723 class geo_rep<T,distributed> : public geo_base_rep<T,distributed> {
724 public:
725 
726  typedef geo_base_rep<T,distributed> base;
727  typedef typename base::size_type size_type;
728  typedef typename base::node_type node_type;
729  typedef typename base::variant_type variant_type;
730  typedef typename base::node_map_type node_map_type;
731  typedef typename base::reference reference;
732  typedef typename base::const_reference const_reference;
733  typedef typename base::iterator iterator;
734  typedef typename base::const_iterator const_iterator;
735  typedef typename base::iterator_by_variant iterator_by_variant;
736  typedef typename base::const_iterator_by_variant const_iterator_by_variant;
737  typedef typename base::coordinate_type coordinate_type;
738 
739 
740  geo_rep ();
741  geo_rep (const geo_rep<T,distributed>&);
742  geo_abstract_rep<T,distributed>* clone() const;
743 
744  void build_from_domain (
745  const domain_indirect_rep<distributed>& indirect,
746  const geo_abstract_rep<T,distributed>& omega,
747  std::map<size_type,size_type>& bgd_ie2dom_ie);
748 
749  geo_rep (
750  const geo_basic<T,distributed>& lambda,
751  const array<point_basic<T>,distributed>& node_list,
752  const boost::array<array<geo_element_auto<heap_allocator<size_type> >,distributed>,
753  reference_element::max_variant>& elt_list);
754 
755 
756  distributor geo_element_ios_ownership (size_type dim) const;
757 
758  size_type ige2ios_dis_ige (size_type dim, size_type ige) const;
759  size_type dis_ige2ios_dis_ige (size_type dim, size_type dis_ige) const;
760  size_type ios_ige2dis_ige (size_type dim, size_type ios_ige) const;
761  const geo_size& ios_sizes() const { return _ios_gs; }
762  void locate (
763  const array<point_basic<T>,distributed>& x,
764  array<size_type, distributed>& dis_ie,
765  bool do_check = false) const;
766  void trace_ray_boundary (
767  const array<point_basic<T>,distributed>& x,
768  const array<point_basic<T>,distributed>& v,
769  array<size_type, distributed>& dis_ie,
770  array<point_basic<T>,distributed>& y,
771  bool do_check = false) const;
772  void trace_move (
773  const array<point_basic<T>,distributed>& x,
774  const array<point_basic<T>,distributed>& v,
775  array<size_type, distributed>& dis_ie,
776  array<point_basic<T>,distributed>& y) const;
777 
778 
779  size_type map_dimension () const { return base::map_dimension(); }
780  size_type size (size_type dim) const { return base::size(dim); }
781  const distributor& geo_element_ownership(size_type dim) const { return base::geo_element_ownership (dim); }
782  const_reference get_geo_element (size_type dim, size_type ige) const { return base::get_geo_element (dim, ige); }
783  reference get_geo_element (size_type dim, size_type ige) { return base::get_geo_element (dim, ige); }
784 
785  const_iterator begin (size_type dim) const { return base::begin(dim); }
786  const_iterator end (size_type dim) const { return base::end (dim); }
787  iterator begin (size_type dim) { return base::begin(dim); }
788  iterator end (size_type dim) { return base::end (dim); }
789 
790  const domain_indirect_basic<distributed>& get_domain_indirect (size_type idom) const { return base::get_domain_indirect (idom); }
791  const domain_indirect_basic<distributed>& get_domain_indirect (const std::string& name) const { return base::get_domain_indirect (name); }
792 
793 
794  size_type size () const { return size (map_dimension()); }
795  const distributor& vertex_ownership() const { return geo_element_ownership(0); }
796  const_reference operator[] (size_type ie) const { return get_geo_element (map_dimension(), ie); }
797 #ifdef TODO
798  reference operator[] (size_type ie) { return get_geo_element (map_dimension(), ie); }
799 #endif // TODO
800 
801 
802  void reset_order (size_type order);
803 
804 // i/o:
805 
806  idiststream& get (idiststream&);
807  odiststream& put (odiststream&) const;
808  void dump (std::string name) const;
809  void load (std::string name, const communicator& comm);
810  bool check(bool verbose) const;
811 
812 
813  void set_ios_permutation (
814  boost::array<size_type,reference_element::max_variant>& loc_ndof_by_variant,
815  array<size_type,distributed>& idof2ios_dis_idof) const;
816 
817 protected:
818  void build_external_entities ();
819  void set_element_side_index (size_type side_dim);
820  void domain_set_side_part1 (
821  const domain_indirect_rep<distributed>& indirect,
822  const geo_abstract_rep<T,distributed>& bgd_omega,
823  size_type sid_dim,
824  array<size_type>& bgd_isid2dom_dis_isid,
825  array<size_type>& dom_isid2bgd_isid,
826  array<size_type>& dom_isid2dom_ios_dis_isid,
827  size_type size_by_variant [reference_element::max_variant]);
828  void domain_set_side_part2 (
829  const domain_indirect_rep<distributed>& indirect,
830  const geo_abstract_rep<T,distributed>& bgd_omega,
831  array<size_type>& bgd_iv2dom_dis_iv,
832  size_type sid_dim,
833  array<size_type>& bgd_isid2dom_dis_isid,
834  array<size_type>& dom_isid2bgd_isid,
835  array<size_type>& dom_isid2dom_ios_dis_isid,
836  size_type size_by_variant [reference_element::max_variant]);
837 
838  void node_renumbering (const distributor& ios_node_ownership);
839 
840  array<size_type> _inod2ios_dis_inod; // permutation for node
841  array<size_type> _ios_inod2dis_inod; // reverse permutation for node
842  array<size_type> _ios_ige2dis_ige[4]; // reverse permutation for geo_element[dim]
843  geo_size _ios_gs;
844  array<size_type> _igev2ios_dis_igev [reference_element::max_variant];
845  array<size_type> _ios_igev2dis_igev [reference_element::max_variant];
846 };
847 template <class T>
848 inline
849 geo_rep<T,distributed>::geo_rep ()
850  : geo_base_rep<T,distributed>(),
851  _inod2ios_dis_inod(),
852  _ios_inod2dis_inod(),
853  _ios_ige2dis_ige(),
854  _ios_gs(),
855  _igev2ios_dis_igev(),
856  _ios_igev2dis_igev()
857 {
858 }
859 template <class T>
860 inline
861 geo_rep<T,distributed>::geo_rep (const geo_rep<T,distributed>& o)
862  : geo_base_rep<T,distributed>(o),
863  _inod2ios_dis_inod(o._inod2ios_dis_inod),
864  _ios_inod2dis_inod(o._ios_inod2dis_inod),
865  _ios_ige2dis_ige (o._ios_ige2dis_ige),
866  _ios_gs (o._ios_gs),
867  _igev2ios_dis_igev(o._igev2ios_dis_igev),
868  _ios_igev2dis_igev(o._ios_igev2dis_igev)
869 {
870 }
871 template <class T>
872 inline
873 geo_abstract_rep<T,distributed>*
874 geo_rep<T,distributed>::clone() const
875 {
876  typedef geo_rep<T,distributed> rep;
877  return new_macro (rep(*this));
878 }
879 #endif // _RHEOLEF_HAVE_MPI
880 
881 template <class T, class M = rheo_default_memory_model>
882 class geo_basic {
883 public:
884  typedef M memory_type;
885 };
894 
895 template<class T, class M> geo_basic<T,M> geo_load (const std::string& name);
896 
897 template <class T>
898 class geo_basic<T,sequential> : public smart_pointer_clone<geo_abstract_rep<T,sequential> > {
899 public:
900 
901 
902  typedef sequential memory_type;
903  typedef geo_abstract_rep<T,sequential> rep;
904  typedef geo_rep<T,sequential> rep_geo_rep;
905  typedef smart_pointer_clone<rep> base;
906  typedef typename rep::size_type size_type;
907  typedef typename rep::node_type node_type;
908  typedef typename rep::variant_type variant_type;
909  typedef typename rep::reference reference;
910  typedef typename rep::const_reference const_reference;
911  typedef typename rep::iterator iterator;
912  typedef typename rep::const_iterator const_iterator;
913  typedef typename rep::iterator_by_variant iterator_by_variant;
914  typedef typename rep::const_iterator_by_variant const_iterator_by_variant;
915  typedef typename rep::coordinate_type coordinate_type;
916 
917 
918  geo_basic ();
919  geo_basic (std::string name, const communicator& comm = communicator());
920  void load (std::string name, const communicator& comm = communicator());
921  geo_basic (const domain_indirect_basic<sequential>& dom, const geo_basic<T,sequential>& omega);
922 
923  geo_basic (
924  const geo_basic<T,sequential>& lambda,
925  const array<point_basic<T>,sequential>& node_list,
926  const boost::array<array<geo_element_auto<heap_allocator<size_type> >,sequential>,
927  reference_element::max_variant>& elt_list)
928  : base (new_macro(rep_geo_rep(lambda,node_list,elt_list))) {}
929 
930 
931  std::string name() const { return base::data().name(); }
932  std::string familyname() const { return base::data().familyname(); }
933  size_type dimension() const { return base::data().dimension(); }
934  size_type map_dimension() const { return base::data().map_dimension(); }
935  size_type serial_number() const { return base::data().serial_number(); }
936  coordinate_type coordinate_system() const { return base::data().coordinate_system(); }
937  std::string coordinate_system_name() const { return space_constant::coordinate_system_name(coordinate_system()); }
938  const basis_basic<T>& get_piola_basis() const { return base::data().get_piola_basis(); }
939  size_type order() const { return base::data().get_piola_basis().degree(); }
940  const node_type& xmin() const { return base::data().xmin(); }
941  const node_type& xmax() const { return base::data().xmax(); }
942  const distributor& geo_element_ownership(size_type dim) const { return base::data().geo_element_ownership(dim); }
943  const geo_size& sizes() const { return base::data().sizes(); }
944  const geo_size& ios_sizes() const { return base::data().ios_sizes(); }
945  const_reference get_geo_element (size_type dim, size_type ige) const { return base::data().get_geo_element (dim, ige); }
946  reference get_geo_element (size_type dim, size_type ige) { return base::data().get_geo_element (dim, ige); }
947  const_reference dis_get_geo_element (size_type dim, size_type dis_ige) const
948  { return get_geo_element (dim, dis_ige); }
949  size_type neighbour (size_type ie, size_type loc_isid) const {
950  return base::data().neighbour (ie, loc_isid); }
951  void neighbour_guard() const { base::data().neighbour_guard(); }
952  size_type n_node() const { return base::data().n_node(); }
953  const node_type& node(size_type inod) const { return base::data().node(inod); }
954  const node_type& dis_node(size_type dis_inod) const { return base::data().dis_node(dis_inod); }
955  void dis_inod (const geo_element& K, std::vector<size_type>& dis_inod) const {
956  return base::data().dis_inod(K,dis_inod); }
957  node_type piola (const geo_element& K, const node_type& hat_x) const { return base::data().piola (K, hat_x); }
958  const array<node_type,sequential>& get_nodes() const { return base::data().get_nodes(); }
959  size_type dis_inod2dis_iv (size_type dis_inod) const { return base::data().dis_inod2dis_iv(dis_inod); }
960 
961  size_type n_domain_indirect () const { return base::data().n_domain_indirect (); }
962  bool have_domain_indirect (const std::string& name) const { return base::data().have_domain_indirect (name); }
963  const domain_indirect_basic<sequential>& get_domain_indirect (size_type i) const {
964  return base::data().get_domain_indirect (i); }
965  const domain_indirect_basic<sequential>& get_domain_indirect (const std::string& name) const {
966  return base::data().get_domain_indirect (name); }
967  void insert_domain_indirect (const domain_indirect_basic<sequential>& dom) const {
968  base::data().insert_domain_indirect (dom); }
969 
970  size_type n_domain () const { return base::data().n_domain_indirect (); }
971  geo_basic<T,sequential> get_domain (size_type i) const;
972  geo_basic<T,sequential> operator[] (const std::string& name) const;
973 
974  size_type seq_locate (
975  const point_basic<T>& x,
976  size_type dis_ie_guest = std::numeric_limits<size_type>::max()) const
977  { return base::data().seq_locate (x, dis_ie_guest); }
978  size_type dis_locate (
979  const point_basic<T>& x,
980  size_type dis_ie_guest = std::numeric_limits<size_type>::max()) const
981  { return base::data().dis_locate (x, dis_ie_guest); }
982  void locate (
983  const array<point_basic<T>, sequential>& x,
984  array<size_type, sequential>& dis_ie) const
985  { return base::data().locate (x, dis_ie); }
986  size_type seq_trace_move (
987  const point_basic<T>& x,
988  const point_basic<T>& v,
989  point_basic<T>& y) const
990  { return base::data().seq_trace_move (x,v,y); }
991  size_type dis_trace_move (
992  const point_basic<T>& x,
993  const point_basic<T>& v,
994  point_basic<T>& y) const
995  { return base::data().dis_trace_move (x,v,y); }
996  void trace_ray_boundary (
997  const array<point_basic<T>,sequential>& x,
998  const array<point_basic<T>,sequential>& v,
999  array<size_type, sequential>& dis_ie,
1000  array<point_basic<T>,sequential>& y) const
1001  { return base::data().trace_ray_boundary (x,v,dis_ie,y); }
1002  void trace_move (
1003  const array<point_basic<T>,sequential>& x,
1004  const array<point_basic<T>,sequential>& v,
1005  array<size_type, sequential>& dis_ie,
1006  array<point_basic<T>,sequential>& y) const
1007  { return base::data().trace_move (x,v,dis_ie,y); }
1008 
1009 
1010  void set_name (std::string name);
1011  void set_dimension (size_type dim);
1012  void set_serial_number (size_type i);
1013  void reset_order (size_type order);
1014  void set_coordinate_system (coordinate_type sys_coord);
1015  void set_coordinate_system (std::string sys_coord_name) { set_coordinate_system (space_constant::coordinate_system(sys_coord_name)); }
1016  void set_nodes (const array<node_type,sequential>& x);
1017  void build_from_data (
1018  const geo_header& hdr,
1019  const array<node_type, sequential>& node,
1020  boost::array<array<geo_element_auto<heap_allocator<size_type> >,sequential, heap_allocator<size_type> >, reference_element::max_variant>&
1021  tmp_geo_element,
1022  bool do_upgrade);
1023 
1024 
1025  const communicator& comm() const { return geo_element_ownership (0).comm(); }
1026  size_type size(size_type dim) const { return base::data().geo_element_ownership(dim).size(); }
1027  size_type dis_size(size_type dim) const { return base::data().geo_element_ownership(dim).dis_size(); }
1028  size_type size() const { return size (map_dimension()); }
1029  size_type dis_size() const { return dis_size (map_dimension()); }
1030  size_type n_vertex() const { return size (0); }
1031  size_type dis_n_vertex() const { return dis_size (0); }
1032  const_reference operator[] (size_type ie) const { return get_geo_element (map_dimension(), ie); }
1033  reference operator[] (size_type ie) { return get_geo_element (map_dimension(), ie); }
1034  const_iterator begin (size_type dim) const { return base::data().begin(dim); }
1035  const_iterator end (size_type dim) const { return base::data().end (dim); }
1036  const_iterator begin () const { return begin(map_dimension()); }
1037  const_iterator end () const { return end (map_dimension()); }
1038 
1039  const_iterator_by_variant begin_by_variant (variant_type variant) const
1040  { return base::data().begin_by_variant (variant); }
1041  const_iterator_by_variant end_by_variant (variant_type variant) const
1042  { return base::data(). end_by_variant (variant); }
1043 
1044  const geo_basic<T,sequential>& get_background_geo() const; // code in geo_domain.h
1045  geo_basic<T,sequential> get_background_domain() const;
1046 
1047 
1048  size_type ige2ios_dis_ige (size_type dim, size_type ige) const { return ige; }
1049  size_type dis_ige2ios_dis_ige (size_type dim, size_type dis_ige) const { return dis_ige; }
1050  size_type ios_ige2dis_ige (size_type dim, size_type ios_ige) const { return ios_ige; }
1051 
1052 
1053  bool operator== (const geo_basic<T,sequential>& omega2) const { return base::data().operator== (omega2.data()); }
1054 
1055 // i/o:
1056 
1057  idiststream& get (idiststream& ips);
1058  odiststream& put (odiststream& ops) const;
1059  void save (std::string filename = "") const;
1060  void dump (std::string name) const { base::data().dump (name); }
1061  bool check (bool verbose = true) const { return base::data().check(verbose); }
1062 };
1063 @endcode
1064 template <class T>
1065 inline
1066 geo_basic<T,sequential>::geo_basic()
1067  : base (new_macro((geo_rep<T,sequential>)))
1068 {
1069 }
1070 template <class T>
1071 inline
1072 geo_basic<T,sequential>::geo_basic (std::string name, const communicator& comm)
1073  : base (0)
1074 {
1075 #ifdef TO_CLEAN
1076  // allocate a new geo_rep<mpi> object ; TODO geo_rep(name,comm)
1077  geo_rep<T,sequential>* ptr = new_macro((geo_rep<T,sequential>));
1078  ptr->load (name);
1079  base::operator= (ptr);
1080 #endif // TO_CLEAN
1081  base::operator= (geo_load<T,sequential>(name));
1082 }
1083 template <class T>
1084 inline
1085 void
1086 geo_basic<T,sequential>::load (std::string name, const communicator& comm)
1087 {
1088 #ifdef TO_CLEAN
1089  // allocate a new geo_rep<mpi> object ; TODO geo_rep(name,comm)
1090  geo_rep<T,sequential>* ptr = new_macro((geo_rep<T,sequential>));
1091  ptr->load (name);
1092  base::operator= (ptr);
1093 #endif // TO_CLEAN
1094  base::operator= (geo_load<T,sequential>(name));
1095 }
1096 #ifdef TO_CLEAN
1097 template <class T>
1098 inline
1099 idiststream&
1100 geo_basic<T,sequential>::get (idiststream& ips)
1101 {
1102  geo_rep<T,sequential>* ptr = new_macro((geo_rep<T,sequential>));
1103  ptr->get (ips);
1104  base::operator= (ptr);
1105  return ips;
1106 }
1107 #endif // TO_CLEAN
1108 #ifdef _RHEOLEF_HAVE_MPI
1109 
1110 template <class T>
1111 class geo_basic<T,distributed> : public smart_pointer_clone<geo_abstract_rep<T,distributed> > {
1112 public:
1113 
1114 
1115  typedef distributed memory_type;
1116  typedef geo_abstract_rep<T,distributed> rep;
1117  typedef geo_rep<T,distributed> rep_geo_rep;
1118  typedef smart_pointer_clone<rep> base;
1119  typedef typename rep::size_type size_type;
1120  typedef typename rep::node_type node_type;
1121  typedef typename rep::variant_type variant_type;
1122  typedef typename rep::node_map_type node_map_type;
1123  typedef typename rep::reference reference;
1124  typedef typename rep::const_reference const_reference;
1125  typedef typename rep::iterator iterator;
1126  typedef typename rep::const_iterator const_iterator;
1127  typedef typename rep::iterator_by_variant iterator_by_variant;
1128  typedef typename rep::const_iterator_by_variant const_iterator_by_variant;
1129  typedef typename rep::coordinate_type coordinate_type;
1130 
1131 
1132  geo_basic ();
1133  geo_basic (std::string name, const communicator& comm = communicator());
1134  void load (std::string name, const communicator& comm = communicator());
1135  geo_basic (const domain_indirect_basic<distributed>& dom, const geo_basic<T,distributed>& omega);
1136 
1137  geo_basic (
1138  const geo_basic<T,distributed>& lambda,
1139  const array<point_basic<T>,distributed>& node_list,
1140  const boost::array<array<geo_element_auto<heap_allocator<size_type> >,distributed>,
1141  reference_element::max_variant>& elt_list)
1142  : base (new_macro(rep_geo_rep(lambda,node_list,elt_list))) {}
1143 
1144 
1145  std::string name() const { return base::data().name(); }
1146  std::string familyname() const { return base::data().familyname(); }
1147  size_type dimension() const { return base::data().dimension(); }
1148  size_type map_dimension() const { return base::data().map_dimension(); }
1149  size_type serial_number() const { return base::data().serial_number(); }
1150  coordinate_type coordinate_system() const { return base::data().coordinate_system(); }
1151  std::string coordinate_system_name() const { return space_constant::coordinate_system_name(coordinate_system()); }
1152  const basis_basic<T>& get_piola_basis() const { return base::data().get_piola_basis(); }
1153  size_type order() const { return base::data().get_piola_basis().degree(); }
1154  const node_type& xmin() const { return base::data().xmin(); }
1155  const node_type& xmax() const { return base::data().xmax(); }
1156  const distributor& geo_element_ownership(size_type dim) const
1157  { return base::data().geo_element_ownership (dim); }
1158  const geo_size& sizes() const { return base::data().sizes(); }
1159  const geo_size& ios_sizes() const { return base::data().ios_sizes(); }
1160  const_reference get_geo_element (size_type dim, size_type ige) const
1161  { return base::data().get_geo_element (dim, ige); }
1162  const_reference dis_get_geo_element (size_type dim, size_type dis_ige) const
1163  { return base::data().dis_get_geo_element (dim, dis_ige); }
1164  size_type neighbour (size_type ie, size_type loc_isid) const {
1165  return base::data().neighbour (ie, loc_isid); }
1166  void neighbour_guard() const { base::data().neighbour_guard(); }
1167  distributor geo_element_ios_ownership (size_type dim) const {
1168  return base::data().geo_element_ios_ownership (dim); }
1169  size_type ige2ios_dis_ige (size_type dim, size_type ige) const {
1170  return base::data().ige2ios_dis_ige (dim,ige); }
1171  size_type dis_ige2ios_dis_ige (size_type dim, size_type dis_ige) const {
1172  return base::data().dis_ige2ios_dis_ige (dim,dis_ige); }
1173  size_type ios_ige2dis_ige (size_type dim, size_type ios_ige) const {
1174  return base::data().ios_ige2dis_ige (dim, ios_ige); }
1175  size_type n_node() const { return base::data().n_node(); }
1176  const node_type& node(size_type inod) const { return base::data().node(inod); }
1177  const node_type& dis_node(size_type dis_inod) const { return base::data().dis_node(dis_inod); }
1178  void dis_inod (const geo_element& K, std::vector<size_type>& dis_inod) const {
1179  return base::data().dis_inod(K,dis_inod); }
1180  node_type piola (const geo_element& K, const node_type& hat_x) const { return base::data().piola (K, hat_x); }
1181  const array<node_type,distributed>& get_nodes() const { return base::data().get_nodes(); }
1182  void set_nodes (const array<node_type,distributed>& x);
1183  void reset_order (size_type order);
1184  size_type dis_inod2dis_iv (size_type dis_inod) const { return base::data().dis_inod2dis_iv(dis_inod); }
1185  void set_coordinate_system (coordinate_type sys_coord);
1186  void set_coordinate_system (std::string sys_coord_name) { set_coordinate_system (space_constant::coordinate_system(sys_coord_name)); }
1187  void set_dimension (size_type dim);
1188  void set_serial_number (size_type i);
1189  void set_name (std::string name);
1190 
1191  size_type n_domain_indirect () const { return base::data().n_domain_indirect (); }
1192  bool have_domain_indirect (const std::string& name) const { return base::data().have_domain_indirect (name); }
1193  const domain_indirect_basic<distributed>& get_domain_indirect (size_type i) const {
1194  return base::data().get_domain_indirect (i); }
1195  const domain_indirect_basic<distributed>& get_domain_indirect (const std::string& name) const {
1196  return base::data().get_domain_indirect (name); }
1197  void insert_domain_indirect (const domain_indirect_basic<distributed>& dom) const {
1198  base::data().insert_domain_indirect (dom); }
1199 
1200  size_type n_domain () const { return base::data().n_domain_indirect (); }
1201  geo_basic<T,distributed> get_domain (size_type i) const;
1202  geo_basic<T,distributed> operator[] (const std::string& name) const;
1203 
1204  size_type seq_locate (
1205  const point_basic<T>& x,
1206  size_type dis_ie_guest = std::numeric_limits<size_type>::max()) const
1207  { return base::data().seq_locate (x, dis_ie_guest); }
1208  size_type dis_locate (
1209  const point_basic<T>& x,
1210  size_type dis_ie_guest = std::numeric_limits<size_type>::max()) const
1211  { return base::data().dis_locate (x, dis_ie_guest); }
1212  void locate (const array<point_basic<T>, distributed>& x, array<size_type, distributed>& dis_ie) const
1213  { return base::data().locate (x, dis_ie); }
1214  size_type seq_trace_move (
1215  const point_basic<T>& x,
1216  const point_basic<T>& v,
1217  point_basic<T>& y) const
1218  { return base::data().seq_trace_move (x,v,y); }
1219  size_type dis_trace_move (
1220  const point_basic<T>& x,
1221  const point_basic<T>& v,
1222  point_basic<T>& y) const
1223  { return base::data().dis_trace_move (x,v,y); }
1224  void trace_ray_boundary (
1225  const array<point_basic<T>,distributed>& x,
1226  const array<point_basic<T>,distributed>& v,
1227  array<size_type, distributed>& dis_ie,
1228  array<point_basic<T>,distributed>& y) const
1229  { return base::data().trace_ray_boundary (x,v,dis_ie,y); }
1230  void trace_move (
1231  const array<point_basic<T>,distributed>& x,
1232  const array<point_basic<T>,distributed>& v,
1233  array<size_type, distributed>& dis_ie,
1234  array<point_basic<T>,distributed>& y) const
1235  { return base::data().trace_move (x,v,dis_ie,y); }
1236 
1237 
1238  size_type size(size_type dim) const { return base::data().geo_element_ownership(dim).size(); }
1239  size_type dis_size(size_type dim) const { return base::data().geo_element_ownership(dim).dis_size(); }
1240  const communicator& comm() const { return geo_element_ownership (0).comm(); }
1241  size_type size() const { return size (map_dimension()); }
1242  size_type dis_size() const { return dis_size (map_dimension()); }
1243  size_type n_vertex() const { return size (0); }
1244  const_reference operator[] (size_type ie) const
1245  { return get_geo_element (map_dimension(), ie); }
1246 
1247  const_iterator begin (size_type dim) const { return base::data().begin(dim); }
1248  const_iterator end (size_type dim) const { return base::data().end (dim); }
1249  const_iterator begin () const { return begin(map_dimension()); }
1250  const_iterator end () const { return end (map_dimension()); }
1251 
1252  const_iterator_by_variant begin_by_variant (variant_type variant) const
1253  { return base::data().begin_by_variant (variant); }
1254  const_iterator_by_variant end_by_variant (variant_type variant) const
1255  { return base::data(). end_by_variant (variant); }
1256 
1257  const geo_basic<T,distributed>& get_background_geo() const; // code in geo_domain.h
1258  geo_basic<T,distributed> get_background_domain() const;
1259 
1260 
1261  bool operator== (const geo_basic<T,distributed>& omega2) const { return base::data().operator== (omega2.data()); }
1262 
1263 // i/o:
1264 
1265  odiststream& put (odiststream& ops) const { return base::data().put (ops); }
1266  idiststream& get (idiststream& ips);
1267  void save (std::string filename = "") const;
1268  bool check (bool verbose = true) const { return base::data().check(verbose); }
1269 
1270 
1271  void set_ios_permutation (
1272  boost::array<size_type,reference_element::max_variant>& loc_ndof_by_variant,
1273  array<size_type,distributed>& idof2ios_dis_idof) const
1274  { base::data().set_ios_permutation (loc_ndof_by_variant, idof2ios_dis_idof); }
1275 };
1276 @endcode
1277 #endif // _RHEOLEF_HAVE_MPI
1278 
1280 typedef geo_basic<Float,rheo_default_memory_model> geo;
1281 
1282 // inlined: geo<T,distributed>
1283 #ifdef _RHEOLEF_HAVE_MPI
1284 template <class T>
1285 inline
1286 geo_basic<T,distributed>::geo_basic()
1287  : base (new_macro((geo_rep<T,distributed>)))
1288 {
1289 }
1290 template <class T>
1291 inline
1292 geo_basic<T,distributed>::geo_basic (std::string name, const communicator& comm)
1293  : base (0)
1294 {
1295 #ifdef TO_CLEAN
1296  // allocate a new geo_rep<mpi> object ; TODO geo_rep(name,comm)
1297  geo_rep<T,distributed>* ptr = new_macro((geo_rep<T,distributed>));
1298  ptr->load (name, comm);
1299  base::operator= (ptr);
1300 #endif // TO_CLEAN
1301  base::operator= (geo_load<T,distributed>(name));
1302 }
1303 template <class T>
1304 inline
1305 void
1306 geo_basic<T,distributed>::load (std::string name, const communicator& comm)
1307 {
1308 #ifdef TO_CLEAN
1309  // allocate a new geo_rep<mpi> object ; TODO geo_rep(name,comm)
1310  geo_rep<T,distributed>* ptr = new_macro((geo_rep<T,distributed>));
1311  ptr->load (name, comm);
1312  base::operator= (ptr);
1313 #endif // TO_CLEAN
1314  base::operator= (geo_load<T,distributed>(name));
1315 }
1316 template <class T>
1317 inline
1318 idiststream&
1319 geo_basic<T,distributed>::get (idiststream& ips)
1320 {
1321  geo_rep<T,distributed>* ptr = new_macro((geo_rep<T,distributed>));
1322  ptr->get (ips);
1323  base::operator= (ptr);
1324  return ips;
1325 }
1326 #endif // _RHEOLEF_HAVE_MPI
1327 
1328 // => check that pointer to geo_abstract_rep points to a geo_rep
1329 
1330 #define _RHEOLEF_save(M) \
1331 template <class T> \
1332 void \
1333 geo_basic<T,M>::save (std::string filename) const \
1334 { \
1335  if (filename == "") filename = name(); \
1336  odiststream out (filename, "geo"); \
1337  put (out); \
1338 }
1339 #define _RHEOLEF_set_nodes(M) \
1340 template <class T> \
1341 void \
1342 geo_basic<T,M>::set_nodes (const array<node_type,M>& x) \
1343 { \
1344  geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
1345  check_macro (ptr != 0, "cannot set_name on geo_domains"); \
1346  ptr->set_nodes(x); \
1347 }
1348 #define _RHEOLEF_reset_order(M) \
1349 template <class T> \
1350 void \
1351 geo_basic<T,M>::reset_order (size_type order) \
1352 { \
1353  geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
1354  check_macro (ptr != 0, "cannot set_name on geo_domains"); \
1355  ptr->reset_order(order); \
1356 }
1357 #define _RHEOLEF_set_coordinate_system(M) \
1358 template <class T> \
1359 void \
1360 geo_basic<T,M>::set_coordinate_system (coordinate_type sys_coord) \
1361 { \
1362  geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
1363  check_macro (ptr != 0, "cannot set_name on geo_domains"); \
1364  ptr->set_coordinate_system(sys_coord); \
1365 }
1366 #define _RHEOLEF_set_dimension(M) \
1367 template <class T> \
1368 void \
1369 geo_basic<T,M>::set_dimension (size_type dim) \
1370 { \
1371  geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
1372  check_macro (ptr != 0, "cannot set_name on geo_domains"); \
1373  ptr->set_dimension(dim); \
1374 }
1375 #define _RHEOLEF_set_serial_number(M) \
1376 template <class T> \
1377 void \
1378 geo_basic<T,M>::set_serial_number (size_type i) \
1379 { \
1380  geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
1381  check_macro (ptr != 0, "cannot set_name on geo_domains"); \
1382  ptr->set_serial_number(i); \
1383 }
1384 #define _RHEOLEF_set_name(M) \
1385 template <class T> \
1386 void \
1387 geo_basic<T,M>::set_name (std::string name) \
1388 { \
1389  geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
1390  check_macro (ptr != 0, "cannot set_name on geo_domains"); \
1391  ptr->set_name(name); \
1392 }
1393 #define _RHEOLEF_build_from_data(M) \
1394 template <class T> \
1395 void \
1396 geo_basic<T,M>::build_from_data ( \
1397  const geo_header& hdr, \
1398  const array<node_type, sequential>& node, \
1399  boost::array<array<geo_element_auto<heap_allocator<size_type> >,sequential, heap_allocator<size_type> >, reference_element::max_variant>& tmp_geo_element, \
1400  bool do_upgrade) \
1401 { \
1402  geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
1403  check_macro (ptr != 0, "cannot set_name on geo_domains"); \
1404  ptr->build_from_data (hdr, node, tmp_geo_element, do_upgrade); \
1405 }
1406 
1407 _RHEOLEF_save(sequential)
1408 _RHEOLEF_set_nodes(sequential)
1409 _RHEOLEF_reset_order(sequential)
1410 _RHEOLEF_set_coordinate_system(sequential)
1411 _RHEOLEF_set_dimension(sequential)
1412 _RHEOLEF_set_serial_number(sequential)
1413 _RHEOLEF_set_name(sequential)
1414 _RHEOLEF_build_from_data(sequential)
1415 #ifdef _RHEOLEF_HAVE_MPI
1416 _RHEOLEF_save(distributed)
1417 _RHEOLEF_set_nodes(distributed)
1418 _RHEOLEF_reset_order(distributed)
1419 _RHEOLEF_set_coordinate_system(distributed)
1420 _RHEOLEF_set_dimension(distributed)
1421 _RHEOLEF_set_serial_number(distributed)
1422 _RHEOLEF_set_name(distributed)
1423 #endif // _RHEOLEF_HAVE_MPI
1424 #undef _RHEOLEF_set_nodes
1425 #undef _RHEOLEF_reset_order
1426 #undef _RHEOLEF_set_coordinate_system
1427 #undef _RHEOLEF_set_dimension
1428 #undef _RHEOLEF_set_name
1429 #undef _RHEOLEF_build_from_data
1430 // inlined: geo<T,M>
1431 template <class T, class M>
1432 inline
1433 idiststream&
1434 operator>> (idiststream& ips, geo_basic<T,M>& omega)
1435 {
1436  return omega.get (ips);
1437 }
1438 template <class T, class M>
1439 inline
1440 odiststream&
1441 operator<< (odiststream& ops, const geo_basic<T,M>& omega)
1442 {
1443  return omega.put (ops);
1444 }
1445 
1446 } // namespace rheolef
1447 #endif // _RHEOLEF_GEO_H
1448