Reference documentation for deal.II version Git b84270a1d4 2020-01-17 16:21:02 -0500
\(\newcommand{\vcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\vcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
tria_objects.h
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 #ifndef dealii_tria_objects_h
17 #define dealii_tria_objects_h
18 
19 #include <deal.II/base/config.h>
20 
21 #include <deal.II/base/exceptions.h>
22 #include <deal.II/base/geometry_info.h>
23 
24 #include <deal.II/grid/tria_object.h>
25 
26 #include <vector>
27 
28 DEAL_II_NAMESPACE_OPEN
29 
30 // TODO: This should all be cleaned up. Currently, only a single
31 // function in the library makes use of the odd specializations, and
32 // this function is Triangulation::execute_refinement() in 3D. I
33 // assume, that the other refinement functions would profit from using
34 // next_free_single_object() and next_free_pair_object, but they seem
35 // to get around it.
36 
37 // TODO: The TriaObjects class contains a std::vector<G>. This is only an
38 // efficient storage scheme if G is relatively well packed, i.e. it's not a
39 // bool and then an integer and then a double, etc. Verify that this is
40 // actually the case.
41 
42 // Forward declarations
43 #ifndef DOXYGEN
44 template <int dim, int spacedim>
45 class Triangulation;
46 template <class Accessor>
47 class TriaRawIterator;
48 template <int, int, int>
49 class TriaAccessor;
50 #endif
51 
52 namespace internal
53 {
54  namespace TriangulationImplementation
55  {
68  template <typename G>
69  class TriaObjects
70  {
71  public:
75  TriaObjects();
76 
81  std::vector<G> cells;
82 
94  std::vector<int> children;
95 
101  std::vector<RefinementCase<G::dimension>> refinement_cases;
102 
111  std::vector<bool> used;
112 
122  std::vector<bool> user_flags;
123 
124 
130  {
131  union
132  {
133  types::boundary_id boundary_id;
134  types::material_id material_id;
135  };
136 
137 
142 
146  static std::size_t
148 
153  template <class Archive>
154  void
155  serialize(Archive &ar, const unsigned int version);
156  };
157 
172  std::vector<BoundaryOrMaterialId> boundary_or_material_id;
173 
178  std::vector<types::manifold_id> manifold_id;
179 
190  void
191  reserve_space(const unsigned int new_objs_in_pairs,
192  const unsigned int new_objs_single = 0);
193 
205  template <int dim, int spacedim>
207  next_free_single_object(const ::Triangulation<dim, spacedim> &tria);
208 
220  template <int dim, int spacedim>
222  next_free_pair_object(const ::Triangulation<dim, spacedim> &tria);
223 
228  template <int dim, int spacedim>
229  typename ::Triangulation<dim, spacedim>::raw_hex_iterator
230  next_free_hex(const ::Triangulation<dim, spacedim> &tria,
231  const unsigned int level);
232 
236  void
237  clear();
238 
253  bool
254  face_orientation(const unsigned int cell, const unsigned int face) const;
255 
256 
260  void *&
261  user_pointer(const unsigned int i);
262 
266  const void *
267  user_pointer(const unsigned int i) const;
268 
272  unsigned int &
273  user_index(const unsigned int i);
274 
278  unsigned int
279  user_index(const unsigned int i) const;
280 
284  void
285  clear_user_data(const unsigned int i);
286 
291  void
292  clear_user_data();
293 
297  void
299 
305  void
306  monitor_memory(const unsigned int true_dimension) const;
307 
312  std::size_t
313  memory_consumption() const;
314 
319  template <class Archive>
320  void
321  serialize(Archive &ar, const unsigned int version);
322 
328  int,
329  int,
330  << "The containers have sizes " << arg1 << " and " << arg2
331  << ", which is not as expected.");
332 
341 
342  protected:
346  unsigned int next_free_single;
347 
351  unsigned int next_free_pair;
352 
357 
361  struct UserData
362  {
363  union
364  {
367  void *p;
370  unsigned int i;
371  };
372 
377  {
378  p = nullptr;
379  }
380 
385  template <class Archive>
386  void
387  serialize(Archive &ar, const unsigned int version);
388  };
389 
394  {
401  };
402 
403 
408  std::vector<UserData> user_data;
409 
416  };
417 
424  class TriaObjectsHex : public TriaObjects<TriaObject<3>>
425  {
426  public:
434  bool
435  face_orientation(const unsigned int cell, const unsigned int face) const;
436 
437 
460  std::vector<bool> face_orientations;
461 
465  std::vector<bool> face_flips;
466 
470  std::vector<bool> face_rotations;
471 
478  void
479  reserve_space(const unsigned int new_objs);
480 
484  void
485  clear();
486 
492  void
493  monitor_memory(const unsigned int true_dimension) const;
494 
499  std::size_t
500  memory_consumption() const;
501 
506  template <class Archive>
507  void
508  serialize(Archive &ar, const unsigned int version);
509  };
510 
511 
518  class TriaObjectsQuad3D : public TriaObjects<TriaObject<2>>
519  {
520  public:
528  bool
529  face_orientation(const unsigned int cell, const unsigned int face) const;
530 
531 
536  std::vector<bool> line_orientations;
537 
545  void
546  reserve_space(const unsigned int new_quads_in_pairs,
547  const unsigned int new_quads_single = 0);
548 
552  void
553  clear();
554 
560  void
561  monitor_memory(const unsigned int true_dimension) const;
562 
567  std::size_t
568  memory_consumption() const;
569 
574  template <class Archive>
575  void
576  serialize(Archive &ar, const unsigned int version);
577  };
578 
579  //----------------------------------------------------------------------//
580 
581 
582  template <typename G>
584  {
585  material_id = numbers::invalid_material_id;
586  }
587 
588 
589 
590  template <typename G>
591  std::size_t
593  {
594  return sizeof(BoundaryOrMaterialId);
595  }
596 
597 
598 
599  template <typename G>
600  template <class Archive>
601  void
603  Archive &ar,
604  const unsigned int /*version*/)
605  {
606  // serialize this
607  // structure by
608  // writing and
609  // reading the larger
610  // of the two values,
611  // in order to make
612  // sure we get all
613  // bits
614  if (sizeof(material_id) > sizeof(boundary_id))
615  ar &material_id;
616  else
617  ar &boundary_id;
618  }
619 
620 
621  template <typename G>
622  inline bool
624  const unsigned int) const
625  {
626  return true;
627  }
628 
629 
630  template <typename G>
631  inline void *&
632  TriaObjects<G>::user_pointer(const unsigned int i)
633  {
637 
638  AssertIndexRange(i, user_data.size());
639  return user_data[i].p;
640  }
641 
642 
643  template <typename G>
644  inline const void *
645  TriaObjects<G>::user_pointer(const unsigned int i) const
646  {
650 
651  AssertIndexRange(i, user_data.size());
652  return user_data[i].p;
653  }
654 
655 
656  template <typename G>
657  inline unsigned int &
658  TriaObjects<G>::user_index(const unsigned int i)
659  {
663 
664  AssertIndexRange(i, user_data.size());
665  return user_data[i].i;
666  }
667 
668 
669  template <typename G>
670  inline void
671  TriaObjects<G>::clear_user_data(const unsigned int i)
672  {
673  AssertIndexRange(i, user_data.size());
674  user_data[i].i = 0;
675  }
676 
677 
678  template <typename G>
680  : next_free_single(numbers::invalid_unsigned_int)
681  , next_free_pair(numbers::invalid_unsigned_int)
684  {}
685 
686 
687  template <typename G>
688  inline unsigned int
689  TriaObjects<G>::user_index(const unsigned int i) const
690  {
691  Assert(user_data_type == data_unknown || user_data_type == data_index,
692  ExcPointerIndexClash());
693  user_data_type = data_index;
694 
695  AssertIndexRange(i, user_data.size());
696  return user_data[i].i;
697  }
698 
699 
700  template <typename G>
701  inline void
703  {
704  user_data_type = data_unknown;
705  for (unsigned int i = 0; i < user_data.size(); ++i)
706  user_data[i].p = nullptr;
707  }
708 
709 
710  template <typename G>
711  inline void
713  {
714  user_flags.assign(user_flags.size(), false);
715  }
716 
717 
718  template <typename G>
719  template <class Archive>
720  void
721  TriaObjects<G>::UserData::serialize(Archive &ar, const unsigned int)
722  {
723  // serialize this as an integer
724  ar &i;
725  }
726 
727 
728 
729  template <typename G>
730  template <class Archive>
731  void
732  TriaObjects<G>::serialize(Archive &ar, const unsigned int)
733  {
734  ar &cells &children;
735  ar & refinement_cases;
736  ar & used;
737  ar & user_flags;
738  ar & boundary_or_material_id;
739  ar & manifold_id;
740  ar &next_free_single &next_free_pair &reverse_order_next_free_single;
741  ar &user_data &user_data_type;
742  }
743 
744 
745  template <class Archive>
746  void
747  TriaObjectsHex::serialize(Archive &ar, const unsigned int version)
748  {
749  this->TriaObjects<TriaObject<3>>::serialize(ar, version);
750 
751  ar &face_orientations &face_flips &face_rotations;
752  }
753 
754 
755  template <class Archive>
756  void
757  TriaObjectsQuad3D::serialize(Archive &ar, const unsigned int version)
758  {
759  this->TriaObjects<TriaObject<2>>::serialize(ar, version);
760 
761  ar &line_orientations;
762  }
763 
764 
765  //----------------------------------------------------------------------//
766 
767  inline bool
768  TriaObjectsHex::face_orientation(const unsigned int cell,
769  const unsigned int face) const
770  {
771  AssertIndexRange(cell,
772  face_orientations.size() /
775 
776  return face_orientations[cell * GeometryInfo<3>::faces_per_cell + face];
777  }
778 
779  //----------------------------------------------------------------------//
780 
781  inline bool
782  TriaObjectsQuad3D::face_orientation(const unsigned int cell,
783  const unsigned int face) const
784  {
785  return line_orientations[cell * GeometryInfo<2>::faces_per_cell + face];
786  }
787 
788 
789  //----------------------------------------------------------------------//
790 
791  template <class G>
792  template <int dim, int spacedim>
795  const ::Triangulation<dim, spacedim> &tria)
796  {
797  // TODO: Think of a way to ensure that we are using the correct
798  // triangulation, i.e. the one containing *this.
799 
800  int pos = next_free_single, last = used.size() - 1;
801  if (!reverse_order_next_free_single)
802  {
803  // first sweep forward, only use really single slots, do not use
804  // pair slots
805  for (; pos < last; ++pos)
806  if (!used[pos])
807  if (used[++pos])
808  {
809  // this was a single slot
810  pos -= 1;
811  break;
812  }
813  if (pos >= last)
814  {
815  reverse_order_next_free_single = true;
816  next_free_single = used.size() - 1;
817  pos = used.size() - 1;
818  }
819  else
820  next_free_single = pos + 1;
821  }
822 
823  if (reverse_order_next_free_single)
824  {
825  // second sweep, use all slots, even
826  // in pairs
827  for (; pos >= 0; --pos)
828  if (!used[pos])
829  break;
830  if (pos > 0)
831  next_free_single = pos - 1;
832  else
833  // no valid single object anymore
834  return ::TriaRawIterator<
836  }
837 
838  return ::TriaRawIterator<
840  }
841 
842 
843 
844  template <class G>
845  template <int dim, int spacedim>
848  const ::Triangulation<dim, spacedim> &tria)
849  {
850  // TODO: Think of a way to ensure that we are using the correct
851  // triangulation, i.e. the one containing *this.
852 
853  int pos = next_free_pair, last = used.size() - 1;
854  for (; pos < last; ++pos)
855  if (!used[pos])
856  if (!used[++pos])
857  {
858  // this was a pair slot
859  pos -= 1;
860  break;
861  }
862  if (pos >= last)
863  // no free slot
864  return ::TriaRawIterator<
866  else
867  next_free_pair = pos + 2;
868 
869  return ::TriaRawIterator<
871  }
872 
873 
874 
875  // declaration of explicit specializations
876 
877  template <>
878  void
879  TriaObjects<TriaObject<2>>::monitor_memory(const unsigned int) const;
880 
881  } // namespace TriangulationImplementation
882 } // namespace internal
883 
884 
885 
886 DEAL_II_NAMESPACE_CLOSE
887 
888 #endif
::TriaRawIterator<::TriaAccessor< G::dimension, dim, spacedim > > next_free_pair_object(const ::Triangulation< dim, spacedim > &tria)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:541
static ::ExceptionBase & ExcMemoryInexact(int arg1, int arg2)
::TriaRawIterator<::TriaAccessor< G::dimension, dim, spacedim > > next_free_single_object(const ::Triangulation< dim, spacedim > &tria)
#define AssertIndexRange(index, range)
Definition: exceptions.h:1641
void monitor_memory(const unsigned int true_dimension) const
unsigned int material_id
Definition: types.h:148
bool face_orientation(const unsigned int cell, const unsigned int face) const
Definition: tria_objects.h:623
unsigned int & user_index(const unsigned int i)
Definition: tria_objects.h:658
std::vector< BoundaryOrMaterialId > boundary_or_material_id
Definition: tria_objects.h:172
#define Assert(cond, exc)
Definition: exceptions.h:1411
#define DeclException0(Exception0)
Definition: exceptions.h:473
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< G::dimension > > refinement_cases
Definition: tria_objects.h:101
std::vector< types::manifold_id > manifold_id
Definition: tria_objects.h:178
typename ::Triangulation< dim, spacedim >::raw_hex_iterator next_free_hex(const ::Triangulation< dim, spacedim > &tria, const unsigned int level)
const types::material_id invalid_material_id
Definition: types.h:217
unsigned int boundary_id
Definition: types.h:125