Reference documentation for deal.II version 9.3.0
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
tria_objects.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2006 - 2021 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 #ifndef dealii_tria_objects_h
17 #define dealii_tria_objects_h
18 
19 #include <deal.II/base/config.h>
20 
24 
25 #include <vector>
26 
28 
29 // Forward declarations
30 #ifndef DOXYGEN
31 template <int dim, int spacedim>
32 class Triangulation;
33 template <class Accessor>
34 class TriaRawIterator;
35 template <int, int, int>
36 class TriaAccessor;
37 #endif
38 
39 namespace internal
40 {
41  namespace TriangulationImplementation
42  {
53  {
54  public:
58  TriaObjects();
59 
63  TriaObjects(const unsigned int structdim);
64 
65  unsigned int structdim;
66 
71  std::vector<int> cells;
72 
76  unsigned int
77  n_objects() const;
78 
87  get_bounding_object_indices(const unsigned int index);
88 
100  std::vector<int> children;
101 
107  std::vector<std::uint8_t> refinement_cases;
108 
117  std::vector<bool> used;
118 
128  std::vector<bool> user_flags;
129 
130 
136  {
137  union
138  {
141  };
142 
143 
148 
152  static std::size_t
154 
160  template <class Archive>
161  void
162  serialize(Archive &ar, const unsigned int version);
163  };
164 
179  std::vector<BoundaryOrMaterialId> boundary_or_material_id;
180 
185  std::vector<types::manifold_id> manifold_id;
186 
198  template <int structdim, int dim, int spacedim>
201 
213  template <int structdim, int dim, int spacedim>
216 
221  template <int dim, int spacedim>
224  const unsigned int level);
225 
229  void *&
230  user_pointer(const unsigned int i);
231 
235  const void *
236  user_pointer(const unsigned int i) const;
237 
241  unsigned int &
242  user_index(const unsigned int i);
243 
247  unsigned int
248  user_index(const unsigned int i) const;
249 
253  void
254  clear_user_data(const unsigned int i);
255 
260  void
261  clear_user_data();
262 
266  void
268 
273  std::size_t
274  memory_consumption() const;
275 
281  template <class Archive>
282  void
283  serialize(Archive &ar, const unsigned int version);
284 
293 
297  unsigned int next_free_single;
298 
302  unsigned int next_free_pair;
303 
308 
312  struct UserData
313  {
314  union
315  {
318  void *p;
321  unsigned int i;
322  };
323 
328  {
329  p = nullptr;
330  }
331 
337  template <class Archive>
338  void
339  serialize(Archive &ar, const unsigned int version);
340  };
341 
346  {
353  };
354 
355 
360  std::vector<UserData> user_data;
361 
368  };
369 
370 
371  //----------------------------------------------------------------------//
372 
373  inline unsigned int
375  {
376  // assume that each cell has the same number of faces
377  const unsigned int faces_per_cell = 2 * this->structdim;
378  return cells.size() / faces_per_cell;
379  }
380 
381 
382 
383  inline ArrayView<int>
385  {
386  // assume that each cell has the same number of faces
387  const unsigned int faces_per_cell = 2 * this->structdim;
388  return ArrayView<int>(cells.data() + index * faces_per_cell,
389  faces_per_cell);
390  }
391 
392 
393 
395  {
397  }
398 
399 
400 
401  inline std::size_t
403  {
404  return sizeof(BoundaryOrMaterialId);
405  }
406 
407 
408 
409  template <class Archive>
410  void
412  const unsigned int /*version*/)
413  {
414  // serialize this
415  // structure by
416  // writing and
417  // reading the larger
418  // of the two values,
419  // in order to make
420  // sure we get all
421  // bits
422  if (sizeof(material_id) > sizeof(boundary_id))
423  ar &material_id;
424  else
425  ar &boundary_id;
426  }
427 
428 
429  inline void *&
430  TriaObjects::user_pointer(const unsigned int i)
431  {
435 
436  AssertIndexRange(i, user_data.size());
437  return user_data[i].p;
438  }
439 
440 
441  inline const void *
442  TriaObjects::user_pointer(const unsigned int i) const
443  {
447 
448  AssertIndexRange(i, user_data.size());
449  return user_data[i].p;
450  }
451 
452 
453  inline unsigned int &
454  TriaObjects::user_index(const unsigned int i)
455  {
459 
460  AssertIndexRange(i, user_data.size());
461  return user_data[i].i;
462  }
463 
464 
465  inline void
466  TriaObjects::clear_user_data(const unsigned int i)
467  {
468  AssertIndexRange(i, user_data.size());
469  user_data[i].i = 0;
470  }
471 
472 
479  {}
480 
481 
482  inline TriaObjects::TriaObjects(const unsigned int structdim)
483  : structdim(structdim)
488  {}
489 
490 
491  inline unsigned int
492  TriaObjects::user_index(const unsigned int i) const
493  {
497 
498  AssertIndexRange(i, user_data.size());
499  return user_data[i].i;
500  }
501 
502 
503  inline void
505  {
507  for (auto &data : user_data)
508  data.p = nullptr;
509  }
510 
511 
512  inline void
514  {
515  user_flags.assign(user_flags.size(), false);
516  }
517 
518 
519  template <class Archive>
520  void
521  TriaObjects::UserData::serialize(Archive &ar, const unsigned int)
522  {
523  // serialize this as an integer
524  ar &i;
525  }
526 
527 
528 
529  template <class Archive>
530  void
531  TriaObjects::serialize(Archive &ar, const unsigned int)
532  {
533  ar &structdim;
534  ar &cells &children;
535  ar & refinement_cases;
536  ar & used;
537  ar & user_flags;
539  ar & manifold_id;
542  }
543 
544 
545  //----------------------------------------------------------------------//
546 
547  template <int structdim_, int dim, int spacedim>
550  const Triangulation<dim, spacedim> &tria)
551  {
552  // TODO: Think of a way to ensure that we are using the correct
553  // triangulation, i.e. the one containing *this.
554 
555  AssertDimension(structdim_, this->structdim);
556 
557  int pos = next_free_single, last = used.size() - 1;
559  {
560  // first sweep forward, only use really single slots, do not use
561  // pair slots
562  for (; pos < last; ++pos)
563  if (!used[pos])
564  if (used[++pos])
565  {
566  // this was a single slot
567  pos -= 1;
568  break;
569  }
570  if (pos >= last)
571  {
573  next_free_single = used.size() - 1;
574  pos = used.size() - 1;
575  }
576  else
577  next_free_single = pos + 1;
578  }
579 
581  {
582  // second sweep, use all slots, even
583  // in pairs
584  for (; pos >= 0; --pos)
585  if (!used[pos])
586  break;
587  if (pos > 0)
588  next_free_single = pos - 1;
589  else
590  // no valid single object anymore
591  return ::TriaRawIterator<
593  }
594 
595  return ::TriaRawIterator<
597  }
598 
599 
600 
601  template <int structdim_, int dim, int spacedim>
604  {
605  // TODO: Think of a way to ensure that we are using the correct
606  // triangulation, i.e. the one containing *this.
607 
608  AssertDimension(structdim_, this->structdim);
609 
610  int pos = next_free_pair, last = used.size() - 1;
611  for (; pos < last; ++pos)
612  if (!used[pos])
613  if (!used[++pos])
614  {
615  // this was a pair slot
616  pos -= 1;
617  break;
618  }
619  if (pos >= last)
620  // no free slot
621  return ::TriaRawIterator<
623  else
624  next_free_pair = pos + 2;
625 
626  return ::TriaRawIterator<
628  }
629  } // namespace TriangulationImplementation
630 } // namespace internal
631 
632 
633 
635 
636 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:196
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1622
::TriaRawIterator<::TriaAccessor< structdim, dim, spacedim > > next_free_single_object(const Triangulation< dim, spacedim > &tria)
ArrayView< int > get_bounding_object_indices(const unsigned int index)
Definition: tria_objects.h:384
unsigned int & user_index(const unsigned int i)
Definition: tria_objects.h:454
#define AssertIndexRange(index, range)
Definition: exceptions.h:1690
void serialize(Archive &ar, const unsigned int version)
Definition: tria_objects.h:521
#define Assert(cond, exc)
Definition: exceptions.h:1465
void serialize(Archive &ar, const unsigned int version)
Definition: tria_objects.h:531
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:395
unsigned int level
Definition: grid_out.cc:4590
::TriaRawIterator<::TriaAccessor< structdim, dim, spacedim > > next_free_pair_object(const Triangulation< dim, spacedim > &tria)
std::vector< BoundaryOrMaterialId > boundary_or_material_id
Definition: tria_objects.h:179
std::vector< types::manifold_id > manifold_id
Definition: tria_objects.h:185
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:394
typename IteratorSelector::raw_hex_iterator raw_hex_iterator
Definition: tria.h:3622
const types::material_id invalid_material_id
Definition: types.h:228
Triangulation< dim, spacedim >::raw_hex_iterator next_free_hex(const Triangulation< dim, spacedim > &tria, const unsigned int level)