Reference documentation for deal.II version Git c86ef8e7a2 2021-02-25 12:12:08 +0100
\(\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.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 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_h
17 #define dealii_tria_h
18 
19 
20 #include <deal.II/base/config.h>
21 
24 #include <deal.II/base/point.h>
27 
28 #include <deal.II/grid/cell_id.h>
32 
33 #include <boost/serialization/map.hpp>
34 #include <boost/serialization/split_member.hpp>
35 #include <boost/serialization/unique_ptr.hpp>
36 #include <boost/serialization/vector.hpp>
37 #include <boost/signals2.hpp>
38 
39 #include <bitset>
40 #include <functional>
41 #include <list>
42 #include <map>
43 #include <memory>
44 #include <numeric>
45 #include <vector>
46 
47 
49 
50 #ifdef signals
51 # error \
52  "The name 'signals' is already defined. You are most likely using the QT library \
53 and using the 'signals' keyword. You can either #include the Qt headers (or any conflicting headers) \
54 *after* the deal.II headers or you can define the 'QT_NO_KEYWORDS' macro and use the 'Q_SIGNALS' macro."
55 #endif
56 
57 // Forward declarations
58 #ifndef DOXYGEN
59 template <int dim, int spacedim>
60 class Manifold;
61 
62 template <int dim>
63 struct CellData;
64 
65 struct SubCellData;
66 
68 {
69  template <int, int>
70  struct Description;
71 }
72 
73 namespace GridTools
74 {
75  template <typename CellIterator>
76  struct PeriodicFacePair;
77 }
78 
79 template <int, int, int>
80 class TriaAccessor;
81 template <int spacedim>
82 class TriaAccessor<0, 1, spacedim>;
83 template <int, int, int>
84 class TriaAccessorBase;
85 
86 namespace internal
87 {
88  namespace TriangulationImplementation
89  {
90  class TriaFaces;
91 
92  class TriaObjects;
93 
94  template <int, int>
95  class Policy;
96 
102  struct Implementation;
103  struct ImplementationMixedMesh;
104  } // namespace TriangulationImplementation
105 
106  namespace TriaAccessorImplementation
107  {
108  struct Implementation;
109  }
110 } // namespace internal
111 
112 namespace hp
113 {
114  template <int dim, int spacedim>
115  class DoFHandler;
116 }
117 #endif
118 
119 
120 /*------------------------------------------------------------------------*/
121 
122 
123 namespace internal
124 {
129  namespace TriangulationImplementation
130  {
142  template <int dim>
143  struct NumberCache
144  {};
145 
157  template <>
158  struct NumberCache<1>
159  {
163  unsigned int n_levels;
164 
168  unsigned int n_lines;
169 
173  std::vector<unsigned int> n_lines_level;
174 
178  unsigned int n_active_lines;
179 
183  std::vector<unsigned int> n_active_lines_level;
184 
188  NumberCache();
189 
194  std::size_t
195  memory_consumption() const;
196 
202  template <class Archive>
203  void
204  serialize(Archive &ar, const unsigned int version);
205  };
206 
207 
220  template <>
221  struct NumberCache<2> : public NumberCache<1>
222  {
226  unsigned int n_quads;
227 
231  std::vector<unsigned int> n_quads_level;
232 
236  unsigned int n_active_quads;
237 
241  std::vector<unsigned int> n_active_quads_level;
242 
246  NumberCache();
247 
252  std::size_t
253  memory_consumption() const;
254 
260  template <class Archive>
261  void
262  serialize(Archive &ar, const unsigned int version);
263  };
264 
265 
279  template <>
280  struct NumberCache<3> : public NumberCache<2>
281  {
285  unsigned int n_hexes;
286 
290  std::vector<unsigned int> n_hexes_level;
291 
295  unsigned int n_active_hexes;
296 
300  std::vector<unsigned int> n_active_hexes_level;
301 
305  NumberCache();
306 
311  std::size_t
312  memory_consumption() const;
313 
319  template <class Archive>
320  void
321  serialize(Archive &ar, const unsigned int version);
322  };
323  } // namespace TriangulationImplementation
324 } // namespace internal
325 
326 
327 /*------------------------------------------------------------------------*/
328 
329 
1126 template <int dim, int spacedim = dim>
1128 {
1129 private:
1134  using IteratorSelector =
1136 
1137 public:
1143  {
1148  none = 0x0,
1192  limit_level_difference_at_vertices = 0x1,
1213  eliminate_unrefined_islands = 0x2,
1229  patch_level_1 = 0x4,
1253  coarsest_level_1 = 0x8,
1278  allow_anisotropic_smoothing = 0x10,
1311  eliminate_refined_inner_islands = 0x100,
1316  eliminate_refined_boundary_islands = 0x200,
1322  do_not_produce_unrefined_islands = 0x400,
1323 
1328  smoothing_on_refinement =
1329  (limit_level_difference_at_vertices | eliminate_unrefined_islands),
1334  smoothing_on_coarsening =
1335  (eliminate_refined_inner_islands | eliminate_refined_boundary_islands |
1336  do_not_produce_unrefined_islands),
1337 
1343  maximum_smoothing = 0xffff ^ allow_anisotropic_smoothing
1344  };
1345 
1362 
1368 
1386 
1400  using face_iterator = TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>;
1401 
1413  using active_face_iterator =
1414  TriaActiveIterator<TriaAccessor<dim - 1, dim, spacedim>>;
1415 
1425 
1438  using active_vertex_iterator =
1440 
1448  using line_iterator = typename IteratorSelector::line_iterator;
1449 
1463  using active_line_iterator = typename IteratorSelector::active_line_iterator;
1464 
1472  using quad_iterator = typename IteratorSelector::quad_iterator;
1473 
1487  using active_quad_iterator = typename IteratorSelector::active_quad_iterator;
1488 
1496  using hex_iterator = typename IteratorSelector::hex_iterator;
1497 
1507  using active_hex_iterator = typename IteratorSelector::active_hex_iterator;
1508 
1528  {
1535  virtual ~DistortedCellList() noexcept override = default;
1536 
1541  std::list<typename Triangulation<dim, spacedim>::cell_iterator>
1543  };
1544 
1548  static const unsigned int dimension = dim;
1549 
1553  static const unsigned int space_dimension = spacedim;
1554 
1569  Triangulation(const MeshSmoothing smooth_grid = none,
1570  const bool check_for_distorted_cells = false);
1571 
1587  Triangulation(const Triangulation<dim, spacedim> &) = delete;
1588 
1596 
1600  Triangulation &
1601  operator=(Triangulation<dim, spacedim> &&tria) noexcept;
1602 
1606  virtual ~Triangulation() override;
1607 
1614  virtual void
1615  clear();
1616 
1622  virtual void
1623  set_mesh_smoothing(const MeshSmoothing mesh_smoothing);
1624 
1628  virtual const MeshSmoothing &
1629  get_mesh_smoothing() const;
1630 
1653  void
1654  set_manifold(const types::manifold_id number,
1655  const Manifold<dim, spacedim> &manifold_object);
1656 
1669  void
1670  reset_manifold(const types::manifold_id manifold_number);
1671 
1683  void
1684  reset_all_manifolds();
1685 
1694  void
1695  set_all_manifold_ids(const types::manifold_id number);
1696 
1705  void
1706  set_all_manifold_ids_on_boundary(const types::manifold_id number);
1707 
1717  void
1718  set_all_manifold_ids_on_boundary(const types::boundary_id b_id,
1719  const types::manifold_id number);
1720 
1732  const Manifold<dim, spacedim> &
1733  get_manifold(const types::manifold_id number) const;
1734 
1747  virtual std::vector<types::boundary_id>
1748  get_boundary_ids() const;
1749 
1762  virtual std::vector<types::manifold_id>
1763  get_manifold_ids() const;
1764 
1791  virtual void
1792  copy_triangulation(const Triangulation<dim, spacedim> &other_tria);
1793 
1841  virtual void
1842  create_triangulation(const std::vector<Point<spacedim>> &vertices,
1843  const std::vector<CellData<dim>> & cells,
1844  const SubCellData & subcelldata);
1845 
1858  virtual void
1861  &construction_data);
1862 
1871  virtual void
1872  create_triangulation_compatibility(
1873  const std::vector<Point<spacedim>> &vertices,
1874  const std::vector<CellData<dim>> & cells,
1875  const SubCellData & subcelldata);
1876 
1883  void
1884  flip_all_direction_flags();
1885 
1897  void
1898  set_all_refine_flags();
1899 
1922  void
1923  refine_global(const unsigned int times = 1);
1924 
1937  void
1938  coarsen_global(const unsigned int times = 1);
1939 
1971  virtual void
1972  execute_coarsening_and_refinement();
1973 
2003  virtual bool
2004  prepare_coarsening_and_refinement();
2005 
2006  /*
2007  * @}
2008  */
2009 
2024  {
2041  CELL_INVALID
2042  };
2043 
2049  template <typename T>
2051  {
2052  using result_type = T;
2053 
2054  template <typename InputIterator>
2055  T
2056  operator()(InputIterator first, InputIterator last) const
2057  {
2058  return std::accumulate(first, last, T());
2059  }
2060  };
2061 
2071  struct Signals
2072  {
2080  boost::signals2::signal<void()> create;
2081 
2089  boost::signals2::signal<void()> pre_refinement;
2090 
2096  boost::signals2::signal<void()> post_refinement;
2097 
2104  boost::signals2::signal<void()> pre_partition;
2105 
2113  boost::signals2::signal<void()> mesh_movement;
2114 
2122  boost::signals2::signal<void(
2123  const typename Triangulation<dim, spacedim>::cell_iterator &cell)>
2125 
2132  boost::signals2::signal<void(
2133  const typename Triangulation<dim, spacedim>::cell_iterator &cell)>
2135 
2142  boost::signals2::signal<void(
2143  const Triangulation<dim, spacedim> &destination_tria)>
2145 
2159  boost::signals2::signal<void()> clear;
2160 
2170  boost::signals2::signal<void()> any_change;
2171 
2197  boost::signals2::signal<unsigned int(const cell_iterator &,
2198  const CellStatus),
2201 
2213  boost::signals2::signal<void()> pre_distributed_refinement;
2214 
2223  boost::signals2::signal<void()> post_distributed_refinement;
2224 
2235  boost::signals2::signal<void()> pre_distributed_repartition;
2236 
2242  boost::signals2::signal<void()> post_distributed_repartition;
2243 
2250  boost::signals2::signal<void()> pre_distributed_save;
2251 
2257  boost::signals2::signal<void()> post_distributed_save;
2258 
2265  boost::signals2::signal<void()> pre_distributed_load;
2266 
2272  boost::signals2::signal<void()> post_distributed_load;
2273  };
2274 
2278  mutable Signals signals;
2279 
2280  /*
2281  * @}
2282  */
2283 
2293  void
2294  save_refine_flags(std::ostream &out) const;
2295 
2299  void
2300  save_refine_flags(std::vector<bool> &v) const;
2301 
2305  void
2306  load_refine_flags(std::istream &in);
2307 
2311  void
2312  load_refine_flags(const std::vector<bool> &v);
2313 
2317  void
2318  save_coarsen_flags(std::ostream &out) const;
2319 
2323  void
2324  save_coarsen_flags(std::vector<bool> &v) const;
2325 
2329  void
2330  load_coarsen_flags(std::istream &out);
2331 
2335  void
2336  load_coarsen_flags(const std::vector<bool> &v);
2337 
2342  bool
2343  get_anisotropic_refinement_flag() const;
2344 
2345  /*
2346  * @}
2347  */
2348 
2358  void
2359  clear_user_flags();
2360 
2366  void
2367  save_user_flags(std::ostream &out) const;
2368 
2374  void
2375  save_user_flags(std::vector<bool> &v) const;
2376 
2381  void
2382  load_user_flags(std::istream &in);
2383 
2388  void
2389  load_user_flags(const std::vector<bool> &v);
2390 
2395  void
2396  clear_user_flags_line();
2397 
2402  void
2403  save_user_flags_line(std::ostream &out) const;
2404 
2410  void
2411  save_user_flags_line(std::vector<bool> &v) const;
2412 
2417  void
2418  load_user_flags_line(std::istream &in);
2419 
2424  void
2425  load_user_flags_line(const std::vector<bool> &v);
2426 
2431  void
2432  clear_user_flags_quad();
2433 
2438  void
2439  save_user_flags_quad(std::ostream &out) const;
2440 
2446  void
2447  save_user_flags_quad(std::vector<bool> &v) const;
2448 
2453  void
2454  load_user_flags_quad(std::istream &in);
2455 
2460  void
2461  load_user_flags_quad(const std::vector<bool> &v);
2462 
2463 
2468  void
2469  clear_user_flags_hex();
2470 
2475  void
2476  save_user_flags_hex(std::ostream &out) const;
2477 
2483  void
2484  save_user_flags_hex(std::vector<bool> &v) const;
2485 
2490  void
2491  load_user_flags_hex(std::istream &in);
2492 
2497  void
2498  load_user_flags_hex(const std::vector<bool> &v);
2499 
2505  void
2506  clear_user_data();
2507 
2513  void
2514  save_user_indices(std::vector<unsigned int> &v) const;
2515 
2520  void
2521  load_user_indices(const std::vector<unsigned int> &v);
2522 
2528  void
2529  save_user_pointers(std::vector<void *> &v) const;
2530 
2535  void
2536  load_user_pointers(const std::vector<void *> &v);
2537 
2543  void
2544  save_user_indices_line(std::vector<unsigned int> &v) const;
2545 
2550  void
2551  load_user_indices_line(const std::vector<unsigned int> &v);
2552 
2558  void
2559  save_user_indices_quad(std::vector<unsigned int> &v) const;
2560 
2565  void
2566  load_user_indices_quad(const std::vector<unsigned int> &v);
2567 
2573  void
2574  save_user_indices_hex(std::vector<unsigned int> &v) const;
2575 
2580  void
2581  load_user_indices_hex(const std::vector<unsigned int> &v);
2587  void
2588  save_user_pointers_line(std::vector<void *> &v) const;
2589 
2594  void
2595  load_user_pointers_line(const std::vector<void *> &v);
2596 
2602  void
2603  save_user_pointers_quad(std::vector<void *> &v) const;
2604 
2609  void
2610  load_user_pointers_quad(const std::vector<void *> &v);
2611 
2617  void
2618  save_user_pointers_hex(std::vector<void *> &v) const;
2619 
2624  void
2625  load_user_pointers_hex(const std::vector<void *> &v);
2626 
2627  /*
2628  * @}
2629  */
2630 
2654  begin(const unsigned int level = 0) const;
2655 
2687  begin_active(const unsigned int level = 0) const;
2688 
2694  end() const;
2695 
2715  end(const unsigned int level) const;
2716 
2737  end_active(const unsigned int level) const;
2738 
2739 
2744  last() const;
2745 
2750  last_active() const;
2751 
2765  create_cell_iterator(const CellId &cell_id) const;
2766 
2782  cell_iterators() const;
2783 
2820  active_cell_iterators() const;
2821 
2838  cell_iterators_on_level(const unsigned int level) const;
2839 
2856  active_cell_iterators_on_level(const unsigned int level) const;
2857 
2858  /*
2859  * @}
2860  */
2861 
2862  /*-------------------------------------------------------------------------*/
2863 
2873  begin_face() const;
2874 
2879  begin_active_face() const;
2880 
2886  end_face() const;
2887 
2907  active_face_iterators() const;
2908 
2909  /*
2910  * @}
2911  */
2912 
2913  /*-------------------------------------------------------------------------*/
2914 
2925  begin_vertex() const;
2926 
2933  begin_active_vertex() const;
2934 
2941  end_vertex() const;
2942 
2943  /*
2944  * @}
2945  */
2946 
2964  unsigned int
2965  n_lines() const;
2966 
2970  unsigned int
2971  n_lines(const unsigned int level) const;
2972 
2976  unsigned int
2977  n_active_lines() const;
2978 
2982  unsigned int
2983  n_active_lines(const unsigned int level) const;
2984 
2988  unsigned int
2989  n_quads() const;
2990 
2994  unsigned int
2995  n_quads(const unsigned int level) const;
2996 
3000  unsigned int
3001  n_active_quads() const;
3002 
3006  unsigned int
3007  n_active_quads(const unsigned int level) const;
3008 
3012  unsigned int
3013  n_hexs() const;
3014 
3019  unsigned int
3020  n_hexs(const unsigned int level) const;
3021 
3025  unsigned int
3026  n_active_hexs() const;
3027 
3032  unsigned int
3033  n_active_hexs(const unsigned int level) const;
3034 
3039  unsigned int
3040  n_cells() const;
3041 
3046  unsigned int
3047  n_cells(const unsigned int level) const;
3048 
3053  unsigned int
3054  n_active_cells() const;
3055 
3063  virtual types::global_cell_index
3064  n_global_active_cells() const;
3065 
3066 
3071  unsigned int
3072  n_active_cells(const unsigned int level) const;
3073 
3079  unsigned int
3080  n_faces() const;
3081 
3087  unsigned int
3088  n_active_faces() const;
3089 
3107  unsigned int
3108  n_levels() const;
3109 
3116  virtual unsigned int
3117  n_global_levels() const;
3118 
3128  virtual bool
3129  has_hanging_nodes() const;
3130 
3138  unsigned int
3139  n_vertices() const;
3140 
3149  const std::vector<Point<spacedim>> &
3150  get_vertices() const;
3151 
3156  unsigned int
3157  n_used_vertices() const;
3158 
3162  bool
3163  vertex_used(const unsigned int index) const;
3164 
3169  const std::vector<bool> &
3170  get_used_vertices() const;
3171 
3183  unsigned int
3184  max_adjacent_cells() const;
3185 
3192  virtual types::subdomain_id
3193  locally_owned_subdomain() const;
3194 
3205  get_triangulation();
3206 
3212  get_triangulation() const;
3213 
3214 
3215  /*
3216  * @}
3217  */
3218 
3233  unsigned int
3234  n_raw_lines() const;
3235 
3245  unsigned int
3246  n_raw_lines(const unsigned int level) const;
3247 
3257  unsigned int
3258  n_raw_quads() const;
3259 
3269  unsigned int
3270  n_raw_quads(const unsigned int level) const;
3271 
3281  unsigned int
3282  n_raw_hexs(const unsigned int level) const;
3283 
3293  unsigned int
3294  n_raw_cells(const unsigned int level) const;
3295 
3307  unsigned int
3308  n_raw_faces() const;
3309 
3310  /*
3311  * @}
3312  */
3313 
3322  virtual std::size_t
3323  memory_consumption() const;
3324 
3334  template <class Archive>
3335  void
3336  save(Archive &ar, const unsigned int version) const;
3337 
3355  template <class Archive>
3356  void
3357  load(Archive &ar, const unsigned int version);
3358 
3359 
3375  virtual void
3376  add_periodicity(
3377  const std::vector<GridTools::PeriodicFacePair<cell_iterator>> &);
3378 
3382  const std::map<
3383  std::pair<cell_iterator, unsigned int>,
3384  std::pair<std::pair<cell_iterator, unsigned int>, std::bitset<3>>> &
3385  get_periodic_face_map() const;
3386 
3391  const std::vector<ReferenceCell> &
3392  get_reference_cells() const;
3393 
3398  bool
3399  all_reference_cells_are_hyper_cube() const;
3400 
3401 #ifdef DOXYGEN
3402 
3407  template <class Archive>
3408  void
3409  serialize(Archive &archive, const unsigned int version);
3410 #else
3411  // This macro defines the serialize() method that is compatible with
3412  // the templated save() and load() method that have been implemented.
3413  BOOST_SERIALIZATION_SPLIT_MEMBER()
3414 #endif
3415 
3426  DeclException2(ExcInvalidLevel,
3427  int,
3428  int,
3429  << "You are requesting information from refinement level "
3430  << arg1
3431  << " of a triangulation, but this triangulation only has "
3432  << arg2 << " refinement levels. The given level " << arg1
3433  << " must be *less* than " << arg2 << ".");
3441  ExcTriangulationNotEmpty,
3442  int,
3443  int,
3444  << "You are trying to perform an operation on a triangulation "
3445  << "that is only allowed if the triangulation is currently empty. "
3446  << "However, it currently stores " << arg1 << " vertices and has "
3447  << "cells on " << arg2 << " levels.");
3453  DeclException0(ExcGridReadError);
3464  DeclException1(ExcEmptyLevel,
3465  int,
3466  << "You tried to do something on level " << arg1
3467  << ", but this level is empty.");
3473  DeclException0(ExcNonOrientableTriangulation);
3474 
3482  DeclException1(ExcBoundaryIdNotFound,
3484  << "The given boundary_id " << arg1
3485  << " is not defined in this Triangulation!");
3486 
3493  ExcInconsistentCoarseningFlags,
3494  "A cell is flagged for coarsening, but either not all of its siblings "
3495  "are active or flagged for coarsening as well. Please clean up all "
3496  "coarsen flags on your triangulation via "
3497  "Triangulation::prepare_coarsening_and_refinement() beforehand!");
3498 
3499  /*
3500  * @}
3501  */
3502 
3503 protected:
3509 
3514  std::vector<ReferenceCell> reference_cells;
3515 
3529  static void
3530  write_bool_vector(const unsigned int magic_number1,
3531  const std::vector<bool> &v,
3532  const unsigned int magic_number2,
3533  std::ostream & out);
3534 
3539  static void
3540  read_bool_vector(const unsigned int magic_number1,
3541  std::vector<bool> &v,
3542  const unsigned int magic_number2,
3543  std::istream & in);
3544 
3549  void
3550  update_periodic_face_map();
3551 
3555  virtual void
3556  update_reference_cells();
3557 
3558 
3559 private:
3564  std::unique_ptr<
3567 
3574  std::vector<GridTools::PeriodicFacePair<cell_iterator>>
3576 
3581  std::map<std::pair<cell_iterator, unsigned int>,
3582  std::pair<std::pair<cell_iterator, unsigned int>, std::bitset<3>>>
3584 
3599  using raw_face_iterator =
3600  TriaRawIterator<TriaAccessor<dim - 1, dim, spacedim>>;
3601  using raw_vertex_iterator =
3603  using raw_line_iterator = typename IteratorSelector::raw_line_iterator;
3604  using raw_quad_iterator = typename IteratorSelector::raw_quad_iterator;
3605  using raw_hex_iterator = typename IteratorSelector::raw_hex_iterator;
3606 
3612  begin_raw(const unsigned int level = 0) const;
3613 
3619  end_raw(const unsigned int level) const;
3620 
3621  /*
3622  * @}
3623  */
3624 
3637  begin_raw_line(const unsigned int level = 0) const;
3638 
3657  begin_line(const unsigned int level = 0) const;
3658 
3677  begin_active_line(const unsigned int level = 0) const;
3678 
3684  end_line() const;
3685 
3686  /*
3687  * @}
3688  */
3689 
3715  begin_raw_quad(const unsigned int level = 0) const;
3716 
3735  begin_quad(const unsigned int level = 0) const;
3736 
3755  begin_active_quad(const unsigned int level = 0) const;
3756 
3762  end_quad() const;
3763 
3764  /*
3765  * @}
3766  */
3767 
3792  begin_raw_hex(const unsigned int level = 0) const;
3793 
3811  hex_iterator
3812  begin_hex(const unsigned int level = 0) const;
3813 
3832  begin_active_hex(const unsigned int level = 0) const;
3833 
3838  hex_iterator
3839  end_hex() const;
3840 
3841  /*
3842  * @}
3843  */
3844 
3845 
3859  void
3860  clear_despite_subscriptions();
3861 
3868  void
3869  reset_active_cell_indices();
3870 
3874  void
3875  reset_global_cell_indices();
3876 
3890  DistortedCellList
3891  execute_refinement();
3892 
3899  void
3900  execute_coarsening();
3901 
3906  void
3907  fix_coarsen_flags();
3908 
3922  virtual unsigned int
3923  coarse_cell_id_to_coarse_cell_index(
3925 
3926 
3939  virtual types::coarse_cell_id
3940  coarse_cell_index_to_coarse_cell_id(
3941  const unsigned int coarse_cell_index) const;
3942 
3947  std::vector<
3948  std::unique_ptr<::internal::TriangulationImplementation::TriaLevel>>
3950 
3956  std::unique_ptr<::internal::TriangulationImplementation::TriaFaces>
3958 
3959 
3963  std::vector<Point<spacedim>> vertices;
3964 
3968  std::vector<bool> vertices_used;
3969 
3974  std::map<types::manifold_id, std::unique_ptr<const Manifold<dim, spacedim>>>
3976 
3981 
3982 
3988 
3999 
4014  std::unique_ptr<std::map<unsigned int, types::boundary_id>>
4016 
4017 
4037  std::unique_ptr<std::map<unsigned int, types::manifold_id>>
4039 
4040  // make a couple of classes friends
4041  template <int, int, int>
4042  friend class TriaAccessorBase;
4043  template <int, int, int>
4044  friend class TriaAccessor;
4045  friend class TriaAccessor<0, 1, spacedim>;
4046 
4047  friend class CellAccessor<dim, spacedim>;
4048 
4049  friend struct ::internal::TriaAccessorImplementation::Implementation;
4050 
4051  friend struct ::internal::TriangulationImplementation::Implementation;
4052  friend struct ::internal::TriangulationImplementation::
4053  ImplementationMixedMesh;
4054 
4055  friend class ::internal::TriangulationImplementation::TriaObjects;
4056 
4057  // explicitly check for sensible template arguments, but not on windows
4058  // because MSVC creates bogus warnings during normal compilation
4059 #ifndef DEAL_II_MSVC
4060  static_assert(dim <= spacedim,
4061  "The dimension <dim> of a Triangulation must be less than or "
4062  "equal to the space dimension <spacedim> in which it lives.");
4063 #endif
4064 };
4065 
4066 
4067 #ifndef DOXYGEN
4068 
4069 
4070 
4071 namespace internal
4072 {
4073  namespace TriangulationImplementation
4074  {
4075  template <class Archive>
4076  void
4077  NumberCache<1>::serialize(Archive &ar, const unsigned int)
4078  {
4079  ar &n_levels;
4080  ar &n_lines &n_lines_level;
4081  ar &n_active_lines &n_active_lines_level;
4082  }
4083 
4084 
4085  template <class Archive>
4086  void
4087  NumberCache<2>::serialize(Archive &ar, const unsigned int version)
4088  {
4089  this->NumberCache<1>::serialize(ar, version);
4090 
4091  ar &n_quads &n_quads_level;
4092  ar &n_active_quads &n_active_quads_level;
4093  }
4094 
4095 
4096  template <class Archive>
4097  void
4098  NumberCache<3>::serialize(Archive &ar, const unsigned int version)
4099  {
4100  this->NumberCache<2>::serialize(ar, version);
4101 
4102  ar &n_hexes &n_hexes_level;
4103  ar &n_active_hexes &n_active_hexes_level;
4104  }
4105 
4106  } // namespace TriangulationImplementation
4107 } // namespace internal
4108 
4109 
4110 template <int dim, int spacedim>
4111 inline bool
4112 Triangulation<dim, spacedim>::vertex_used(const unsigned int index) const
4113 {
4114  AssertIndexRange(index, vertices_used.size());
4115  return vertices_used[index];
4116 }
4117 
4118 
4119 
4120 template <int dim, int spacedim>
4121 inline unsigned int
4123 {
4124  return number_cache.n_levels;
4125 }
4126 
4127 template <int dim, int spacedim>
4128 inline unsigned int
4130 {
4131  return number_cache.n_levels;
4132 }
4133 
4134 
4135 template <int dim, int spacedim>
4136 inline unsigned int
4138 {
4139  return vertices.size();
4140 }
4141 
4142 
4143 
4144 template <int dim, int spacedim>
4145 inline const std::vector<Point<spacedim>> &
4147 {
4148  return vertices;
4149 }
4150 
4151 
4152 template <int dim, int spacedim>
4153 template <class Archive>
4154 void
4155 Triangulation<dim, spacedim>::save(Archive &ar, const unsigned int) const
4156 {
4157  // as discussed in the documentation, do not store the signals as
4158  // well as boundary and manifold description but everything else
4159  ar &smooth_grid;
4160 
4161  unsigned int n_levels = levels.size();
4162  ar & n_levels;
4163  for (const auto &level : levels)
4164  ar &level;
4165 
4166  // boost dereferences a nullptr when serializing a nullptr
4167  // at least up to 1.65.1. This causes problems with clang-5.
4168  // Therefore, work around it.
4169  bool faces_is_nullptr = (faces.get() == nullptr);
4170  ar & faces_is_nullptr;
4171  if (!faces_is_nullptr)
4172  ar &faces;
4173 
4174  ar &vertices;
4175  ar &vertices_used;
4176 
4177  ar &anisotropic_refinement;
4178  ar &number_cache;
4179 
4180  ar &check_for_distorted_cells;
4181 
4182  if (dim == 1)
4183  {
4184  ar &vertex_to_boundary_id_map_1d;
4185  ar &vertex_to_manifold_id_map_1d;
4186  }
4187 }
4188 
4189 
4190 
4191 template <int dim, int spacedim>
4192 template <class Archive>
4193 void
4194 Triangulation<dim, spacedim>::load(Archive &ar, const unsigned int)
4195 {
4196  // clear previous content. this also calls the respective signal
4197  clear();
4198 
4199  // as discussed in the documentation, do not store the signals as
4200  // well as boundary and manifold description but everything else
4201  ar &smooth_grid;
4202 
4203  unsigned int size;
4204  ar & size;
4205  levels.resize(size);
4206  for (auto &level_ : levels)
4207  {
4208  std::unique_ptr<internal::TriangulationImplementation::TriaLevel> level;
4209  ar & level;
4210  level_ = std::move(level);
4211  }
4212 
4213  // Workaround for nullptr, see in save().
4214  bool faces_is_nullptr = true;
4215  ar & faces_is_nullptr;
4216  if (!faces_is_nullptr)
4217  ar &faces;
4218 
4219  ar &vertices;
4220  ar &vertices_used;
4221 
4222  ar &anisotropic_refinement;
4223  ar &number_cache;
4224 
4225  // the levels do not serialize the active_cell_indices because
4226  // they are easy enough to rebuild upon re-loading data. do
4227  // this here. don't forget to first resize the fields appropriately
4228  {
4229  for (auto &level : levels)
4230  level->active_cell_indices.resize(level->refine_flags.size());
4231  reset_active_cell_indices();
4232  }
4233 
4234 
4235  bool my_check_for_distorted_cells;
4236  ar & my_check_for_distorted_cells;
4237 
4238  Assert(my_check_for_distorted_cells == check_for_distorted_cells,
4239  ExcMessage("The triangulation loaded into here must have the "
4240  "same setting with regard to reporting distorted "
4241  "cell as the one previously stored."));
4242 
4243  if (dim == 1)
4244  {
4245  ar &vertex_to_boundary_id_map_1d;
4246  ar &vertex_to_manifold_id_map_1d;
4247  }
4248 
4249  // trigger the create signal to indicate
4250  // that new content has been imported into
4251  // the triangulation
4252  signals.create();
4253 }
4254 
4255 
4256 
4257 template <int dim, int spacedim>
4258 inline unsigned int
4261 {
4262  return coarse_cell_id;
4263 }
4264 
4265 
4266 
4267 template <int dim, int spacedim>
4268 inline types::coarse_cell_id
4270  const unsigned int coarse_cell_index) const
4271 {
4272  return coarse_cell_index;
4273 }
4274 
4275 
4276 
4277 /* -------------- declaration of explicit specializations ------------- */
4278 
4279 template <>
4280 unsigned int
4282 template <>
4283 unsigned int
4284 Triangulation<1, 1>::n_quads(const unsigned int level) const;
4285 template <>
4286 unsigned int
4287 Triangulation<1, 1>::n_raw_quads(const unsigned int level) const;
4288 template <>
4289 unsigned int
4290 Triangulation<2, 2>::n_raw_quads(const unsigned int level) const;
4291 template <>
4292 unsigned int
4293 Triangulation<1, 1>::n_raw_hexs(const unsigned int level) const;
4294 template <>
4295 unsigned int
4296 Triangulation<1, 1>::n_active_quads(const unsigned int level) const;
4297 template <>
4298 unsigned int
4300 template <>
4301 unsigned int
4303 
4304 
4305 // -------------------------------------------------------------------
4306 // -- Explicit specializations for codimension one grids
4307 
4308 
4309 template <>
4310 unsigned int
4312 template <>
4313 unsigned int
4314 Triangulation<1, 2>::n_quads(const unsigned int level) const;
4315 template <>
4316 unsigned int
4317 Triangulation<1, 2>::n_raw_quads(const unsigned int level) const;
4318 template <>
4319 unsigned int
4320 Triangulation<2, 3>::n_raw_quads(const unsigned int level) const;
4321 template <>
4322 unsigned int
4323 Triangulation<1, 2>::n_raw_hexs(const unsigned int level) const;
4324 template <>
4325 unsigned int
4326 Triangulation<1, 2>::n_active_quads(const unsigned int level) const;
4327 template <>
4328 unsigned int
4330 template <>
4331 unsigned int
4333 
4334 // -------------------------------------------------------------------
4335 // -- Explicit specializations for codimension two grids
4336 
4337 
4338 template <>
4339 unsigned int
4341 template <>
4342 unsigned int
4343 Triangulation<1, 3>::n_quads(const unsigned int level) const;
4344 template <>
4345 unsigned int
4346 Triangulation<1, 3>::n_raw_quads(const unsigned int level) const;
4347 template <>
4348 unsigned int
4349 Triangulation<2, 3>::n_raw_quads(const unsigned int level) const;
4350 template <>
4351 unsigned int
4352 Triangulation<1, 3>::n_raw_hexs(const unsigned int level) const;
4353 template <>
4354 unsigned int
4355 Triangulation<1, 3>::n_active_quads(const unsigned int level) const;
4356 template <>
4357 unsigned int
4359 template <>
4360 unsigned int
4362 
4363 
4364 #endif // DOXYGEN
4365 
4367 
4368 // Include tria_accessor.h here, so that it is possible for an end
4369 // user to use the iterators of Triangulation<dim> directly without
4370 // the need to include tria_accessor.h separately. (Otherwise the
4371 // iterators are an 'opaque' or 'incomplete' type.)
4373 
4374 #endif
boost::signals2::signal< void(const Triangulation< dim, spacedim > &destination_tria)> copy
Definition: tria.h:2144
::internal::TriangulationImplementation::NumberCache< dim > number_cache
Definition: tria.h:3998
virtual unsigned int coarse_cell_id_to_coarse_cell_index(const types::coarse_cell_id coarse_cell_id) const
boost::signals2::signal< void()> post_refinement
Definition: tria.h:2096
unsigned int n_vertices() const
boost::signals2::signal< void(const typename Triangulation< dim, spacedim >::cell_iterator &cell)> pre_coarsening_on_cell
Definition: tria.h:2124
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:538
boost::signals2::signal< void()> post_distributed_repartition
Definition: tria.h:2242
typename IteratorSelector::line_iterator line_iterator
Definition: tria.h:1448
T operator()(InputIterator first, InputIterator last) const
Definition: tria.h:2056
MeshSmoothing smooth_grid
Definition: tria.h:3508
typename IteratorSelector::raw_line_iterator raw_line_iterator
Definition: tria.h:3603
std::vector< Point< spacedim > > vertices
Definition: tria.h:3963
std::map< types::manifold_id, std::unique_ptr< const Manifold< dim, spacedim > > > manifold
Definition: tria.h:3975
typename IteratorSelector::hex_iterator hex_iterator
Definition: tria.h:1496
typename IteratorSelector::active_line_iterator active_line_iterator
Definition: tria.h:1463
boost::signals2::signal< void()> pre_distributed_refinement
Definition: tria.h:2213
bool anisotropic_refinement
Definition: tria.h:3980
#define AssertIndexRange(index, range)
Definition: exceptions.h:1691
std::unique_ptr< ::internal::TriangulationImplementation::Policy< dim, spacedim > > policy
Definition: tria.h:3566
std::vector< GridTools::PeriodicFacePair< cell_iterator > > periodic_face_pairs_level_0
Definition: tria.h:3575
std::vector< ReferenceCell > reference_cells
Definition: tria.h:3514
static ::ExceptionBase & ExcFacesHaveNoLevel()
unsigned int n_levels() const
void create_triangulation(Triangulation< dim, dim > &tria, const AdditionalData &additional_data=AdditionalData())
static ::ExceptionBase & ExcMessage(std::string arg1)
boost::signals2::signal< void()> pre_refinement
Definition: tria.h:2089
const bool check_for_distorted_cells
Definition: tria.h:3987
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
static const char T
unsigned int n_raw_quads() const
Definition: tria.cc:13102
#define Assert(cond, exc)
Definition: exceptions.h:1466
Signals signals
Definition: tria.h:2278
std::vector< std::unique_ptr<::internal::TriangulationImplementation::TriaLevel > > levels
Definition: tria.h:3949
boost::signals2::signal< void()> clear
Definition: tria.h:2159
unsigned int n_raw_hexs(const unsigned int level) const
Definition: tria.cc:13158
unsigned int max_adjacent_cells() const
Definition: tria.cc:13275
void load(Archive &ar, const unsigned int version)
std::list< typename Triangulation< dim, spacedim >::cell_iterator > distorted_cells
Definition: tria.h:1542
unsigned int n_quads() const
Definition: tria.cc:13054
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:493
std::unique_ptr< std::map< unsigned int, types::manifold_id > > vertex_to_manifold_id_map_1d
Definition: tria.h:4038
std::vector< bool > vertices_used
Definition: tria.h:3968
#define DeclException0(Exception0)
Definition: exceptions.h:470
virtual types::coarse_cell_id coarse_cell_index_to_coarse_cell_id(const unsigned int coarse_cell_index) const
boost::signals2::signal< void(const typename Triangulation< dim, spacedim >::cell_iterator &cell)> post_refinement_on_cell
Definition: tria.h:2134
typename IteratorSelector::active_quad_iterator active_quad_iterator
Definition: tria.h:1487
const std::vector< Point< spacedim > > & get_vertices() const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:394
boost::signals2::signal< void()> mesh_movement
Definition: tria.h:2113
unsigned int level
Definition: grid_out.cc:4578
typename IteratorSelector::raw_quad_iterator raw_quad_iterator
Definition: tria.h:3604
VectorType::value_type * end(VectorType &V)
std::unique_ptr<::internal::TriangulationImplementation::TriaFaces > faces
Definition: tria.h:3957
Point< 3 > vertices[4]
boost::signals2::signal< void()> pre_distributed_save
Definition: tria.h:2250
void save(Archive &ar, const unsigned int version) const
unsigned int n_active_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:12647
typename IteratorSelector::quad_iterator quad_iterator
Definition: tria.h:1472
TriangulationBase< dim, spacedim > Triangulation
Definition: tria_base.h:395
boost::signals2::signal< void()> create
Definition: tria.h:2080
std::unique_ptr< std::map< unsigned int, types::boundary_id > > vertex_to_boundary_id_map_1d
Definition: tria.h:4015
Definition: hp.h:117
unsigned int n_active_quads() const
Definition: tria.cc:13121
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:12640
Definition: cell_id.h:70
Point< 2 > first
Definition: grid_out.cc:4575
boost::signals2::signal< void()> pre_partition
Definition: tria.h:2104
typename IteratorSelector::active_hex_iterator active_hex_iterator
Definition: tria.h:1507
boost::signals2::signal< void()> post_distributed_save
Definition: tria.h:2257
std::vector< unsigned int > n_active_hexes_level
Definition: tria.h:300
std::map< std::pair< cell_iterator, unsigned int >, std::pair< std::pair< cell_iterator, unsigned int >, std::bitset< 3 > > > periodic_face_map
Definition: tria.h:3583
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:393
virtual unsigned int n_global_levels() const
VectorType::value_type * begin(VectorType &V)
typename IteratorSelector::raw_hex_iterator raw_hex_iterator
Definition: tria.h:3605
boost::signals2::signal< unsigned int(const cell_iterator &, const CellStatus), CellWeightSum< unsigned int > > cell_weight
Definition: tria.h:2200
boost::signals2::signal< void()> pre_distributed_load
Definition: tria.h:2265
std::vector< unsigned int > n_active_lines_level
Definition: tria.h:183
std::vector< unsigned int > n_active_quads_level
Definition: tria.h:241
boost::signals2::signal< void()> post_distributed_load
Definition: tria.h:2272
global_cell_index coarse_cell_id
Definition: types.h:114
boost::signals2::signal< void()> pre_distributed_repartition
Definition: tria.h:2235
bool vertex_used(const unsigned int index) const
boost::signals2::signal< void()> any_change
Definition: tria.h:2170
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
boost::signals2::signal< void()> post_distributed_refinement
Definition: tria.h:2223