Reference documentation for deal.II version Git 840703ecf3 2019-10-14 17:15:20 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
tria_objects.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2006 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #include <deal.II/base/memory_consumption.h>
17 
18 #include <deal.II/grid/tria.h>
19 #include <deal.II/grid/tria_accessor.h>
20 #include <deal.II/grid/tria_iterator.h>
21 #include <deal.II/grid/tria_objects.h>
22 
23 #include <algorithm>
24 #include <functional>
25 
26 
27 
28 DEAL_II_NAMESPACE_OPEN
29 
30 namespace internal
31 {
32  namespace TriangulationImplementation
33  {
34  template <class G>
35  void
36  TriaObjects<G>::reserve_space(const unsigned int new_objects_in_pairs,
37  const unsigned int new_objects_single)
38  {
39  Assert(new_objects_in_pairs % 2 == 0, ExcInternalError());
40 
41  next_free_single = 0;
42  next_free_pair = 0;
43  reverse_order_next_free_single = false;
44 
45  // count the number of objects, of unused single objects and of
46  // unused pairs of objects
47  unsigned int n_objects = 0;
48  unsigned int n_unused_pairs = 0;
49  unsigned int n_unused_singles = 0;
50  for (unsigned int i = 0; i < used.size(); ++i)
51  {
52  if (used[i])
53  ++n_objects;
54  else if (i + 1 < used.size())
55  {
56  if (used[i + 1])
57  {
58  ++n_unused_singles;
59  if (next_free_single == 0)
60  next_free_single = i;
61  }
62  else
63  {
64  ++n_unused_pairs;
65  if (next_free_pair == 0)
66  next_free_pair = i;
67  ++i;
68  }
69  }
70  else
71  ++n_unused_singles;
72  }
73  Assert(n_objects + 2 * n_unused_pairs + n_unused_singles == used.size(),
75 
76  // how many single objects are needed in addition to
77  // n_unused_objects?
78  const int additional_single_objects =
79  new_objects_single - n_unused_singles;
80 
81  unsigned int new_size =
82  used.size() + new_objects_in_pairs - 2 * n_unused_pairs;
83  if (additional_single_objects > 0)
84  new_size += additional_single_objects;
85 
86  // only allocate space if necessary
87  if (new_size > cells.size())
88  {
89  cells.reserve(new_size);
90  cells.insert(cells.end(), new_size - cells.size(), G());
91 
92  used.reserve(new_size);
93  used.insert(used.end(), new_size - used.size(), false);
94 
95  user_flags.reserve(new_size);
96  user_flags.insert(user_flags.end(),
97  new_size - user_flags.size(),
98  false);
99 
100  const unsigned int factor =
102  children.reserve(factor * new_size);
103  children.insert(children.end(),
104  factor * new_size - children.size(),
105  -1);
106 
107  if (G::dimension > 1)
108  {
109  refinement_cases.reserve(new_size);
110  refinement_cases.insert(
111  refinement_cases.end(),
112  new_size - refinement_cases.size(),
114  }
115 
116  // first reserve, then resize. Otherwise the std library can decide to
117  // allocate more entries.
118  boundary_or_material_id.reserve(new_size);
119  boundary_or_material_id.resize(new_size);
120 
121  user_data.reserve(new_size);
122  user_data.resize(new_size);
123 
124  manifold_id.reserve(new_size);
125  manifold_id.insert(manifold_id.end(),
126  new_size - manifold_id.size(),
128  }
129 
130  if (n_unused_singles == 0)
131  {
132  next_free_single = new_size - 1;
133  reverse_order_next_free_single = true;
134  }
135  }
136 
137 
138  template <>
139  template <int dim, int spacedim>
140  typename ::Triangulation<dim, spacedim>::raw_hex_iterator
141  TriaObjects<TriaObject<3>>::next_free_hex(
142  const ::Triangulation<dim, spacedim> &tria,
143  const unsigned int level)
144  {
145  // TODO: Think of a way to ensure that we are using the correct
146  // triangulation, i.e. the one containing *this.
147 
148  int pos = next_free_pair, last = used.size() - 1;
149  for (; pos < last; ++pos)
150  if (!used[pos])
151  {
152  // this should be a pair slot
153  Assert(!used[pos + 1], ExcInternalError());
154  break;
155  }
156  if (pos >= last)
157  // no free slot
158  return tria.end_hex();
159  else
160  next_free_pair = pos + 2;
161 
162  return
163  typename ::Triangulation<dim, spacedim>::raw_hex_iterator(&tria,
164  level,
165  pos);
166  }
167 
168 
169  void
170  TriaObjectsHex::reserve_space(const unsigned int new_hexes)
171  {
172  const unsigned int new_size =
173  new_hexes + std::count(used.begin(), used.end(), true);
174 
175  // see above...
176  if (new_size > cells.size())
177  {
178  cells.reserve(new_size);
179  cells.insert(cells.end(), new_size - cells.size(), TriaObject<3>());
180 
181  used.reserve(new_size);
182  used.insert(used.end(), new_size - used.size(), false);
183 
184  user_flags.reserve(new_size);
185  user_flags.insert(user_flags.end(),
186  new_size - user_flags.size(),
187  false);
188 
189  children.reserve(4 * new_size);
190  children.insert(children.end(), 4 * new_size - children.size(), -1);
191 
192  // for the following fields, we know exactly how many elements
193  // we need, so first reserve then resize (resize itself, at least
194  // with some compiler libraries, appears to round up the size it
195  // actually reserves)
196  boundary_or_material_id.reserve(new_size);
197  boundary_or_material_id.resize(new_size);
198 
199  manifold_id.reserve(new_size);
200  manifold_id.insert(manifold_id.end(),
201  new_size - manifold_id.size(),
203 
204  user_data.reserve(new_size);
205  user_data.resize(new_size);
206 
209  new_size * GeometryInfo<3>::faces_per_cell -
210  face_orientations.size(),
211  true);
212 
213  refinement_cases.reserve(new_size);
214  refinement_cases.insert(refinement_cases.end(),
215  new_size - refinement_cases.size(),
217 
218  face_flips.reserve(new_size * GeometryInfo<3>::faces_per_cell);
219  face_flips.insert(face_flips.end(),
220  new_size * GeometryInfo<3>::faces_per_cell -
221  face_flips.size(),
222  false);
223  face_rotations.reserve(new_size * GeometryInfo<3>::faces_per_cell);
224  face_rotations.insert(face_rotations.end(),
225  new_size * GeometryInfo<3>::faces_per_cell -
226  face_rotations.size(),
227  false);
228  }
230  }
231 
232 
233  void
234  TriaObjectsQuad3D::reserve_space(const unsigned int new_quads_in_pairs,
235  const unsigned int new_quads_single)
236  {
237  Assert(new_quads_in_pairs % 2 == 0, ExcInternalError());
238 
239  next_free_single = 0;
240  next_free_pair = 0;
242 
243  // count the number of objects, of unused single objects and of
244  // unused pairs of objects
245  unsigned int n_quads = 0;
246  unsigned int n_unused_pairs = 0;
247  unsigned int n_unused_singles = 0;
248  for (unsigned int i = 0; i < used.size(); ++i)
249  {
250  if (used[i])
251  ++n_quads;
252  else if (i + 1 < used.size())
253  {
254  if (used[i + 1])
255  {
256  ++n_unused_singles;
257  if (next_free_single == 0)
258  next_free_single = i;
259  }
260  else
261  {
262  ++n_unused_pairs;
263  if (next_free_pair == 0)
264  next_free_pair = i;
265  ++i;
266  }
267  }
268  else
269  ++n_unused_singles;
270  }
271  Assert(n_quads + 2 * n_unused_pairs + n_unused_singles == used.size(),
272  ExcInternalError());
273 
274  // how many single quads are needed in addition to n_unused_quads?
275  const int additional_single_quads = new_quads_single - n_unused_singles;
276 
277  unsigned int new_size =
278  used.size() + new_quads_in_pairs - 2 * n_unused_pairs;
279  if (additional_single_quads > 0)
280  new_size += additional_single_quads;
281 
282  // see above...
283  if (new_size > cells.size())
284  {
285  // reseve space for the base class
286  TriaObjects<TriaObject<2>>::reserve_space(new_quads_in_pairs,
287  new_quads_single);
288  // reserve the field of the derived class
289  line_orientations.reserve(new_size * GeometryInfo<2>::lines_per_cell);
290  line_orientations.insert(line_orientations.end(),
291  new_size * GeometryInfo<2>::lines_per_cell -
292  line_orientations.size(),
293  true);
294  }
295 
296  if (n_unused_singles == 0)
297  {
298  next_free_single = new_size - 1;
300  }
301  }
302 
303 
304  template <>
305  void
306  TriaObjects<TriaObject<1>>::monitor_memory(const unsigned int) const
307  {
308  Assert(cells.size() == used.size(),
309  ExcMemoryInexact(cells.size(), used.size()));
310  Assert(cells.size() == user_flags.size(),
311  ExcMemoryInexact(cells.size(), user_flags.size()));
312  Assert(cells.size() == children.size(),
313  ExcMemoryInexact(cells.size(), children.size()));
314  Assert(cells.size() == boundary_or_material_id.size(),
316  Assert(cells.size() == manifold_id.size(),
317  ExcMemoryInexact(cells.size(), manifold_id.size()));
318  Assert(cells.size() == user_data.size(),
319  ExcMemoryInexact(cells.size(), user_data.size()));
320  }
321 
322 
323  template <>
324  void
325  TriaObjects<TriaObject<2>>::monitor_memory(const unsigned int) const
326  {
327  Assert(cells.size() == used.size(),
328  ExcMemoryInexact(cells.size(), used.size()));
329  Assert(cells.size() == user_flags.size(),
330  ExcMemoryInexact(cells.size(), user_flags.size()));
331  Assert(2 * cells.size() == children.size(),
332  ExcMemoryInexact(cells.size(), children.size()));
333  Assert(cells.size() == refinement_cases.size(),
334  ExcMemoryInexact(cells.size(), refinement_cases.size()));
335  Assert(cells.size() == boundary_or_material_id.size(),
337  Assert(cells.size() == manifold_id.size(),
338  ExcMemoryInexact(cells.size(), manifold_id.size()));
339  Assert(cells.size() == user_data.size(),
340  ExcMemoryInexact(cells.size(), user_data.size()));
341  }
342 
343 
344  void
345  TriaObjectsHex::monitor_memory(const unsigned int) const
346  {
347  Assert(cells.size() == used.size(),
348  ExcMemoryInexact(cells.size(), used.size()));
349  Assert(cells.size() == user_flags.size(),
350  ExcMemoryInexact(cells.size(), user_flags.size()));
351  Assert(4 * cells.size() == children.size(),
352  ExcMemoryInexact(cells.size(), children.size()));
353  Assert(cells.size() == boundary_or_material_id.size(),
355  Assert(cells.size() == manifold_id.size(),
356  ExcMemoryInexact(cells.size(), manifold_id.size()));
357  Assert(cells.size() == user_data.size(),
358  ExcMemoryInexact(cells.size(), user_data.size()));
360  face_orientations.size(),
362  face_orientations.size()));
364  face_flips.size(),
366  face_flips.size()));
368  face_rotations.size(),
370  face_rotations.size()));
371  }
372 
373 
374  void
375  TriaObjectsQuad3D::monitor_memory(const unsigned int) const
376  {
377  // check that we have not allocated too much memory. note that bool
378  // vectors allocate their memory in chunks of whole integers, so they
379  // may over-allocate by up to as many elements as an integer has bits
381  line_orientations.size(),
383  line_orientations.size()));
385  }
386 
387 
388  template <typename G>
389  void
391  {
392  cells.clear();
393  children.clear();
394  refinement_cases.clear();
395  used.clear();
396  user_flags.clear();
397  boundary_or_material_id.clear();
398  manifold_id.clear();
399  user_data.clear();
401  }
402 
403 
404  void
406  {
408  face_orientations.clear();
409  face_flips.clear();
410  face_rotations.clear();
411  }
412 
413 
414  void
416  {
418  line_orientations.clear();
419  }
420 
421 
422  template <typename G>
423  std::size_t
425  {
433  user_data.capacity() * sizeof(UserData) + sizeof(user_data));
434  }
435 
436 
437  std::size_t
439  {
444  }
445 
446 
447  std::size_t
449  {
450  return (MemoryConsumption::memory_consumption(line_orientations) +
452  }
453 
454 
455 
456  // explicit instantiations
457  template class TriaObjects<TriaObject<1>>;
458  template class TriaObjects<TriaObject<2>>;
459 
460 #include "tria_objects.inst"
461  } // namespace TriangulationImplementation
462 } // namespace internal
463 
464 DEAL_II_NAMESPACE_CLOSE
void reserve_space(const unsigned int new_objs)
const types::manifold_id flat_manifold_id
Definition: types.h:267
static ::ExceptionBase & ExcMemoryInexact(int arg1, int arg2)
void monitor_memory(const unsigned int true_dimension) const
void monitor_memory(const unsigned int true_dimension) const
#define Assert(cond, exc)
Definition: exceptions.h:1407
void reserve_space(const unsigned int new_objs_in_pairs, const unsigned int new_objs_single=0)
Definition: tria_objects.cc:36
std::vector< RefinementCase< TriaObject< 3 > ::dimension > > refinement_cases
Definition: tria_objects.h:101
void reserve_space(const unsigned int new_quads_in_pairs, const unsigned int new_quads_single=0)
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
static ::ExceptionBase & ExcInternalError()