GetFEM  5.4.2
getfem_im_data.cc
1 /*===========================================================================
2 
3  Copyright (C) 2012-2020 Liang Jin Lim.
4 
5  This file is a part of GetFEM
6 
7  GetFEM is free software; you can redistribute it and/or modify it
8  under the terms of the GNU Lesser General Public License as published
9  by the Free Software Foundation; either version 3 of the License, or
10  (at your option) any later version along with the GCC Runtime Library
11  Exception either version 3.1 or (at your option) any later version.
12  This program is distributed in the hope that it will be useful, but
13  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14  or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public
15  License and GCC Runtime Library Exception for more details.
16  You should have received a copy of the GNU Lesser General Public License
17  along with this program; if not, write to the Free Software Foundation,
18  Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA.
19 
20 ===========================================================================*/
21 
22 #include "getfem/getfem_im_data.h"
23 #include "getfem/getfem_omp.h"
24 
25 namespace getfem
26 {
27 
29  bgeot::multi_index tsize,
30  size_type filtered_region_,
31  bgeot::multi_index actual_tsize)
32  : im_(mim), region_(filtered_region_),
33  nb_int_pts_intern(0), nb_int_pts_onfaces(0),
34  nb_filtered_int_pts_intern(0), nb_filtered_int_pts_onfaces(0),
35  locks_()
36  {
37  set_tensor_size(tsize);
38  set_actual_tensor_size(actual_tsize);
39  add_dependency(im_);
41  }
42 
43  im_data::im_data(const getfem::mesh_im& mim, size_type filtered_region_)
44  : im_(mim), region_(filtered_region_),
45  nb_int_pts_intern(0), nb_int_pts_onfaces(0),
46  nb_filtered_int_pts_intern(0), nb_filtered_int_pts_onfaces(0),
47  locks_()
48  {
49  tensor_size_.resize(1);
50  tensor_size_[0] = 1;
51  actual_tensor_size_ = tensor_size_;
52  nb_tensor_elem_ = 1;
53  add_dependency(im_);
55  }
56 
57  size_type im_data::nb_index(bool use_filter) const {
58  context_check();
59  if (use_filter)
60  return nb_filtered_int_pts_intern + nb_filtered_int_pts_onfaces;
61  else
62  return nb_int_pts_intern + nb_int_pts_onfaces;
63  }
64 
65  size_type im_data::nb_points_of_element(size_type cv, bool use_filter) const {
66  context_check();
67  if (cv < convexes.size()) {
68  size_type nb_int_pts(0);
69  if (use_filter) {
70  for (short_type f=0, nb_faces=nb_faces_of_element(cv);
71  f < nb_faces; ++f)
72  if (convexes[cv].first_int_pt_onface_fid[f] != size_type(-1))
73  nb_int_pts += convexes[cv].nb_int_pts_onface[f];
74  if (convexes[cv].first_int_pt_fid != size_type(-1)) // convex in filtered_region()
75  nb_int_pts += convexes[cv].nb_int_pts;
76  }
77  else {
78  for (auto nb_pts : convexes[cv].nb_int_pts_onface)
79  nb_int_pts += nb_pts;
80  if (nb_int_pts_intern > 0) // has_convexes
81  nb_int_pts += convexes[cv].nb_int_pts;
82  }
83  return nb_int_pts;
84  }
85  return 0;
86  }
87 
89  bool use_filter) const {
90  context_check();
91  if (cv < convexes.size()) {
92  if (f == short_type(-1)) {
93  if (!use_filter || convexes[cv].first_int_pt_fid != size_type(-1))
94  return convexes[cv].nb_int_pts;
95  }
96  else if (f < convexes[cv].nb_int_pts_onface.size()) {
97  if (!use_filter || convexes[cv].first_int_pt_onface_fid[f] != size_type(-1))
98  return convexes[cv].nb_int_pts_onface[f];
99  }
100  }
101  return 0;
102  }
103 
105  context_check();
106  if (cv < convexes.size())
107  return short_type(convexes[cv].first_int_pt_onface_id.size());
108  return 0;
109  }
110 
112  bool use_filter) const {
113  context_check();
114  if (cv < convexes.size()) {
115  if (i < convexes[cv].nb_int_pts) { // internal point
116  size_type int_pt_id = use_filter ? convexes[cv].first_int_pt_fid
117  : convexes[cv].first_int_pt_id;
118  if (int_pt_id != size_type(-1))
119  return int_pt_id + i;
120  }
121  else if (nb_faces_of_element(cv) != 0) {
122  const getfem::papprox_integration pim = approx_int_method_of_element(cv);
123  for (short_type f=0, nb_faces=pim->nb_convex_faces();
124  f < nb_faces; ++f)
125  if (i < pim->repart()[f+1]) {
126  size_type int_pt_id = use_filter ? convexes[cv].first_int_pt_onface_fid[f]
127  : convexes[cv].first_int_pt_onface_id[f];
128  if (int_pt_id != size_type(-1))
129  return int_pt_id + i - pim->ind_first_point_on_face(f);
130  break;
131  }
132  }
133  }
134  return size_type(-1);
135  }
136 
138  { return index_of_point(cv, i, true); }
139 
141  bool use_filter) const {
142  context_check();
143  if (cv < convexes.size()) {
144  if (f == short_type(-1)) {
145  return use_filter ? convexes[cv].first_int_pt_fid
146  : convexes[cv].first_int_pt_id;
147  }
148  else if (f < nb_faces_of_element(cv)) {
149  return use_filter ? convexes[cv].first_int_pt_onface_fid[f]
150  : convexes[cv].first_int_pt_onface_id[f];
151  }
152  }
153  return size_type(-1);
154  }
155 
157  { return index_of_first_point(cv, f, true); }
158 
159  dal::bit_vector im_data::convex_index(bool use_filter) const {
160  context_check();
161  dal::bit_vector ind = im_.convex_index();
162  if (use_filter && filtered_region() != size_type(-1))
163  ind &= linked_mesh().region(filtered_region()).index();
164  return ind;
165  }
166 
168  region_ = rg;
170  }
171 
173  local_guard lock = locks_.get_lock();
174  bool no_region = (filtered_region() == size_type(-1));
175  const getfem::mesh_region &rg = no_region
178  bool has_convexes = (no_region || !rg.is_only_faces());
179  bool has_faces = (!no_region && !rg.is_only_convexes());
180 
181  size_type nb_cv = im_.convex_index().last_true() + 1;
182  convexes.clear();
183  convexes.resize(nb_cv);
184 
185  for (dal::bv_visitor cv(im_.convex_index()); !cv.finished(); ++cv)
186  convexes[cv].nb_int_pts
187  = im_.int_method_of_element(cv)->approx_method()->nb_points_on_convex();
188 
189  nb_int_pts_intern = 0;
190  nb_filtered_int_pts_intern = 0;
191  if (has_convexes)
192  // The global indexing of integration points follows the indexing of convexes
193  for (dal::bv_visitor cv(im_.convex_index()); !cv.finished(); ++cv) {
194  convexes[cv].first_int_pt_id = nb_int_pts_intern;
195  nb_int_pts_intern += convexes[cv].nb_int_pts;
196  if (no_region || rg.is_in(cv)) {
197  convexes[cv].first_int_pt_fid = nb_filtered_int_pts_intern;
198  nb_filtered_int_pts_intern += convexes[cv].nb_int_pts;
199  }
200  }
201 
202  nb_int_pts_onfaces = 0;
203  nb_filtered_int_pts_onfaces = 0;
204  if (has_faces) {
205  for (dal::bv_visitor cv(im_.convex_index()); !cv.finished(); ++cv) {
206  short_type nb_faces = im_.linked_mesh().nb_faces_of_convex(cv);
207  convexes[cv].first_int_pt_onface_id.assign(nb_faces, size_type(-1));
208  convexes[cv].first_int_pt_onface_fid.assign(nb_faces, size_type(-1));
209  convexes[cv].nb_int_pts_onface.assign(nb_faces, size_type(-1));
210  const getfem::papprox_integration pim(approx_int_method_of_element(cv));
211  for (short_type f = 0; f < nb_faces; ++f) {
212  convexes[cv].first_int_pt_onface_id[f] = nb_int_pts_intern +
213  nb_int_pts_onfaces;
214  size_type nb_pts = pim->nb_points_on_face(f);
215  nb_int_pts_onfaces += nb_pts;
216  if (rg.is_in(cv, f)) {
217  convexes[cv].first_int_pt_onface_fid[f] = nb_filtered_int_pts_intern +
218  nb_filtered_int_pts_onfaces;
219  nb_filtered_int_pts_onfaces += nb_pts;
220  }
221  convexes[cv].nb_int_pts_onface[f] = nb_pts;
222  }
223  }
224  }
225  v_num_ = act_counter();
226  touch();
227  }
228 
230  return nb_tensor_elem_;
231  }
232 
233  void im_data::set_tensor_size(const bgeot::multi_index& tsize) {
234  if (actual_tensor_size_ == tensor_size_) actual_tensor_size_ = tsize;
235  tensor_size_ = tsize;
236  nb_tensor_elem_ = tensor_size_.total_size();
237  }
238 
239  void im_data::set_actual_tensor_size(const bgeot::multi_index& tsize) {
240  actual_tensor_size_ = tsize.empty() ? tensor_size_ : tsize;
241  }
242 
243  bool is_equivalent_with_vector(const bgeot::multi_index &sizes, size_type vector_size) {
244  bool checked = false;
245  size_type size = 1;
246  for (size_type i = 0; i < sizes.size(); ++i) {
247  if (sizes[i] > 1 && checked) return false;
248  if (sizes[i] > 1) {
249  checked = true;
250  size = sizes[i];
251  if (size != vector_size) return false;
252  }
253  }
254  return (vector_size == size);
255  }
256 
257  bool is_equivalent_with_matrix(const bgeot::multi_index &sizes, size_type nrows, size_type ncols) {
258  if (nrows == 1 || ncols == 1) {
259  return is_equivalent_with_vector(sizes, nrows + ncols - 1);
260  }
261  size_type tensor_row = 1;
262  size_type tensor_col = 1;
263  bool first_checked = false;
264  bool second_checked = false;
265  for (size_type i = 0; i < sizes.size(); ++i) {
266  if (sizes[i] > 1 && !first_checked) {
267  first_checked = true;
268  tensor_row = sizes[i];
269  if (tensor_row != nrows) return false;
270  } else if (sizes[i] > 1 && !second_checked) {
271  second_checked = true;
272  tensor_col = sizes[i];
273  if (tensor_col != ncols) return false;
274  }
275  else if (sizes[i] > 1 && first_checked && second_checked) return false;
276  }
277  return (nrows == tensor_row && ncols == tensor_col);
278  }
279 }
getfem::im_data::index_of_point
size_type index_of_point(size_type cv, size_type i, bool use_filter=false) const
Returns the index of an integration point with no filtering.
Definition: getfem_im_data.cc:111
getfem::im_data::linked_mesh
const mesh & linked_mesh() const
linked mesh
Definition: getfem_im_data.h:159
getfem::im_data::update_from_context
void update_from_context() const
called automatically when there is a change in dependencies
Definition: getfem_im_data.cc:172
getfem::im_data::nb_points_of_element
size_type nb_points_of_element(size_type cv, bool use_filter=false) const
Total number of points in element cv.
Definition: getfem_im_data.cc:65
getfem::mesh::region
const mesh_region region(size_type id) const
Return the region of index 'id'.
Definition: getfem_mesh.h:421
getfem::im_data::nb_tensor_elem
size_type nb_tensor_elem() const
sum of tensor elements, M(3,3) will have 3*3=9 elements
Definition: getfem_im_data.cc:229
getfem::im_data::filtered_index_of_first_point
size_type filtered_index_of_first_point(size_type cv, short_type f=short_type(-1)) const
Returns the index of the first integration point with filtering.
Definition: getfem_im_data.cc:156
getfem::mesh_im::convex_index
const dal::bit_vector & convex_index(void) const
Get the set of convexes where an integration method has been assigned.
Definition: getfem_mesh_im.h:73
bgeot::size_type
size_t size_type
used as the common size type in the library
Definition: bgeot_poly.h:49
getfem::mesh_im
Describe an integration method linked to a mesh.
Definition: getfem_mesh_im.h:47
getfem::im_data::set_region
void set_region(size_type region)
set filtered region id
Definition: getfem_im_data.cc:167
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
getfem::mesh_region::is_only_faces
bool is_only_faces() const
Return true if the region do contain only convex faces.
Definition: getfem_mesh_region.cc:350
getfem_omp.h
Tools for multithreaded, OpenMP and Boost based parallelization.
getfem::mesh_region::index
const dal::bit_vector & index() const
Index of the region convexes, or the convexes from the partition on the current thread.
Definition: getfem_mesh_region.cc:237
getfem::im_data::filtered_region
size_type filtered_region() const
return filtered region id
Definition: getfem_im_data.h:98
getfem::im_data::convex_index
dal::bit_vector convex_index(bool use_filter=false) const
List of convexes.
Definition: getfem_im_data.cc:159
getfem::im_data::index_of_first_point
size_type index_of_first_point(size_type cv, short_type f=short_type(-1), bool use_filter=false) const
Returns the index of the first integration point with no filtering.
Definition: getfem_im_data.cc:140
getfem_im_data.h
Provides indexing of integration points for mesh_im.
getfem::im_data::filtered_index_of_point
size_type filtered_index_of_point(size_type cv, size_type i) const
Returns the index of an integration point with filtering.
Definition: getfem_im_data.cc:137
getfem::mesh_im::linked_mesh
const mesh & linked_mesh() const
Give a reference to the linked mesh of type mesh.
Definition: getfem_mesh_im.h:79
getfem::mesh_region
structure used to hold a set of convexes and/or convex faces.
Definition: getfem_mesh_region.h:55
getfem::im_data::nb_faces_of_element
short_type nb_faces_of_element(size_type cv) const
Number of (active) faces in element cv.
Definition: getfem_im_data.cc:104
getfem::mesh_im::int_method_of_element
virtual pintegration_method int_method_of_element(size_type cv) const
return the integration method associated with an element (in no integration is associated,...
Definition: getfem_mesh_im.h:117
getfem::mesh_region::is_only_convexes
bool is_only_convexes() const
Return true if the region do not contain any convex face.
Definition: getfem_mesh_region.cc:345
getfem::is_equivalent_with_vector
bool is_equivalent_with_vector(const bgeot::multi_index &sizes, size_type vector_size)
check if a given tensor size is equivalent to a vector
Definition: getfem_im_data.cc:243
getfem::context_dependencies::context_check
bool context_check() const
return true if update_from_context was called
Definition: getfem_context.h:126
getfem::is_equivalent_with_matrix
bool is_equivalent_with_matrix(const bgeot::multi_index &sizes, size_type nrows, size_type ncols)
check if a given tensor size is equivalent to a matrix
Definition: getfem_im_data.cc:257
getfem::im_data::im_data
im_data(const mesh_im &mim_, bgeot::multi_index tensor_size, size_type filtered_region_=size_type(-1), bgeot::multi_index actual_tensor_size={})
Constructor.
Definition: getfem_im_data.cc:28
getfem::im_data::nb_index
size_type nb_index(bool use_filter=false) const
Total numbers of index (integration points)
Definition: getfem_im_data.cc:57
bgeot::mesh_structure::nb_faces_of_convex
short_type nb_faces_of_convex(size_type ic) const
Return the number of faces of convex ic.
Definition: bgeot_mesh_structure.h:118