Reference documentation for deal.II version Git 08727cc441 2020-07-02 15:45:42 -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(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 
159  template <class Archive>
160  void
161  serialize(Archive &ar, const unsigned int version);
162  };
163 
178  std::vector<BoundaryOrMaterialId> boundary_or_material_id;
179 
184  std::vector<types::manifold_id> manifold_id;
185 
197  template <int structdim, int dim, int spacedim>
200 
212  template <int structdim, int dim, int spacedim>
215 
220  template <int dim, int spacedim>
223  const unsigned int level);
224 
228  void *&
229  user_pointer(const unsigned int i);
230 
234  const void *
235  user_pointer(const unsigned int i) const;
236 
240  unsigned int &
241  user_index(const unsigned int i);
242 
246  unsigned int
247  user_index(const unsigned int i) const;
248 
252  void
253  clear_user_data(const unsigned int i);
254 
259  void
260  clear_user_data();
261 
265  void
267 
272  std::size_t
273  memory_consumption() const;
274 
279  template <class Archive>
280  void
281  serialize(Archive &ar, const unsigned int version);
282 
291 
295  unsigned int next_free_single;
296 
300  unsigned int next_free_pair;
301 
306 
310  struct UserData
311  {
312  union
313  {
316  void *p;
319  unsigned int i;
320  };
321 
326  {
327  p = nullptr;
328  }
329 
334  template <class Archive>
335  void
336  serialize(Archive &ar, const unsigned int version);
337  };
338 
343  {
350  };
351 
352 
357  std::vector<UserData> user_data;
358 
365  };
366 
367 
368  //----------------------------------------------------------------------//
369 
370  inline unsigned int
372  {
373  // assume that each cell has the same number of faces
374 
375  unsigned int faces_per_cell = 1;
376 
377  if (this->structdim == 1)
378  faces_per_cell = GeometryInfo<1>::faces_per_cell;
379  else if (this->structdim == 2)
380  faces_per_cell = GeometryInfo<2>::faces_per_cell;
381  else if (this->structdim == 3)
382  faces_per_cell = GeometryInfo<3>::faces_per_cell;
383  else
384  AssertThrow(false, ExcNotImplemented());
385 
386  return cells.size() / faces_per_cell;
387  }
388 
389 
390 
391  inline ArrayView<int>
393  {
394  // assume that each cell has the same number of faces
395 
396  unsigned int faces_per_cell = 1;
397 
398  if (this->structdim == 1)
399  faces_per_cell = GeometryInfo<1>::faces_per_cell;
400  else if (this->structdim == 2)
401  faces_per_cell = GeometryInfo<2>::faces_per_cell;
402  else if (this->structdim == 3)
403  faces_per_cell = GeometryInfo<3>::faces_per_cell;
404  else
405  AssertThrow(false, ExcNotImplemented());
406 
407  return ArrayView<int>(cells.data() + index * faces_per_cell,
408  faces_per_cell);
409  }
410 
411 
412 
414  {
416  }
417 
418 
419 
420  inline std::size_t
422  {
423  return sizeof(BoundaryOrMaterialId);
424  }
425 
426 
427 
428  template <class Archive>
429  void
431  const unsigned int /*version*/)
432  {
433  // serialize this
434  // structure by
435  // writing and
436  // reading the larger
437  // of the two values,
438  // in order to make
439  // sure we get all
440  // bits
441  if (sizeof(material_id) > sizeof(boundary_id))
442  ar &material_id;
443  else
444  ar &boundary_id;
445  }
446 
447 
448  inline void *&
449  TriaObjects::user_pointer(const unsigned int i)
450  {
454 
455  AssertIndexRange(i, user_data.size());
456  return user_data[i].p;
457  }
458 
459 
460  inline const void *
461  TriaObjects::user_pointer(const unsigned int i) const
462  {
466 
467  AssertIndexRange(i, user_data.size());
468  return user_data[i].p;
469  }
470 
471 
472  inline unsigned int &
473  TriaObjects::user_index(const unsigned int i)
474  {
478 
479  AssertIndexRange(i, user_data.size());
480  return user_data[i].i;
481  }
482 
483 
484  inline void
485  TriaObjects::clear_user_data(const unsigned int i)
486  {
487  AssertIndexRange(i, user_data.size());
488  user_data[i].i = 0;
489  }
490 
491 
493  : structdim(-1)
498  {}
499 
500 
501  inline TriaObjects::TriaObjects(const unsigned int structdim)
502  : structdim(structdim)
507  {}
508 
509 
510  inline unsigned int
511  TriaObjects::user_index(const unsigned int i) const
512  {
516 
517  AssertIndexRange(i, user_data.size());
518  return user_data[i].i;
519  }
520 
521 
522  inline void
524  {
526  for (auto &data : user_data)
527  data.p = nullptr;
528  }
529 
530 
531  inline void
533  {
534  user_flags.assign(user_flags.size(), false);
535  }
536 
537 
538  template <class Archive>
539  void
540  TriaObjects::UserData::serialize(Archive &ar, const unsigned int)
541  {
542  // serialize this as an integer
543  ar &i;
544  }
545 
546 
547 
548  template <class Archive>
549  void
550  TriaObjects::serialize(Archive &ar, const unsigned int)
551  {
552  ar &structdim;
553  ar &cells &children;
554  ar & refinement_cases;
555  ar & used;
556  ar & user_flags;
558  ar & manifold_id;
561  }
562 
563 
564  //----------------------------------------------------------------------//
565 
566  template <int structdim_, int dim, int spacedim>
569  const Triangulation<dim, spacedim> &tria)
570  {
571  // TODO: Think of a way to ensure that we are using the correct
572  // triangulation, i.e. the one containing *this.
573 
574  AssertDimension(structdim_, this->structdim);
575 
576  int pos = next_free_single, last = used.size() - 1;
578  {
579  // first sweep forward, only use really single slots, do not use
580  // pair slots
581  for (; pos < last; ++pos)
582  if (!used[pos])
583  if (used[++pos])
584  {
585  // this was a single slot
586  pos -= 1;
587  break;
588  }
589  if (pos >= last)
590  {
592  next_free_single = used.size() - 1;
593  pos = used.size() - 1;
594  }
595  else
596  next_free_single = pos + 1;
597  }
598 
600  {
601  // second sweep, use all slots, even
602  // in pairs
603  for (; pos >= 0; --pos)
604  if (!used[pos])
605  break;
606  if (pos > 0)
607  next_free_single = pos - 1;
608  else
609  // no valid single object anymore
610  return ::TriaRawIterator<
612  }
613 
614  return ::TriaRawIterator<
616  }
617 
618 
619 
620  template <int structdim_, int dim, int spacedim>
623  {
624  // TODO: Think of a way to ensure that we are using the correct
625  // triangulation, i.e. the one containing *this.
626 
627  AssertDimension(structdim_, this->structdim);
628 
629  int pos = next_free_pair, last = used.size() - 1;
630  for (; pos < last; ++pos)
631  if (!used[pos])
632  if (!used[++pos])
633  {
634  // this was a pair slot
635  pos -= 1;
636  break;
637  }
638  if (pos >= last)
639  // no free slot
640  return ::TriaRawIterator<
642  else
643  next_free_pair = pos + 2;
644 
645  return ::TriaRawIterator<
647  }
648  } // namespace TriangulationImplementation
649 } // namespace internal
650 
651 
652 
654 
655 #endif
static const unsigned int invalid_unsigned_int
Definition: types.h:191
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1560
::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:392
unsigned int & user_index(const unsigned int i)
Definition: tria_objects.h:473
#define AssertIndexRange(index, range)
Definition: exceptions.h:1628
void serialize(Archive &ar, const unsigned int version)
Definition: tria_objects.h:540
#define AssertThrow(cond, exc)
Definition: exceptions.h:1513
#define Assert(cond, exc)
Definition: exceptions.h:1403
void serialize(Archive &ar, const unsigned int version)
Definition: tria_objects.h:550
#define DeclException0(Exception0)
Definition: exceptions.h:470
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
unsigned int level
Definition: grid_out.cc:4355
::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:178
std::vector< types::manifold_id > manifold_id
Definition: tria_objects.h:184
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
typename IteratorSelector::raw_hex_iterator raw_hex_iterator
Definition: tria.h:3490
static ::ExceptionBase & ExcNotImplemented()
const types::material_id invalid_material_id
Definition: types.h:223
Triangulation< dim, spacedim >::raw_hex_iterator next_free_hex(const Triangulation< dim, spacedim > &tria, const unsigned int level)