GetFEM  5.4.2
getfem_mesh_region.h
Go to the documentation of this file.
1 /* -*- c++ -*- (enables emacs c++ mode) */
2 /*===========================================================================
3 
4  Copyright (C) 2005-2020 Yves Renard, Julien Pommier
5 
6  This file is a part of GetFEM
7 
8  GetFEM is free software; you can redistribute it and/or modify it
9  under the terms of the GNU Lesser General Public License as published
10  by the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version along with the GCC Runtime Library
12  Exception either version 3.1 or (at your option) any later version.
13  This program is distributed in the hope that it will be useful, but
14  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
15  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
16  License and GCC Runtime Library Exception for more details.
17  You should have received a copy of the GNU Lesser General Public License
18  along with this program; if not, write to the Free Software Foundation,
19  Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
20 
21  As a special exception, you may use this file as it is a part of a free
22  software library without restriction. Specifically, if other files
23  instantiate templates or use macros or inline functions from this file,
24  or you compile this file and link it with other files to produce an
25  executable, this file does not by itself cause the resulting executable
26  to be covered by the GNU Lesser General Public License. This exception
27  does not however invalidate any other reasons why the executable file
28  might be covered by the GNU Lesser General Public License.
29 
30 ===========================================================================*/
31 
32 /**@file getfem_mesh_region.h
33 @author Yves Renard, Julien Pommier
34 @date 2005.
35 @brief region objects (set of convexes and/or convex faces)
36 */
37 
38 #pragma once
39 
40 #include <atomic>
41 #include <bitset>
42 #include <iostream>
43 #include <map>
44 
45 #include "dal_bit_vector.h"
46 #include "bgeot_convex_structure.h"
47 #include "getfem_config.h"
48 
49 namespace getfem {
50  class mesh;
51 
52  /** structure used to hold a set of convexes and/or convex faces.
53  @see mesh::region
54  */
55  class APIDECL mesh_region {
56  public:
57  using face_bitset = std::bitset<MAX_FACES_PER_CV+1>;
58  using map_t = std::map<size_type, face_bitset>;
59 
60  private:
61 
62  using const_iterator = map_t::const_iterator;
63 
64  struct impl {
65  mutable map_t m;
66  mutable omp_distribute<dal::bit_vector> index_;
67  mutable dal::bit_vector serial_index_;
68  };
69  std::shared_ptr<impl> p; /* the real region data */
70 
71  size_type id_; /* used temporarily when the
72  mesh_region(size_type)
73  constructor is used */
74 
75  size_type type_; //optional type of the region
76 
77  std::atomic_bool partitioning_allowed; /* specifies that in multithreaded code only a
78  partition of the region is visible in index()
79  and size() methods, as well as during iteration
80  with mr_visitor */
81 
82  mesh *parent_mesh; /* used for mesh_region "extracted" from
83  a mesh (to provide feedback) */
84 
85 
86  //cashing iterators for partitions
87  mutable omp_distribute<const_iterator> itbegin;
88  mutable omp_distribute<const_iterator> itend;
89 
90  //flags for all the cashes
91  mutable omp_distribute<bool> index_updated;
92  mutable omp_distribute<bool> partitions_updated;
93  mutable bool serial_index_updated;
94 
95  void mark_region_changed() const;
96 
97  void update_index() const;
98 
99  void update_partition_iterators() const;
100 
101  impl &wp() { return *p.get(); }
102  const impl &rp() const { return *p.get(); }
103  void clean();
104  /** tells the owner mesh that the region is valid */
105  void touch_parent_mesh();
106 
107  /**when running while multithreaded, gives the iterator
108  for the beginning of the region partition for the current thread*/
109  const_iterator partition_begin() const;
110 
111  /**when running while multithreaded, gives the iterator
112  for the end of the region partition for the current thread*/
113  const_iterator partition_end() const;
114 
115  /**begin iterator of the region depending if its partitioned or not*/
116  const_iterator begin() const;
117 
118  /**end iterator of the region depending if its partitioned or not*/
119  const_iterator end() const;
120 
121  /**number of region entries before partitioning*/
122  size_type unpartitioned_size() const;
123 
124  public:
125  mesh_region(const mesh_region &other);
126  mesh_region();
127  /** a mesh_region can be built from a integer parameter
128  (a region number in a mesh),
129  but it won't be usable until 'from_mesh(m)' has been called
130  Note that these regions are read-only, this constructor is
131  mostly used for backward-compatibility.
132  */
133  mesh_region(size_type id__);
134 
135  /** internal constructor. You should used m.region(id) instead. */
136  mesh_region(mesh& m, size_type id__, size_type type = size_type(-1));
137  /** build a mesh_region from a convex list stored in a bit_vector. */
138  mesh_region(const dal::bit_vector &bv);
139 
140  /** provide a default value for the mesh_region parameters of assembly
141  procedures etc. */
143  return mesh_region(size_type(-1));
144  }
145  /** return the intersection of two mesh regions */
146  static mesh_region intersection(const mesh_region& a,
147  const mesh_region& b);
148  /** return the union of two mesh_regions */
149  static mesh_region merge(const mesh_region &a,
150  const mesh_region &b);
151  /** subtract the second region from the first one */
152  static mesh_region subtract(const mesh_region &a,
153  const mesh_region &b);
154  /** Test if the region is a boundary of a list of faces of elements of
155  region `rg`. Return 0 if not, -1 if only partially, 1 if the region
156  contains only some faces which are all faces of elements of `rg`. */
157  int region_is_faces_of(const getfem::mesh& m1,
158  const mesh_region &rg2,
159  const getfem::mesh& m2) const;
160 
161  size_type id() const { return id_; }
162 
163  size_type get_type() const { return type_; }
164 
165  void set_type(size_type type) { type_ = type; }
166 
167  /** In multithreaded part of the program makes only a partition of the
168  region visible in the index() and size() operations, as well as during
169  iterations with mr_visitor. This is a default behaviour*/
170  void allow_partitioning();
171 
172  /** Return the bounding box [Pmin - Pmax] of the mesh_region. */
173  void bounding_box(base_node& Pmin, base_node& Pmax) const;
174 
175  /** Disregard partitioning, which means being able to see the whole region
176  in multithreaded code. Can be used, for instance, for contact problems
177  where master region is partitioned, while the slave region is not*/
178  void prohibit_partitioning();
179 
180  bool is_partitioning_allowed() const;
181 
182  /** Extract the next region number
183  that does not yet exists in the mesh*/
184  static size_type free_region_id(const getfem::mesh& m);
185 
186  /** For regions which have been built with just a number 'id',
187  from_mesh(m) sets the current region to 'm.region(id)'.
188  (works only once)
189  */
190  const mesh_region& from_mesh(const mesh &m) const;
191 
192  mesh_region& operator=(const mesh_region &mr);
193 
194  bool compare(const mesh &m1, const mesh_region &mr, const mesh &m2) const;
195 
196  face_bitset operator[](size_t cv) const;
197 
198  /** Index of the region convexes, or the convexes from the partition on the
199  current thread. */
200  const dal::bit_vector& index() const;
201  void add(const dal::bit_vector &bv);
202  void add(size_type cv, short_type f = short_type(-1));
203  void sup(size_type cv, short_type f = short_type(-1));
204  void sup_all(size_type cv);
205  void clear();
206  void swap_convex(size_type cv1, size_type cv2);
207  bool is_in(size_type cv, short_type f = short_type(-1)) const;
208  bool is_in(size_type cv, short_type f, const mesh &m) const;
209 
210  /** Region size, or the size of the region partition on the current
211  thread if the region is partitioned*/
212  size_type size() const;
213 
214  /** Number of convexes in the region, or on the partition on the current
215  thread*/
216  size_type nb_convex() const { return index().card();}
217  bool is_empty() const;
218  /** Return true if the region do contain only convex faces */
219  bool is_only_faces() const;
220  bool is_boundary() const { return is_only_faces(); }
221  /** Return true if the region do not contain any convex face */
222  bool is_only_convexes() const;
223  face_bitset faces_of_convex(size_type cv) const;
224  face_bitset and_mask() const;
225  face_bitset or_mask() const;
226  void error_if_not_faces() const;
227  void error_if_not_convexes() const;
228  void error_if_not_homogeneous() const;
229  const mesh *get_parent_mesh(void) const { return parent_mesh; }
230  void set_parent_mesh(mesh *pm) { parent_mesh = pm; }
231 
232  /** "iterator" class for regions. Usage similar to bv_visitor:
233  for (mr_visitor i(region); !i.finished(); ++i) {
234  ...
235  }
236  */
237  class visitor {
238 
239  typedef mesh_region::map_t::const_iterator const_iterator;
240  bool whole_mesh;
241  dal::bit_const_iterator itb, iteb;
242  const_iterator it, ite;
243  face_bitset c;
244  size_type cv_;
245  short_type f_;
246  bool finished_;
247 #if GETFEM_PARA_LEVEL > 1
248  std::unique_ptr<mesh_region> mpi_rg;
249 #endif
250  void init(const mesh_region &s);
251  void init(const dal::bit_vector &s);
252 
253  public:
254  visitor(const mesh_region &s);
255  visitor(const mesh_region &s, const mesh &m,
256  bool intersect_with_mpi = false);
257  size_type cv() const { return cv_; }
258  size_type is_face() const { return f_ != 0; }
259  short_type f() const { return short_type(f_-1); }
260 
261  bool next();
262  bool operator++() { return next(); }
263 
264  bool finished() const { return finished_; }
265 
266  bool next_face(){
267  if (whole_mesh) return false;
268  if (c.none()) return false;
269  do { ++f_; } while (!c.test(f_));
270  c.set(f_,0);
271  return true;
272  }
273  };
274 
275  friend std::ostream & operator <<(std::ostream &os, const mesh_region &w);
276  };
277 
279 
280  /** Dummy mesh_region for default parameter of functions. */
282 
283 } /* end of namespace getfem. */
getfem_config.h
defines and typedefs for namespace getfem
bgeot::size_type
size_t size_type
used as the common size type in the library
Definition: bgeot_poly.h:49
gmm::clear
void clear(L &l)
clear (fill with zeros) a vector or matrix.
Definition: gmm_blas.h:59
getfem::dummy_mesh_region
const mesh_region & dummy_mesh_region()
Dummy mesh_region for default parameter of functions.
Definition: getfem_mesh_region.cc:617
getfem::mesh_region::all_convexes
static mesh_region all_convexes()
provide a default value for the mesh_region parameters of assembly procedures etc.
Definition: getfem_mesh_region.h:142
bgeot::short_type
gmm::uint16_type short_type
used as the common short type integer in the library
Definition: bgeot_config.h:72
getfem
GEneric Tool for Finite Element Methods.
Definition: getfem_accumulated_distro.h:46
dal_bit_vector.h
Provide a dynamic bit container.
getfem::mesh_region::visitor
"iterator" class for regions.
Definition: getfem_mesh_region.h:237
getfem::omp_distribute< dal::bit_vector >
bgeot_convex_structure.h
Definition of convex structures.
getfem::mesh_region::nb_convex
size_type nb_convex() const
Number of convexes in the region, or on the partition on the current thread.
Definition: getfem_mesh_region.h:216
getfem::mesh
Describe a mesh (collection of convexes (elements) and points).
Definition: getfem_mesh.h:95
getfem::mesh_region
structure used to hold a set of convexes and/or convex faces.
Definition: getfem_mesh_region.h:55
gmm::add
void add(const L1 &l1, L2 &l2)
*‍/
Definition: gmm_blas.h:1268