Reference documentation for deal.II version Git 6cf3cbfd16 2020-07-12 08:18:15 -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.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 
31 
32 #include <boost/serialization/map.hpp>
33 #include <boost/serialization/split_member.hpp>
34 #include <boost/serialization/unique_ptr.hpp>
35 #include <boost/serialization/vector.hpp>
36 #include <boost/signals2.hpp>
37 
38 #include <bitset>
39 #include <functional>
40 #include <list>
41 #include <map>
42 #include <memory>
43 #include <numeric>
44 #include <vector>
45 
46 
48 
49 // Forward declarations
50 #ifndef DOXYGEN
51 template <int dim, int spacedim>
52 class Manifold;
53 
54 template <int dim>
55 struct CellData;
56 
57 struct SubCellData;
58 
60 {
61  template <int, int>
62  struct Description;
63 }
64 
65 namespace GridTools
66 {
67  template <typename CellIterator>
68  struct PeriodicFacePair;
69 }
70 
71 template <int, int, int>
72 class TriaAccessor;
73 template <int spacedim>
74 class TriaAccessor<0, 1, spacedim>;
75 template <int, int, int>
76 class TriaAccessorBase;
77 
78 namespace internal
79 {
80  namespace TriangulationImplementation
81  {
82  class TriaFaces;
83 
84  class TriaObjects;
85 
91  struct Implementation;
92  } // namespace TriangulationImplementation
93 
94  namespace TriaAccessorImplementation
95  {
96  struct Implementation;
97  }
98 } // namespace internal
99 
100 namespace hp
101 {
102  template <int dim, int spacedim>
103  class DoFHandler;
104 }
105 #endif
106 
107 
108 /*------------------------------------------------------------------------*/
109 
110 
111 namespace internal
112 {
117  namespace TriangulationImplementation
118  {
130  template <int dim>
131  struct NumberCache
132  {};
133 
145  template <>
146  struct NumberCache<1>
147  {
151  unsigned int n_levels;
152 
156  unsigned int n_lines;
157 
161  std::vector<unsigned int> n_lines_level;
162 
166  unsigned int n_active_lines;
167 
171  std::vector<unsigned int> n_active_lines_level;
172 
176  NumberCache();
177 
182  std::size_t
183  memory_consumption() const;
184 
189  template <class Archive>
190  void
191  serialize(Archive &ar, const unsigned int version);
192  };
193 
194 
207  template <>
208  struct NumberCache<2> : public NumberCache<1>
209  {
213  unsigned int n_quads;
214 
218  std::vector<unsigned int> n_quads_level;
219 
223  unsigned int n_active_quads;
224 
228  std::vector<unsigned int> n_active_quads_level;
229 
233  NumberCache();
234 
239  std::size_t
240  memory_consumption() const;
241 
246  template <class Archive>
247  void
248  serialize(Archive &ar, const unsigned int version);
249  };
250 
251 
265  template <>
266  struct NumberCache<3> : public NumberCache<2>
267  {
271  unsigned int n_hexes;
272 
276  std::vector<unsigned int> n_hexes_level;
277 
281  unsigned int n_active_hexes;
282 
286  std::vector<unsigned int> n_active_hexes_level;
287 
291  NumberCache();
292 
297  std::size_t
298  memory_consumption() const;
299 
304  template <class Archive>
305  void
306  serialize(Archive &ar, const unsigned int version);
307  };
308  } // namespace TriangulationImplementation
309 } // namespace internal
310 
311 
312 /*------------------------------------------------------------------------*/
313 
314 
1096 template <int dim, int spacedim = dim>
1098 {
1099 private:
1104  using IteratorSelector =
1106 
1107 public:
1113  {
1118  none = 0x0,
1162  limit_level_difference_at_vertices = 0x1,
1183  eliminate_unrefined_islands = 0x2,
1199  patch_level_1 = 0x4,
1223  coarsest_level_1 = 0x8,
1248  allow_anisotropic_smoothing = 0x10,
1281  eliminate_refined_inner_islands = 0x100,
1286  eliminate_refined_boundary_islands = 0x200,
1292  do_not_produce_unrefined_islands = 0x400,
1293 
1298  smoothing_on_refinement =
1299  (limit_level_difference_at_vertices | eliminate_unrefined_islands),
1304  smoothing_on_coarsening =
1305  (eliminate_refined_inner_islands | eliminate_refined_boundary_islands |
1306  do_not_produce_unrefined_islands),
1307 
1313  maximum_smoothing = 0xffff ^ allow_anisotropic_smoothing
1314  };
1315 
1332 
1338 
1356 
1370  using face_iterator = TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>;
1371 
1383  using active_face_iterator =
1384  TriaActiveIterator<TriaAccessor<dim - 1, dim, spacedim>>;
1385 
1395 
1408  using active_vertex_iterator =
1410 
1418  using line_iterator = typename IteratorSelector::line_iterator;
1419 
1433  using active_line_iterator = typename IteratorSelector::active_line_iterator;
1434 
1442  using quad_iterator = typename IteratorSelector::quad_iterator;
1443 
1457  using active_quad_iterator = typename IteratorSelector::active_quad_iterator;
1458 
1466  using hex_iterator = typename IteratorSelector::hex_iterator;
1467 
1477  using active_hex_iterator = typename IteratorSelector::active_hex_iterator;
1478 
1498  {
1505  virtual ~DistortedCellList() noexcept override = default;
1506 
1511  std::list<typename Triangulation<dim, spacedim>::cell_iterator>
1513  };
1514 
1518  static const unsigned int dimension = dim;
1519 
1523  static const unsigned int space_dimension = spacedim;
1524 
1539  Triangulation(const MeshSmoothing smooth_grid = none,
1540  const bool check_for_distorted_cells = false);
1541 
1557  Triangulation(const Triangulation<dim, spacedim> &) = delete;
1558 
1566 
1570  Triangulation &
1571  operator=(Triangulation<dim, spacedim> &&tria) noexcept;
1572 
1576  virtual ~Triangulation() override;
1577 
1584  virtual void
1585  clear();
1586 
1592  virtual void
1593  set_mesh_smoothing(const MeshSmoothing mesh_smoothing);
1594 
1598  virtual const MeshSmoothing &
1599  get_mesh_smoothing() const;
1600 
1623  void
1624  set_manifold(const types::manifold_id number,
1625  const Manifold<dim, spacedim> &manifold_object);
1626 
1639  void
1640  reset_manifold(const types::manifold_id manifold_number);
1641 
1653  void
1654  reset_all_manifolds();
1655 
1664  void
1665  set_all_manifold_ids(const types::manifold_id number);
1666 
1675  void
1676  set_all_manifold_ids_on_boundary(const types::manifold_id number);
1677 
1687  void
1688  set_all_manifold_ids_on_boundary(const types::boundary_id b_id,
1689  const types::manifold_id number);
1690 
1700  const Manifold<dim, spacedim> &
1701  get_manifold(const types::manifold_id number) const;
1702 
1715  virtual std::vector<types::boundary_id>
1716  get_boundary_ids() const;
1717 
1730  virtual std::vector<types::manifold_id>
1731  get_manifold_ids() const;
1732 
1759  virtual void
1760  copy_triangulation(const Triangulation<dim, spacedim> &other_tria);
1761 
1807  virtual void
1808  create_triangulation(const std::vector<Point<spacedim>> &vertices,
1809  const std::vector<CellData<dim>> & cells,
1810  const SubCellData & subcelldata);
1811 
1821  virtual void
1824  &construction_data);
1825 
1834  virtual void
1835  create_triangulation_compatibility(
1836  const std::vector<Point<spacedim>> &vertices,
1837  const std::vector<CellData<dim>> & cells,
1838  const SubCellData & subcelldata);
1839 
1846  void
1847  flip_all_direction_flags();
1848 
1860  void
1861  set_all_refine_flags();
1862 
1877  void
1878  refine_global(const unsigned int times = 1);
1879 
1911  virtual void
1912  execute_coarsening_and_refinement();
1913 
1943  virtual bool
1944  prepare_coarsening_and_refinement();
1945 
1946  /*
1947  * @}
1948  */
1949 
1964  {
1981  CELL_INVALID
1982  };
1983 
1989  template <typename T>
1991  {
1992  using result_type = T;
1993 
1994  template <typename InputIterator>
1995  T
1996  operator()(InputIterator first, InputIterator last) const
1997  {
1998  return std::accumulate(first, last, T());
1999  }
2000  };
2001 
2011  struct Signals
2012  {
2020  boost::signals2::signal<void()> create;
2021 
2029  boost::signals2::signal<void()> pre_refinement;
2030 
2036  boost::signals2::signal<void()> post_refinement;
2037 
2044  boost::signals2::signal<void()> pre_partition;
2045 
2053  boost::signals2::signal<void()> mesh_movement;
2054 
2062  boost::signals2::signal<void(
2063  const typename Triangulation<dim, spacedim>::cell_iterator &cell)>
2065 
2072  boost::signals2::signal<void(
2073  const typename Triangulation<dim, spacedim>::cell_iterator &cell)>
2075 
2082  boost::signals2::signal<void(
2083  const Triangulation<dim, spacedim> &destination_tria)>
2085 
2099  boost::signals2::signal<void()> clear;
2100 
2110  boost::signals2::signal<void()> any_change;
2111 
2135  boost::signals2::signal<unsigned int(const cell_iterator &,
2136  const CellStatus),
2139 
2151  boost::signals2::signal<void()> pre_distributed_refinement;
2152 
2161  boost::signals2::signal<void()> post_distributed_refinement;
2162 
2173  boost::signals2::signal<void()> pre_distributed_repartition;
2174 
2180  boost::signals2::signal<void()> post_distributed_repartition;
2181 
2188  boost::signals2::signal<void()> pre_distributed_save;
2189 
2195  boost::signals2::signal<void()> post_distributed_save;
2196 
2203  boost::signals2::signal<void()> pre_distributed_load;
2204 
2210  boost::signals2::signal<void()> post_distributed_load;
2211  };
2212 
2216  mutable Signals signals;
2217 
2218  /*
2219  * @}
2220  */
2221 
2231  void
2232  save_refine_flags(std::ostream &out) const;
2233 
2237  void
2238  save_refine_flags(std::vector<bool> &v) const;
2239 
2243  void
2244  load_refine_flags(std::istream &in);
2245 
2249  void
2250  load_refine_flags(const std::vector<bool> &v);
2251 
2255  void
2256  save_coarsen_flags(std::ostream &out) const;
2257 
2261  void
2262  save_coarsen_flags(std::vector<bool> &v) const;
2263 
2267  void
2268  load_coarsen_flags(std::istream &out);
2269 
2273  void
2274  load_coarsen_flags(const std::vector<bool> &v);
2275 
2280  bool
2281  get_anisotropic_refinement_flag() const;
2282 
2283  /*
2284  * @}
2285  */
2286 
2296  void
2297  clear_user_flags();
2298 
2304  void
2305  save_user_flags(std::ostream &out) const;
2306 
2312  void
2313  save_user_flags(std::vector<bool> &v) const;
2314 
2319  void
2320  load_user_flags(std::istream &in);
2321 
2326  void
2327  load_user_flags(const std::vector<bool> &v);
2328 
2333  void
2334  clear_user_flags_line();
2335 
2340  void
2341  save_user_flags_line(std::ostream &out) const;
2342 
2348  void
2349  save_user_flags_line(std::vector<bool> &v) const;
2350 
2355  void
2356  load_user_flags_line(std::istream &in);
2357 
2362  void
2363  load_user_flags_line(const std::vector<bool> &v);
2364 
2369  void
2370  clear_user_flags_quad();
2371 
2376  void
2377  save_user_flags_quad(std::ostream &out) const;
2378 
2384  void
2385  save_user_flags_quad(std::vector<bool> &v) const;
2386 
2391  void
2392  load_user_flags_quad(std::istream &in);
2393 
2398  void
2399  load_user_flags_quad(const std::vector<bool> &v);
2400 
2401 
2406  void
2407  clear_user_flags_hex();
2408 
2413  void
2414  save_user_flags_hex(std::ostream &out) const;
2415 
2421  void
2422  save_user_flags_hex(std::vector<bool> &v) const;
2423 
2428  void
2429  load_user_flags_hex(std::istream &in);
2430 
2435  void
2436  load_user_flags_hex(const std::vector<bool> &v);
2437 
2443  void
2444  clear_user_data();
2445 
2451  void
2452  save_user_indices(std::vector<unsigned int> &v) const;
2453 
2458  void
2459  load_user_indices(const std::vector<unsigned int> &v);
2460 
2466  void
2467  save_user_pointers(std::vector<void *> &v) const;
2468 
2473  void
2474  load_user_pointers(const std::vector<void *> &v);
2475 
2481  void
2482  save_user_indices_line(std::vector<unsigned int> &v) const;
2483 
2488  void
2489  load_user_indices_line(const std::vector<unsigned int> &v);
2490 
2496  void
2497  save_user_indices_quad(std::vector<unsigned int> &v) const;
2498 
2503  void
2504  load_user_indices_quad(const std::vector<unsigned int> &v);
2505 
2511  void
2512  save_user_indices_hex(std::vector<unsigned int> &v) const;
2513 
2518  void
2519  load_user_indices_hex(const std::vector<unsigned int> &v);
2525  void
2526  save_user_pointers_line(std::vector<void *> &v) const;
2527 
2532  void
2533  load_user_pointers_line(const std::vector<void *> &v);
2534 
2540  void
2541  save_user_pointers_quad(std::vector<void *> &v) const;
2542 
2547  void
2548  load_user_pointers_quad(const std::vector<void *> &v);
2549 
2555  void
2556  save_user_pointers_hex(std::vector<void *> &v) const;
2557 
2562  void
2563  load_user_pointers_hex(const std::vector<void *> &v);
2564 
2565  /*
2566  * @}
2567  */
2568 
2592  begin(const unsigned int level = 0) const;
2593 
2625  begin_active(const unsigned int level = 0) const;
2626 
2632  end() const;
2633 
2653  end(const unsigned int level) const;
2654 
2675  end_active(const unsigned int level) const;
2676 
2677 
2682  last() const;
2683 
2688  last_active() const;
2689 
2705  cell_iterators() const;
2706 
2743  active_cell_iterators() const;
2744 
2761  cell_iterators_on_level(const unsigned int level) const;
2762 
2779  active_cell_iterators_on_level(const unsigned int level) const;
2780 
2781  /*
2782  * @}
2783  */
2784 
2785  /*-------------------------------------------------------------------------*/
2786 
2796  begin_face() const;
2797 
2802  begin_active_face() const;
2803 
2809  end_face() const;
2810 
2830  active_face_iterators() const;
2831 
2832  /*
2833  * @}
2834  */
2835 
2836  /*-------------------------------------------------------------------------*/
2837 
2848  begin_vertex() const;
2849 
2856  begin_active_vertex() const;
2857 
2864  end_vertex() const;
2865 
2866  /*
2867  * @}
2868  */
2869 
2887  unsigned int
2888  n_lines() const;
2889 
2893  unsigned int
2894  n_lines(const unsigned int level) const;
2895 
2899  unsigned int
2900  n_active_lines() const;
2901 
2905  unsigned int
2906  n_active_lines(const unsigned int level) const;
2907 
2911  unsigned int
2912  n_quads() const;
2913 
2917  unsigned int
2918  n_quads(const unsigned int level) const;
2919 
2923  unsigned int
2924  n_active_quads() const;
2925 
2929  unsigned int
2930  n_active_quads(const unsigned int level) const;
2931 
2935  unsigned int
2936  n_hexs() const;
2937 
2942  unsigned int
2943  n_hexs(const unsigned int level) const;
2944 
2948  unsigned int
2949  n_active_hexs() const;
2950 
2955  unsigned int
2956  n_active_hexs(const unsigned int level) const;
2957 
2962  unsigned int
2963  n_cells() const;
2964 
2969  unsigned int
2970  n_cells(const unsigned int level) const;
2971 
2976  unsigned int
2977  n_active_cells() const;
2978 
2986  virtual types::global_cell_index
2987  n_global_active_cells() const;
2988 
2989 
2994  unsigned int
2995  n_active_cells(const unsigned int level) const;
2996 
3002  unsigned int
3003  n_faces() const;
3004 
3010  unsigned int
3011  n_active_faces() const;
3012 
3030  unsigned int
3031  n_levels() const;
3032 
3039  virtual unsigned int
3040  n_global_levels() const;
3041 
3051  virtual bool
3052  has_hanging_nodes() const;
3053 
3061  unsigned int
3062  n_vertices() const;
3063 
3072  const std::vector<Point<spacedim>> &
3073  get_vertices() const;
3074 
3079  unsigned int
3080  n_used_vertices() const;
3081 
3085  bool
3086  vertex_used(const unsigned int index) const;
3087 
3092  const std::vector<bool> &
3093  get_used_vertices() const;
3094 
3106  unsigned int
3107  max_adjacent_cells() const;
3108 
3115  virtual types::subdomain_id
3116  locally_owned_subdomain() const;
3117 
3128  get_triangulation();
3129 
3135  get_triangulation() const;
3136 
3137 
3138  /*
3139  * @}
3140  */
3141 
3156  unsigned int
3157  n_raw_lines() const;
3158 
3168  unsigned int
3169  n_raw_lines(const unsigned int level) const;
3170 
3180  unsigned int
3181  n_raw_quads() const;
3182 
3192  unsigned int
3193  n_raw_quads(const unsigned int level) const;
3194 
3204  unsigned int
3205  n_raw_hexs(const unsigned int level) const;
3206 
3216  unsigned int
3217  n_raw_cells(const unsigned int level) const;
3218 
3230  unsigned int
3231  n_raw_faces() const;
3232 
3233  /*
3234  * @}
3235  */
3236 
3245  virtual std::size_t
3246  memory_consumption() const;
3247 
3256  template <class Archive>
3257  void
3258  save(Archive &ar, const unsigned int version) const;
3259 
3275  template <class Archive>
3276  void
3277  load(Archive &ar, const unsigned int version);
3278 
3279 
3295  virtual void
3296  add_periodicity(
3297  const std::vector<GridTools::PeriodicFacePair<cell_iterator>> &);
3298 
3302  const std::map<
3303  std::pair<cell_iterator, unsigned int>,
3304  std::pair<std::pair<cell_iterator, unsigned int>, std::bitset<3>>> &
3305  get_periodic_face_map() const;
3306 
3307 #ifdef DOXYGEN
3308 
3312  template <class Archive>
3313  void
3314  serialize(Archive &archive, const unsigned int version);
3315 #else
3316  // This macro defines the serialize() method that is compatible with
3317  // the templated save() and load() method that have been implemented.
3318  BOOST_SERIALIZATION_SPLIT_MEMBER()
3319 #endif
3320 
3331  DeclException2(ExcInvalidLevel,
3332  int,
3333  int,
3334  << "You are requesting information from refinement level "
3335  << arg1
3336  << " of a triangulation, but this triangulation only has "
3337  << arg2 << " refinement levels. The given level " << arg1
3338  << " must be *less* than " << arg2 << ".");
3346  ExcTriangulationNotEmpty,
3347  int,
3348  int,
3349  << "You are trying to perform an operation on a triangulation "
3350  << "that is only allowed if the triangulation is currently empty. "
3351  << "However, it currently stores " << arg1 << " vertices and has "
3352  << "cells on " << arg2 << " levels.");
3358  DeclException0(ExcGridReadError);
3369  DeclException1(ExcEmptyLevel,
3370  int,
3371  << "You tried to do something on level " << arg1
3372  << ", but this level is empty.");
3378  DeclException0(ExcNonOrientableTriangulation);
3379 
3387  DeclException1(ExcBoundaryIdNotFound,
3389  << "The given boundary_id " << arg1
3390  << " is not defined in this Triangulation!");
3391 
3398  ExcInconsistentCoarseningFlags,
3399  "A cell is flagged for coarsening, but either not all of its siblings "
3400  "are active or flagged for coarsening as well. Please clean up all "
3401  "coarsen flags on your triangulation via "
3402  "Triangulation::prepare_coarsening_and_refinement() beforehand!");
3403 
3404  /*
3405  * @}
3406  */
3407 
3408 protected:
3414 
3428  static void
3429  write_bool_vector(const unsigned int magic_number1,
3430  const std::vector<bool> &v,
3431  const unsigned int magic_number2,
3432  std::ostream & out);
3433 
3438  static void
3439  read_bool_vector(const unsigned int magic_number1,
3440  std::vector<bool> &v,
3441  const unsigned int magic_number2,
3442  std::istream & in);
3443 
3448  void
3449  update_periodic_face_map();
3450 
3451 
3452 private:
3459  std::vector<GridTools::PeriodicFacePair<cell_iterator>>
3461 
3466  std::map<std::pair<cell_iterator, unsigned int>,
3467  std::pair<std::pair<cell_iterator, unsigned int>, std::bitset<3>>>
3469 
3484  using raw_face_iterator =
3485  TriaRawIterator<TriaAccessor<dim - 1, dim, spacedim>>;
3486  using raw_vertex_iterator =
3488  using raw_line_iterator = typename IteratorSelector::raw_line_iterator;
3489  using raw_quad_iterator = typename IteratorSelector::raw_quad_iterator;
3490  using raw_hex_iterator = typename IteratorSelector::raw_hex_iterator;
3491 
3497  begin_raw(const unsigned int level = 0) const;
3498 
3504  end_raw(const unsigned int level) const;
3505 
3506  /*
3507  * @}
3508  */
3509 
3522  begin_raw_line(const unsigned int level = 0) const;
3523 
3542  begin_line(const unsigned int level = 0) const;
3543 
3562  begin_active_line(const unsigned int level = 0) const;
3563 
3569  end_line() const;
3570 
3571  /*
3572  * @}
3573  */
3574 
3600  begin_raw_quad(const unsigned int level = 0) const;
3601 
3620  begin_quad(const unsigned int level = 0) const;
3621 
3640  begin_active_quad(const unsigned int level = 0) const;
3641 
3647  end_quad() const;
3648 
3649  /*
3650  * @}
3651  */
3652 
3677  begin_raw_hex(const unsigned int level = 0) const;
3678 
3696  hex_iterator
3697  begin_hex(const unsigned int level = 0) const;
3698 
3717  begin_active_hex(const unsigned int level = 0) const;
3718 
3723  hex_iterator
3724  end_hex() const;
3725 
3726  /*
3727  * @}
3728  */
3729 
3730 
3744  void
3745  clear_despite_subscriptions();
3746 
3753  void
3754  reset_active_cell_indices();
3755 
3769  DistortedCellList
3770  execute_refinement();
3771 
3778  void
3779  execute_coarsening();
3780 
3785  void
3786  fix_coarsen_flags();
3787 
3801  virtual unsigned int
3802  coarse_cell_id_to_coarse_cell_index(
3804 
3805 
3818  virtual types::coarse_cell_id
3819  coarse_cell_index_to_coarse_cell_id(
3820  const unsigned int coarse_cell_index) const;
3821 
3826  std::vector<
3827  std::unique_ptr<::internal::TriangulationImplementation::TriaLevel>>
3829 
3835  std::unique_ptr<::internal::TriangulationImplementation::TriaFaces>
3837 
3838 
3842  std::vector<Point<spacedim>> vertices;
3843 
3847  std::vector<bool> vertices_used;
3848 
3853  std::map<types::manifold_id, std::unique_ptr<const Manifold<dim, spacedim>>>
3855 
3860 
3861 
3867 
3878 
3893  std::unique_ptr<std::map<unsigned int, types::boundary_id>>
3895 
3896 
3916  std::unique_ptr<std::map<unsigned int, types::manifold_id>>
3918 
3919  // make a couple of classes friends
3920  template <int, int, int>
3921  friend class TriaAccessorBase;
3922  template <int, int, int>
3923  friend class TriaAccessor;
3924  friend class TriaAccessor<0, 1, spacedim>;
3925 
3926  friend class CellAccessor<dim, spacedim>;
3927 
3928  friend struct ::internal::TriaAccessorImplementation::Implementation;
3929 
3930  friend class hp::DoFHandler<dim, spacedim>;
3931 
3932  friend struct ::internal::TriangulationImplementation::Implementation;
3933 
3934  friend class ::internal::TriangulationImplementation::TriaObjects;
3935 
3936  friend class CellId;
3937 
3938  // explicitly check for sensible template arguments, but not on windows
3939  // because MSVC creates bogus warnings during normal compilation
3940 #ifndef DEAL_II_MSVC
3941  static_assert(dim <= spacedim,
3942  "The dimension <dim> of a Triangulation must be less than or "
3943  "equal to the space dimension <spacedim> in which it lives.");
3944 #endif
3945 };
3946 
3947 
3948 #ifndef DOXYGEN
3949 
3950 
3951 
3952 namespace internal
3953 {
3954  namespace TriangulationImplementation
3955  {
3956  template <class Archive>
3957  void
3958  NumberCache<1>::serialize(Archive &ar, const unsigned int)
3959  {
3960  ar &n_levels;
3961  ar &n_lines &n_lines_level;
3962  ar &n_active_lines &n_active_lines_level;
3963  }
3964 
3965 
3966  template <class Archive>
3967  void
3968  NumberCache<2>::serialize(Archive &ar, const unsigned int version)
3969  {
3970  this->NumberCache<1>::serialize(ar, version);
3971 
3972  ar &n_quads &n_quads_level;
3973  ar &n_active_quads &n_active_quads_level;
3974  }
3975 
3976 
3977  template <class Archive>
3978  void
3979  NumberCache<3>::serialize(Archive &ar, const unsigned int version)
3980  {
3981  this->NumberCache<2>::serialize(ar, version);
3982 
3983  ar &n_hexes &n_hexes_level;
3984  ar &n_active_hexes &n_active_hexes_level;
3985  }
3986 
3987  } // namespace TriangulationImplementation
3988 } // namespace internal
3989 
3990 
3991 template <int dim, int spacedim>
3992 inline bool
3993 Triangulation<dim, spacedim>::vertex_used(const unsigned int index) const
3994 {
3995  AssertIndexRange(index, vertices_used.size());
3996  return vertices_used[index];
3997 }
3998 
3999 
4000 
4001 template <int dim, int spacedim>
4002 inline unsigned int
4004 {
4005  return number_cache.n_levels;
4006 }
4007 
4008 template <int dim, int spacedim>
4009 inline unsigned int
4011 {
4012  return number_cache.n_levels;
4013 }
4014 
4015 
4016 template <int dim, int spacedim>
4017 inline unsigned int
4019 {
4020  return vertices.size();
4021 }
4022 
4023 
4024 
4025 template <int dim, int spacedim>
4026 inline const std::vector<Point<spacedim>> &
4028 {
4029  return vertices;
4030 }
4031 
4032 
4033 template <int dim, int spacedim>
4034 template <class Archive>
4035 void
4036 Triangulation<dim, spacedim>::save(Archive &ar, const unsigned int) const
4037 {
4038  // as discussed in the documentation, do not store the signals as
4039  // well as boundary and manifold description but everything else
4040  ar &smooth_grid;
4041 
4042  unsigned int n_levels = levels.size();
4043  ar & n_levels;
4044  for (const auto &level : levels)
4045  ar &level;
4046 
4047  // boost dereferences a nullptr when serializing a nullptr
4048  // at least up to 1.65.1. This causes problems with clang-5.
4049  // Therefore, work around it.
4050  bool faces_is_nullptr = (faces.get() == nullptr);
4051  ar & faces_is_nullptr;
4052  if (!faces_is_nullptr)
4053  ar &faces;
4054 
4055  ar &vertices;
4056  ar &vertices_used;
4057 
4058  ar &anisotropic_refinement;
4059  ar &number_cache;
4060 
4061  ar &check_for_distorted_cells;
4062 
4063  if (dim == 1)
4064  {
4065  ar &vertex_to_boundary_id_map_1d;
4066  ar &vertex_to_manifold_id_map_1d;
4067  }
4068 }
4069 
4070 
4071 
4072 template <int dim, int spacedim>
4073 template <class Archive>
4074 void
4075 Triangulation<dim, spacedim>::load(Archive &ar, const unsigned int)
4076 {
4077  // clear previous content. this also calls the respective signal
4078  clear();
4079 
4080  // as discussed in the documentation, do not store the signals as
4081  // well as boundary and manifold description but everything else
4082  ar &smooth_grid;
4083 
4084  unsigned int size;
4085  ar & size;
4086  levels.resize(size);
4087  for (auto &level_ : levels)
4088  {
4089  std::unique_ptr<internal::TriangulationImplementation::TriaLevel> level;
4090  ar & level;
4091  level_ = std::move(level);
4092  }
4093 
4094  // Workaround for nullptr, see in save().
4095  bool faces_is_nullptr = true;
4096  ar & faces_is_nullptr;
4097  if (!faces_is_nullptr)
4098  ar &faces;
4099 
4100  ar &vertices;
4101  ar &vertices_used;
4102 
4103  ar &anisotropic_refinement;
4104  ar &number_cache;
4105 
4106  // the levels do not serialize the active_cell_indices because
4107  // they are easy enough to rebuild upon re-loading data. do
4108  // this here. don't forget to first resize the fields appropriately
4109  {
4110  for (auto &level : levels)
4111  level->active_cell_indices.resize(level->refine_flags.size());
4112  reset_active_cell_indices();
4113  }
4114 
4115 
4116  bool my_check_for_distorted_cells;
4117  ar & my_check_for_distorted_cells;
4118 
4119  Assert(my_check_for_distorted_cells == check_for_distorted_cells,
4120  ExcMessage("The triangulation loaded into here must have the "
4121  "same setting with regard to reporting distorted "
4122  "cell as the one previously stored."));
4123 
4124  if (dim == 1)
4125  {
4126  ar &vertex_to_boundary_id_map_1d;
4127  ar &vertex_to_manifold_id_map_1d;
4128  }
4129 
4130  // trigger the create signal to indicate
4131  // that new content has been imported into
4132  // the triangulation
4133  signals.create();
4134 }
4135 
4136 
4137 
4138 template <int dim, int spacedim>
4139 inline unsigned int
4142 {
4143  return coarse_cell_id;
4144 }
4145 
4146 
4147 
4148 template <int dim, int spacedim>
4149 inline types::coarse_cell_id
4151  const unsigned int coarse_cell_index) const
4152 {
4153  return coarse_cell_index;
4154 }
4155 
4156 
4157 
4158 /* -------------- declaration of explicit specializations ------------- */
4159 
4160 template <>
4161 unsigned int
4163 template <>
4164 unsigned int
4165 Triangulation<1, 1>::n_quads(const unsigned int level) const;
4166 template <>
4167 unsigned int
4168 Triangulation<1, 1>::n_raw_quads(const unsigned int level) const;
4169 template <>
4170 unsigned int
4171 Triangulation<2, 2>::n_raw_quads(const unsigned int level) const;
4172 template <>
4173 unsigned int
4174 Triangulation<1, 1>::n_raw_hexs(const unsigned int level) const;
4175 template <>
4176 unsigned int
4177 Triangulation<1, 1>::n_active_quads(const unsigned int level) const;
4178 template <>
4179 unsigned int
4181 template <>
4182 unsigned int
4184 
4185 
4186 // -------------------------------------------------------------------
4187 // -- Explicit specializations for codimension one grids
4188 
4189 
4190 template <>
4191 unsigned int
4193 template <>
4194 unsigned int
4195 Triangulation<1, 2>::n_quads(const unsigned int level) const;
4196 template <>
4197 unsigned int
4198 Triangulation<1, 2>::n_raw_quads(const unsigned int level) const;
4199 template <>
4200 unsigned int
4201 Triangulation<2, 3>::n_raw_quads(const unsigned int level) const;
4202 template <>
4203 unsigned int
4204 Triangulation<1, 2>::n_raw_hexs(const unsigned int level) const;
4205 template <>
4206 unsigned int
4207 Triangulation<1, 2>::n_active_quads(const unsigned int level) const;
4208 template <>
4209 unsigned int
4211 template <>
4212 unsigned int
4214 
4215 // -------------------------------------------------------------------
4216 // -- Explicit specializations for codimension two grids
4217 
4218 
4219 template <>
4220 unsigned int
4222 template <>
4223 unsigned int
4224 Triangulation<1, 3>::n_quads(const unsigned int level) const;
4225 template <>
4226 unsigned int
4227 Triangulation<1, 3>::n_raw_quads(const unsigned int level) const;
4228 template <>
4229 unsigned int
4230 Triangulation<2, 3>::n_raw_quads(const unsigned int level) const;
4231 template <>
4232 unsigned int
4233 Triangulation<1, 3>::n_raw_hexs(const unsigned int level) const;
4234 template <>
4235 unsigned int
4236 Triangulation<1, 3>::n_active_quads(const unsigned int level) const;
4237 template <>
4238 unsigned int
4240 template <>
4241 unsigned int
4243 
4244 
4245 #endif // DOXYGEN
4246 
4248 
4249 // Include tria_accessor.h here, so that it is possible for an end
4250 // user to use the iterators of Triangulation<dim> directly without
4251 // the need to include tria_accessor.h separately. (Otherwise the
4252 // iterators are an 'opaque' or 'incomplete' type.)
4254 
4255 #endif
boost::signals2::signal< void(const Triangulation< dim, spacedim > &destination_tria)> copy
Definition: tria.h:2084
::internal::TriangulationImplementation::NumberCache< dim > number_cache
Definition: tria.h:3877
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:2036
unsigned int n_vertices() const
boost::signals2::signal< void(const typename Triangulation< dim, spacedim >::cell_iterator &cell)> pre_coarsening_on_cell
Definition: tria.h:2064
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:538
boost::signals2::signal< void()> post_distributed_repartition
Definition: tria.h:2180
typename IteratorSelector::line_iterator line_iterator
Definition: tria.h:1418
T operator()(InputIterator first, InputIterator last) const
Definition: tria.h:1996
std::vector< GridTools::PeriodicFacePair< cell_iterator > > periodic_face_pairs_level_0
Definition: tria.h:3460
typename IteratorSelector::raw_line_iterator raw_line_iterator
Definition: tria.h:3488
std::vector< Point< spacedim > > vertices
Definition: tria.h:3842
std::map< types::manifold_id, std::unique_ptr< const Manifold< dim, spacedim > > > manifold
Definition: tria.h:3854
typename IteratorSelector::hex_iterator hex_iterator
Definition: tria.h:1466
typename IteratorSelector::active_line_iterator active_line_iterator
Definition: tria.h:1433
boost::signals2::signal< void()> pre_distributed_refinement
Definition: tria.h:2151
bool anisotropic_refinement
Definition: tria.h:3859
#define AssertIndexRange(index, range)
Definition: exceptions.h:1628
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:2029
const bool check_for_distorted_cells
Definition: tria.h:3866
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:515
static const char T
unsigned int n_raw_quads() const
Definition: tria.cc:13311
#define Assert(cond, exc)
Definition: exceptions.h:1403
Signals signals
Definition: tria.h:2216
std::vector< std::unique_ptr<::internal::TriangulationImplementation::TriaLevel > > levels
Definition: tria.h:3828
boost::signals2::signal< void()> clear
Definition: tria.h:2099
unsigned int n_raw_hexs(const unsigned int level) const
Definition: tria.cc:13367
unsigned int max_adjacent_cells() const
Definition: tria.cc:13484
void load(Archive &ar, const unsigned int version)
std::list< typename Triangulation< dim, spacedim >::cell_iterator > distorted_cells
Definition: tria.h:1512
unsigned int n_quads() const
Definition: tria.cc:13263
#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:3917
std::vector< bool > vertices_used
Definition: tria.h:3847
#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:2074
typename IteratorSelector::active_quad_iterator active_quad_iterator
Definition: tria.h:1457
const std::vector< Point< spacedim > > & get_vertices() const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
boost::signals2::signal< void()> mesh_movement
Definition: tria.h:2053
unsigned int level
Definition: grid_out.cc:4339
typename IteratorSelector::raw_quad_iterator raw_quad_iterator
Definition: tria.h:3489
VectorType::value_type * end(VectorType &V)
std::unique_ptr<::internal::TriangulationImplementation::TriaFaces > faces
Definition: tria.h:3836
Point< 3 > vertices[4]
boost::signals2::signal< void()> pre_distributed_save
Definition: tria.h:2188
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:3468
void save(Archive &ar, const unsigned int version) const
unsigned int n_active_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:12856
typename IteratorSelector::quad_iterator quad_iterator
Definition: tria.h:1442
TriangulationBase< dim, spacedim > Triangulation
Definition: tria_base.h:302
boost::signals2::signal< void()> create
Definition: tria.h:2020
std::unique_ptr< std::map< unsigned int, types::boundary_id > > vertex_to_boundary_id_map_1d
Definition: tria.h:3894
Definition: hp.h:117
unsigned int n_active_quads() const
Definition: tria.cc:13330
unsigned int n_cells(const internal::TriangulationImplementation::NumberCache< 1 > &c)
Definition: tria.cc:12849
Definition: cell_id.h:69
Point< 2 > first
Definition: grid_out.cc:4336
boost::signals2::signal< void()> pre_partition
Definition: tria.h:2044
typename IteratorSelector::active_hex_iterator active_hex_iterator
Definition: tria.h:1477
boost::signals2::signal< void()> post_distributed_save
Definition: tria.h:2195
std::vector< unsigned int > n_active_hexes_level
Definition: tria.h:286
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
virtual unsigned int n_global_levels() const
VectorType::value_type * begin(VectorType &V)
typename IteratorSelector::raw_hex_iterator raw_hex_iterator
Definition: tria.h:3490
MeshSmoothing smooth_grid
Definition: tria.h:3413
boost::signals2::signal< unsigned int(const cell_iterator &, const CellStatus), CellWeightSum< unsigned int > > cell_weight
Definition: tria.h:2138
boost::signals2::signal< void()> pre_distributed_load
Definition: tria.h:2203
std::vector< unsigned int > n_active_lines_level
Definition: tria.h:171
std::vector< unsigned int > n_active_quads_level
Definition: tria.h:228
boost::signals2::signal< void()> post_distributed_load
Definition: tria.h:2210
global_cell_index coarse_cell_id
Definition: types.h:114
boost::signals2::signal< void()> pre_distributed_repartition
Definition: tria.h:2173
bool vertex_used(const unsigned int index) const
boost::signals2::signal< void()> any_change
Definition: tria.h:2110
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:2161