1 #include "rheolef/geo.h"
22 #include "rheolef/geo_domain.h"
23 #include "rheolef/space_numbering.h"
24 #include "rheolef/rheostream.h"
25 #include "rheolef/iorheo.h"
35 template <
class T,
class M>
44 _have_connectivity(false),
45 _have_neighbour(false),
55 _tracer_ray_boundary(),
59 template <
class T,
class M>
63 _version (o._version),
64 _serial_number (o._serial_number),
65 _geo_element (o._geo_element),
67 _domains (o._domains),
68 _have_connectivity (o._have_connectivity),
69 _have_neighbour (o._have_neighbour),
71 _dimension (o._dimension),
72 _sys_coord (o._sys_coord),
77 _piola_basis (o._piola_basis),
78 _locator (o._locator),
79 _tracer_ray_boundary (o._tracer_ray_boundary),
80 _nearestor (o._nearestor)
82 trace_macro (
"*** PHYSICAL COPY OF GEO_BASE_REP ***");
93 trace_macro (
"*** PHYSICAL COPY OF GEO_REP<seq> ***");
99 trace_macro (
"*** CLONE OF GEO_REP<seq> ***");
101 return new_macro (rep(*
this));
110 #ifdef _RHEOLEF_HAVE_MPI
114 _inod2ios_dis_inod(),
115 _ios_inod2dis_inod(),
118 _igev2ios_dis_igev(),
123 geo_rep<T,distributed>::geo_rep (
const geo_rep<T,distributed>& o)
124 : geo_base_rep<
T,distributed>(o),
125 _inod2ios_dis_inod(o._inod2ios_dis_inod),
126 _ios_inod2dis_inod(o._ios_inod2dis_inod),
127 _ios_ige2dis_ige (o._ios_ige2dis_ige),
129 _igev2ios_dis_igev(o._igev2ios_dis_igev),
130 _ios_igev2dis_igev(o._ios_igev2dis_igev)
132 trace_macro (
"*** PHYSICAL COPY OF GEO_REP<dis> ***");
135 geo_abstract_rep<T,distributed>*
136 geo_rep<T,distributed>::clone()
const
138 trace_macro (
"*** CLONE OF GEO_REP<dis> ***");
139 typedef geo_rep<T,distributed> rep;
140 return new_macro (rep(*
this));
146 template<
class T,
class M>
149 template<
class T,
class M>
161 is_domain (
const std::string&
name, std::string& bgd_name, std::string& dom_name)
163 size_t n =
name.length();
164 if (
n <= 2 ||
name[
n-1] !=
']')
return false;
166 for (; i > 0 &&
name[i] !=
'[' &&
name[i] !=
'/'; i--) ;
168 if (i == 0 || i ==
n-2 ||
name[i] ==
'/')
return false;
170 bgd_name =
name.substr(0,i);
171 dom_name =
name.substr(i+1);
172 dom_name = dom_name.substr(0,dom_name.length()-1);
175 template<
class T,
class M>
181 std::string bgd_name, dom_name;
182 if (is_domain(root_name,bgd_name,dom_name)) {
184 return omega[dom_name];
192 #ifdef DEBUG_LOADED_TABLE_TRACE
196 omega.base::operator= (base((*iter).second,
typename base::internal()));
203 omega.geo_basic<
T,
M>::base::operator= (ptr);
206 #ifdef DEBUG_LOADED_TABLE_TRACE
217 template <
class T,
class M>
229 new_gamma.geo_basic<
T,
M>::base::operator= (geo_ptr);
238 #define _RHEOLEF_save(M) \
241 geo_basic<T,M>::save (std::string filename) const \
243 if (filename == "") filename = name(); \
244 odiststream out (filename, "geo"); \
247 #define _RHEOLEF_set_name(M) \
250 geo_basic<T,M>::set_name (std::string new_name) \
252 geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
253 check_macro (ptr != 0, "cannot set_name on geo_domains"); \
254 ptr->set_name(new_name); \
255 geo_base_rep<T,M>::loaded_map().insert (make_pair(name(), base::get_count()));\
257 #define _RHEOLEF_set_serial_number(M) \
260 geo_basic<T,M>::set_serial_number (size_type i) \
262 geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
263 check_macro (ptr != 0, "cannot set_serial_number on geo_domains"); \
264 ptr->set_serial_number(i); \
265 geo_base_rep<T,M>::loaded_map().insert (make_pair(name(), base::get_count()));\
267 #define _RHEOLEF_reset_order(M) \
270 geo_basic<T,M>::reset_order (size_type order) \
272 geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
273 check_macro (ptr != 0, "cannot reset_order on geo_domains"); \
274 ptr->reset_order(order); \
276 #define _RHEOLEF_set_nodes(M) \
279 geo_basic<T,M>::set_nodes (const disarray<node_type,M>& x) \
281 geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
282 check_macro (ptr != 0, "cannot set_nodes on geo_domains"); \
285 #define _RHEOLEF_set_coordinate_system(M) \
288 geo_basic<T,M>::set_coordinate_system (coordinate_type sys_coord) \
290 geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
291 check_macro (ptr != 0, "cannot set_coordinate_system on geo_domains"); \
292 ptr->set_coordinate_system(sys_coord); \
294 #define _RHEOLEF_set_dimension(M) \
297 geo_basic<T,M>::set_dimension (size_type dim) \
299 geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
300 check_macro (ptr != 0, "cannot set_dimension on geo_domains"); \
301 ptr->set_dimension(dim); \
303 #define _RHEOLEF_build_by_subdividing(M) \
306 geo_basic<T,M>::build_by_subdividing ( \
307 const geo_basic<T,M>& omega, \
310 geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
311 check_macro (ptr != 0, "cannot build_by_subdividing on geo_domains"); \
312 ptr->build_by_subdividing (omega, k); \
314 #define _RHEOLEF_build_from_data(M) \
317 geo_basic<T,M>::build_from_data ( \
318 const geo_header& hdr, \
319 const disarray<node_type, M>& node, \
320 std::array<disarray<geo_element_auto<>,M>, reference_element::max_variant>& tmp_geo_element, \
323 geo_rep<T,M>* ptr = dynamic_cast<geo_rep<T,M>*>(base::pointer()); \
324 check_macro (ptr != 0, "cannot build_from_data on geo_domains"); \
325 ptr->build_from_data (hdr, node, tmp_geo_element, do_upgrade); \
337 #ifdef _RHEOLEF_HAVE_MPI
350 #undef _RHEOLEF_set_nodes
351 #undef _RHEOLEF_reset_order
352 #undef _RHEOLEF_set_coordinate_system
353 #undef _RHEOLEF_set_dimension
354 #undef _RHEOLEF_set_name
355 #undef _RHEOLEF_build_from_data
356 #undef _RHEOLEF_build_by_subdividing
362 #define _RHEOLEF_reset_order(M) \
365 geo_rep<T,M>::reset_order (size_type new_order) \
367 if (new_order == base::_piola_basis.degree()) return; \
368 base::_piola_basis.reset_family_index (new_order); \
369 size_type dis_nnod = space_numbering::dis_ndof (base::_piola_basis, base::_gs, base::_gs._map_dimension); \
370 size_type nnod = space_numbering::ndof (base::_piola_basis, base::_gs, base::_gs._map_dimension); \
371 base::_gs.node_ownership = distributor (nnod, base::comm(), nnod); \
372 disarray<point_basic<T>, M> new_node (base::_gs.node_ownership); \
373 for (size_type iv = 0, nv = base::_gs.ownership_by_dimension[0].size(); iv < nv; iv++) { \
374 new_node [iv] = base::_node [iv]; \
376 base::_node = new_node; \
377 build_external_entities (); \
380 #ifdef _RHEOLEF_HAVE_MPI
383 #undef _RHEOLEF_reset_order
387 template <
class T,
class M>
391 if (_name ==
"*nogeo*" || _name ==
"unnamed") {
394 #ifdef DEBUG_LOADED_TABLE_TRACE
406 template <
class T,
class M>
410 if (_serial_number == 0)
return _name;
412 sprintf (buffer,
"%.3d",
int(_serial_number));
413 return _name +
"-" + buffer;
415 template <
class T,
class M>
420 iter = _domains.begin(), last = _domains.end(); iter != last; iter++) {
422 if (
name == dom.name())
return true;
426 template <
class T,
class M>
431 iter = _domains.begin(), last = _domains.end(); iter != last; iter++) {
433 if (
name == dom.name())
return dom;
435 error_macro (
"undefined domain \""<<
name<<
"\" in mesh \"" << _name <<
"\"");
436 return *(_domains.begin());
438 template <
class T,
class M>
443 iter = _domains.begin(), last = _domains.end(); iter != last; iter++) {
445 if (dom.name() == curr_dom.name()) {
452 _domains.push_back (dom);
454 template <
class T,
class M>
461 sz += _geo_element [
variant].size();
465 template <
class T,
class M>
472 sz += _geo_element [
variant].dis_size();
476 template <
class T,
class M>
480 std::vector<size_type>&
dis_inod)
const
484 template <
class T,
class M>
489 for (
size_type j = 0; j < _dimension; j++) {
490 _xmin[j] = std::numeric_limits<T>::max();
491 _xmax[j] = -std::numeric_limits<T>::max();
493 for (
size_type j = _dimension+1; j < 3; j++) {
494 _xmin [j] = _xmax [j] =
T(0);
498 for (
size_type j = 0 ; j < _dimension; j++) {
499 _xmin[j] = min(x[j], _xmin[j]);
500 _xmax[j] = max(x[j], _xmax[j]);
504 #ifdef _RHEOLEF_HAVE_MPI
506 for (
size_type j = 0 ; j < _dimension; j++) {
507 _xmin[j] = mpi::all_reduce (comm(), _xmin[j], mpi::minimum<T>());
508 _xmax[j] = mpi::all_reduce (comm(), _xmax[j], mpi::maximum<T>());
514 _hmin = std::numeric_limits<T>::max();
516 for (
size_t iedg = 0, nedg = geo_element_ownership(1).size(); iedg < nedg; ++iedg) {
520 T hloc =
dist(x0,x1);
521 _hmin = min(_hmin, hloc);
522 _hmax = max(_hmax, hloc);
524 #ifdef _RHEOLEF_HAVE_MPI
526 _hmin = mpi::all_reduce (comm(), _hmin, mpi::minimum<T>());
527 _hmax = mpi::all_reduce (comm(), _hmax, mpi::maximum<T>());
531 template <
class T,
class M>
538 last_ige += _geo_element [
variant].size();
539 if (ige < last_ige)
return _geo_element [
variant] [ige - first_ige];
540 first_ige = last_ige;
542 error_macro (
"geo_element index " << ige <<
" out of range [0:"<< last_ige <<
"[");
543 return _geo_element [0][0];
545 template <
class T,
class M>
552 last_ige += _geo_element [
variant].size();
553 if (ige < last_ige)
return _geo_element [
variant] [ige - first_ige];
554 first_ige = last_ige;
556 error_macro (
"geo_element index " << ige <<
" out of range [0:"<< last_ige <<
"[");
557 return _geo_element [0][0];
563 template <
class T,
class M>
568 return _gs.dis_inod2dis_iv (
dis_inod);
570 template <
class T,
class M>
575 return _gs.dis_iv2dis_inod (dis_iv);
580 template <
class T,
class M>
585 if (_gs.ownership_by_dimension[
dim].is_owned (dis_ige)) {
586 size_type first_dis_ige = _gs.ownership_by_dimension[
dim].first_index();
588 return get_geo_element (
dim, ige);
599 template <
class T,
class M>
603 if (omega.have_domain_indirect (
"boundary"))
return;
604 omega.neighbour_guard();
608 std::vector<size_type> isid_list;
609 for (
size_type isid = 0, nsid = omega.size(sid_dim); isid < nsid; isid++) {
610 const geo_element& S = omega.get_geo_element (sid_dim, isid);
611 if (S.
master(1) != std::numeric_limits<size_type>::max())
continue;
612 isid_list.push_back (isid);
614 communicator comm = omega.sizes().ownership_by_dimension[sid_dim].comm();
616 omega.insert_domain_indirect (isid_dom);
618 template <
class T,
class M>
622 if (omega.have_domain_indirect (
"internal_sides"))
return;
623 omega.neighbour_guard();
627 std::vector<size_type> isid_list;
628 for (
size_type isid = 0, nsid = omega.size(sid_dim); isid < nsid; isid++) {
629 const geo_element& S = omega.get_geo_element (sid_dim, isid);
630 if (S.
master(1) == std::numeric_limits<size_type>::max())
continue;
631 isid_list.push_back (isid);
633 communicator comm = omega.sizes().ownership_by_dimension[sid_dim].comm();
635 isid_dom.set_broken (
true);
636 omega.insert_domain_indirect (isid_dom);
638 template <
class T,
class M>
642 if (omega.have_domain_indirect (
"sides"))
return;
643 omega.neighbour_guard();
647 std::vector<size_type> isid_list;
648 for (
size_type isid = 0, nsid = omega.size(sid_dim); isid < nsid; isid++) {
650 isid_list.push_back (isid);
652 communicator comm = omega.sizes().ownership_by_dimension[sid_dim].comm();
654 isid_dom.set_broken (
true);
655 omega.insert_domain_indirect (isid_dom);
660 #define _RHEOLEF_instanciation(T,M) \
661 template class geo_base_rep<T,M>; \
662 template class geo_rep<T,M>; \
663 template class geo_basic<T,M>; \
664 template geo_basic<T,M> geo_load (const std::string& name); \
665 template geo_basic<T,M> compact (const geo_basic<T,M>&); \
666 template void boundary_guard (const geo_basic<T,M>&); \
667 template void internal_sides_guard (const geo_basic<T,M>&); \
668 template void sides_guard (const geo_basic<T,M>&);
671 #ifdef _RHEOLEF_HAVE_MPI
field::size_type size_type
see the Float page for the full documentation
abstract base interface class
std::map< size_type, geo_element_auto<> > geo_element_map_type
geo_element_hack::size_type size_type
base class for M=sequential or distributed meshes representations
size_type dis_size() const
void dis_inod(const geo_element &K, std::vector< size_type > &dis_inod) const
size_type dis_iv2dis_inod(size_type dis_iv) const
void insert_domain_indirect(const domain_indirect_basic< M > &dom) const
static loaded_map_t _loaded_map
base::size_type size_type
std::unordered_map< std::string, void * > loaded_map_t
const_reference get_geo_element(size_type dim, size_type ige) const
base::const_reference const_reference
const_reference dis_get_geo_element(size_type dim, size_type dis_ige) const
size_type dis_inod2dis_iv(size_type dis_inod) const
static loaded_map_t & loaded_map()
base::const_iterator const_iterator
const domain_indirect_basic< M > & get_domain_indirect(size_type i) const
base::reference reference
bool have_domain_indirect(const std::string &name) const
generic mesh with rerefence counting
see the geo_element page for the full documentation
size_type master(bool i) const
sequential mesh representation
static variant_type last_variant_by_dimension(size_type dim)
static variant_type first_variant_by_dimension(size_type dim)
static iorheo::force_initialization dummy
#define dis_warning_macro(message)
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
#define _RHEOLEF_set_name(M)
#define _RHEOLEF_build_by_subdividing(M)
#define _RHEOLEF_set_dimension(M)
#define _RHEOLEF_set_coordinate_system(M)
#define _RHEOLEF_build_from_data(M)
#define _RHEOLEF_reset_order(M)
#define _RHEOLEF_set_serial_number(M)
#define _RHEOLEF_set_nodes(M)
void dis_idof(const basis_basic< T > &b, const geo_size &gs, const geo_element &K, typename std::vector< size_type >::iterator dis_idof_tab)
void dis_inod(const basis_basic< T > &b, const geo_size &gs, const geo_element &K, typename std::vector< size_type >::iterator dis_inod_tab)
size_type nnod(const basis_basic< T > &b, const geo_size &gs, size_type map_dim)
This file is part of Rheolef.
T dist(const point_basic< T > &x, const point_basic< T > &y)
string delete_suffix(const string &name, const string &suffix)
delete_suffix: see the rheostream page for the full documentation
string get_basename(const string &name)
get_basename: see the rheostream page for the full documentation
void boundary_guard(const geo_basic< T, M > &omega)
_RHEOLEF_instanciation(Float) _RHEOLEF_instanciation_evaluate(Float
geo_basic< T, M > geo_load(const std::string &filename)
sequential mesh with reference counting
void sides_guard(const geo_basic< T, M > &omega)
_RHEOLEF_save(sequential) _RHEOLEF_set_name(sequential) _RHEOLEF_set_serial_number(sequential) _RHEOLEF_reset_order(sequential) _RHEOLEF_set_nodes(sequential) _RHEOLEF_set_coordinate_system(sequential) _RHEOLEF_set_dimension(sequential) _RHEOLEF_build_by_subdividing(sequential) _RHEOLEF_build_from_data(sequential) _RHEOLEF_reset_order(sequential) template< class T
void internal_sides_guard(const geo_basic< T, M > &omega)
geo_basic< T, M > compact(const geo_basic< T, M > &gamma)