Reference documentation for deal.II version Git 78a8940608 2021-04-17 09:24:19 -0400
\(\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 - 2020 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 
378  unsigned int faces_per_cell = 1;
379 
380  if (this->structdim == 1)
381  faces_per_cell = GeometryInfo<1>::faces_per_cell;
382  else if (this->structdim == 2)
383  faces_per_cell = GeometryInfo<2>::faces_per_cell;
384  else if (this->structdim == 3)
385  faces_per_cell = GeometryInfo<3>::faces_per_cell;
386  else
387  AssertThrow(false, ExcNotImplemented());
388 
389  return cells.size() / faces_per_cell;
390  }
391 
392 
393 
394  inline ArrayView<int>
396  {
397  // assume that each cell has the same number of faces
398 
399  unsigned int faces_per_cell = 1;
400 
401  if (this->structdim == 1)
402  faces_per_cell = GeometryInfo<1>::faces_per_cell;
403  else if (this->structdim == 2)
404  faces_per_cell = GeometryInfo<2>::faces_per_cell;
405  else if (this->structdim == 3)
406  faces_per_cell = GeometryInfo<3>::faces_per_cell;
407  else
408  AssertThrow(false, ExcNotImplemented());
409 
410  return ArrayView<int>(cells.data() + index * faces_per_cell,
411  faces_per_cell);
412  }
413 
414 
415 
417  {
419  }
420 
421 
422 
423  inline std::size_t
425  {
426  return sizeof(BoundaryOrMaterialId);
427  }
428 
429 
430 
431  template <class Archive>
432  void
434  const unsigned int /*version*/)
435  {
436  // serialize this
437  // structure by
438  // writing and
439  // reading the larger
440  // of the two values,
441  // in order to make
442  // sure we get all
443  // bits
444  if (sizeof(material_id) > sizeof(boundary_id))
445  ar &material_id;
446  else
447  ar &boundary_id;
448  }
449 
450 
451  inline void *&
452  TriaObjects::user_pointer(const unsigned int i)
453  {
457 
458  AssertIndexRange(i, user_data.size());
459  return user_data[i].p;
460  }
461 
462 
463  inline const void *
464  TriaObjects::user_pointer(const unsigned int i) const
465  {
469 
470  AssertIndexRange(i, user_data.size());
471  return user_data[i].p;
472  }
473 
474 
475  inline unsigned int &
476  TriaObjects::user_index(const unsigned int i)
477  {
481 
482  AssertIndexRange(i, user_data.size());
483  return user_data[i].i;
484  }
485 
486 
487  inline void
488  TriaObjects::clear_user_data(const unsigned int i)
489  {
490  AssertIndexRange(i, user_data.size());
491  user_data[i].i = 0;
492  }
493 
494 
501  {}
502 
503 
504  inline TriaObjects::TriaObjects(const unsigned int structdim)
505  : structdim(structdim)
510  {}
511 
512 
513  inline unsigned int
514  TriaObjects::user_index(const unsigned int i) const
515  {
519 
520  AssertIndexRange(i, user_data.size());
521  return user_data[i].i;
522  }
523 
524 
525  inline void
527  {
529  for (auto &data : user_data)
530  data.p = nullptr;
531  }
532 
533 
534  inline void
536  {
537  user_flags.assign(user_flags.size(), false);
538  }
539 
540 
541  template <class Archive>
542  void
543  TriaObjects::UserData::serialize(Archive &ar, const unsigned int)
544  {
545  // serialize this as an integer
546  ar &i;
547  }
548 
549 
550 
551  template <class Archive>
552  void
553  TriaObjects::serialize(Archive &ar, const unsigned int)
554  {
555  ar &structdim;
556  ar &cells &children;
557  ar & refinement_cases;
558  ar & used;
559  ar & user_flags;
561  ar & manifold_id;
564  }
565 
566 
567  //----------------------------------------------------------------------//
568 
569  template <int structdim_, int dim, int spacedim>
572  const Triangulation<dim, spacedim> &tria)
573  {
574  // TODO: Think of a way to ensure that we are using the correct
575  // triangulation, i.e. the one containing *this.
576 
577  AssertDimension(structdim_, this->structdim);
578 
579  int pos = next_free_single, last = used.size() - 1;
581  {
582  // first sweep forward, only use really single slots, do not use
583  // pair slots
584  for (; pos < last; ++pos)
585  if (!used[pos])
586  if (used[++pos])
587  {
588  // this was a single slot
589  pos -= 1;
590  break;
591  }
592  if (pos >= last)
593  {
595  next_free_single = used.size() - 1;
596  pos = used.size() - 1;
597  }
598  else
599  next_free_single = pos + 1;
600  }
601 
603  {
604  // second sweep, use all slots, even
605  // in pairs
606  for (; pos >= 0; --pos)
607  if (!used[pos])
608  break;
609  if (pos > 0)
610  next_free_single = pos - 1;
611  else
612  // no valid single object anymore
613  return ::TriaRawIterator<
615  }
616 
617  return ::TriaRawIterator<
619  }
620 
621 
622 
623  template <int structdim_, int dim, int spacedim>
626  {
627  // TODO: Think of a way to ensure that we are using the correct
628  // triangulation, i.e. the one containing *this.
629 
630  AssertDimension(structdim_, this->structdim);
631 
632  int pos = next_free_pair, last = used.size() - 1;
633  for (; pos < last; ++pos)
634  if (!used[pos])
635  if (!used[++pos])
636  {
637  // this was a pair slot
638  pos -= 1;
639  break;
640  }
641  if (pos >= last)
642  // no free slot
643  return ::TriaRawIterator<
645  else
646  next_free_pair = pos + 2;
647 
648  return ::TriaRawIterator<
650  }
651  } // namespace TriangulationImplementation
652 } // namespace internal
653 
654 
655 
657 
658 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:196
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1623
::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:395
unsigned int & user_index(const unsigned int i)
Definition: tria_objects.h:476
#define AssertIndexRange(index, range)
Definition: exceptions.h:1691
void serialize(Archive &ar, const unsigned int version)
Definition: tria_objects.h:543
#define AssertThrow(cond, exc)
Definition: exceptions.h:1576
#define Assert(cond, exc)
Definition: exceptions.h:1466
void serialize(Archive &ar, const unsigned int version)
Definition: tria_objects.h:553
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:396
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:395
typename IteratorSelector::raw_hex_iterator raw_hex_iterator
Definition: tria.h:3622
static ::ExceptionBase & ExcNotImplemented()
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)