Reference documentation for deal.II version Git 840703ecf3 2019-10-14 17:15:20 -0400
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
tria.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2019 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 
22 #include <deal.II/base/geometry_info.h>
23 #include <deal.II/base/iterator_range.h>
24 #include <deal.II/base/point.h>
25 #include <deal.II/base/smartpointer.h>
26 #include <deal.II/base/subscriptor.h>
27 
28 #include <deal.II/grid/tria_iterator_selector.h>
29 
30 #include <boost/serialization/map.hpp>
31 #include <boost/serialization/split_member.hpp>
32 #include <boost/serialization/unique_ptr.hpp>
33 #include <boost/serialization/vector.hpp>
34 #include <boost/signals2.hpp>
35 
36 #include <bitset>
37 #include <functional>
38 #include <list>
39 #include <map>
40 #include <memory>
41 #include <numeric>
42 #include <vector>
43 
44 
45 DEAL_II_NAMESPACE_OPEN
46 
47 // Forward declarations
48 #ifndef DOXYGEN
49 template <int dim, int spacedim>
50 class Manifold;
51 
52 namespace GridTools
53 {
54  template <typename CellIterator>
55  struct PeriodicFacePair;
56 }
57 
58 template <int, int, int>
59 class TriaAccessor;
60 template <int spacedim>
61 class TriaAccessor<0, 1, spacedim>;
62 template <int, int, int>
63 class TriaAccessorBase;
64 
65 namespace internal
66 {
67  namespace TriangulationImplementation
68  {
69  template <int dim>
70  class TriaLevel;
71  template <int dim>
72  class TriaFaces;
73 
74  template <typename>
75  class TriaObjects;
76 
82  struct Implementation;
83  } // namespace TriangulationImplementation
84 
85  namespace TriaAccessorImplementation
86  {
87  struct Implementation;
88  }
89 } // namespace internal
90 
91 namespace hp
92 {
93  template <int dim, int spacedim>
94  class DoFHandler;
95 }
96 #endif
97 
98 /*------------------------------------------------------------------------*/
99 
135 template <int structdim>
136 struct CellData
137 {
144 
152  union
153  {
163 
174  };
175 
185 
193  CellData();
194 
198  bool
199  operator==(const CellData<structdim> &other) const;
200 
204  template <class Archive>
205  void
206  serialize(Archive &ar, const unsigned int version);
207 
208  static_assert(structdim > 0,
209  "The class CellData can only be used for structdim>0.");
210 };
211 
212 
213 
269 {
276  std::vector<CellData<1>> boundary_lines;
277 
284  std::vector<CellData<2>> boundary_quads;
285 
294  bool
295  check_consistency(const unsigned int dim) const;
296 };
297 
298 
299 /*------------------------------------------------------------------------*/
300 
301 
302 namespace internal
303 {
308  namespace TriangulationImplementation
309  {
323  template <int dim>
324  struct NumberCache
325  {};
326 
340  template <>
341  struct NumberCache<1>
342  {
346  unsigned int n_levels;
347 
351  unsigned int n_lines;
352 
356  std::vector<unsigned int> n_lines_level;
357 
361  unsigned int n_active_lines;
362 
366  std::vector<unsigned int> n_active_lines_level;
367 
371  NumberCache();
372 
377  std::size_t
378  memory_consumption() const;
379 
384  template <class Archive>
385  void
386  serialize(Archive &ar, const unsigned int version);
387  };
388 
389 
404  template <>
405  struct NumberCache<2> : public NumberCache<1>
406  {
410  unsigned int n_quads;
411 
415  std::vector<unsigned int> n_quads_level;
416 
420  unsigned int n_active_quads;
421 
425  std::vector<unsigned int> n_active_quads_level;
426 
430  NumberCache();
431 
436  std::size_t
437  memory_consumption() const;
438 
443  template <class Archive>
444  void
445  serialize(Archive &ar, const unsigned int version);
446  };
447 
448 
464  template <>
465  struct NumberCache<3> : public NumberCache<2>
466  {
470  unsigned int n_hexes;
471 
475  std::vector<unsigned int> n_hexes_level;
476 
480  unsigned int n_active_hexes;
481 
485  std::vector<unsigned int> n_active_hexes_level;
486 
490  NumberCache();
491 
496  std::size_t
497  memory_consumption() const;
498 
503  template <class Archive>
504  void
505  serialize(Archive &ar, const unsigned int version);
506  };
507  } // namespace TriangulationImplementation
508 } // namespace internal
509 
510 
511 /*------------------------------------------------------------------------*/
512 
513 
1294 template <int dim, int spacedim = dim>
1296 {
1297 private:
1302  using IteratorSelector =
1303  ::internal::TriangulationImplementation::Iterators<dim, spacedim>;
1304 
1305 public:
1311  {
1316  none = 0x0,
1360  limit_level_difference_at_vertices = 0x1,
1381  eliminate_unrefined_islands = 0x2,
1397  patch_level_1 = 0x4,
1421  coarsest_level_1 = 0x8,
1446  allow_anisotropic_smoothing = 0x10,
1479  eliminate_refined_inner_islands = 0x100,
1484  eliminate_refined_boundary_islands = 0x200,
1490  do_not_produce_unrefined_islands = 0x400,
1491 
1496  smoothing_on_refinement =
1497  (limit_level_difference_at_vertices | eliminate_unrefined_islands),
1502  smoothing_on_coarsening =
1503  (eliminate_refined_inner_islands | eliminate_refined_boundary_islands |
1504  do_not_produce_unrefined_islands),
1505 
1511  maximum_smoothing = 0xffff ^ allow_anisotropic_smoothing
1512  };
1513 
1530 
1548 
1562  using face_iterator = TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>;
1563 
1575  using active_face_iterator =
1576  TriaActiveIterator<TriaAccessor<dim - 1, dim, spacedim>>;
1577 
1587 
1600  using active_vertex_iterator =
1602 
1610  using line_iterator = typename IteratorSelector::line_iterator;
1611 
1625  using active_line_iterator = typename IteratorSelector::active_line_iterator;
1626 
1634  using quad_iterator = typename IteratorSelector::quad_iterator;
1635 
1649  using active_quad_iterator = typename IteratorSelector::active_quad_iterator;
1650 
1658  using hex_iterator = typename IteratorSelector::hex_iterator;
1659 
1669  using active_hex_iterator = typename IteratorSelector::active_hex_iterator;
1670 
1690  {
1697  virtual ~DistortedCellList() noexcept override = default;
1698 
1703  std::list<typename Triangulation<dim, spacedim>::cell_iterator>
1705  };
1706 
1710  static const unsigned int dimension = dim;
1711 
1715  static const unsigned int space_dimension = spacedim;
1716 
1731  Triangulation(const MeshSmoothing smooth_grid = none,
1732  const bool check_for_distorted_cells = false);
1733 
1749  Triangulation(const Triangulation<dim, spacedim> &) = delete;
1750 
1758 
1762  Triangulation &
1763  operator=(Triangulation<dim, spacedim> &&tria) noexcept;
1764 
1768  virtual ~Triangulation() override;
1769 
1776  virtual void
1777  clear();
1778 
1784  virtual void
1785  set_mesh_smoothing(const MeshSmoothing mesh_smoothing);
1786 
1790  virtual const MeshSmoothing &
1791  get_mesh_smoothing() const;
1792 
1815  void
1816  set_manifold(const types::manifold_id number,
1817  const Manifold<dim, spacedim> &manifold_object);
1818 
1819 
1834  DEAL_II_DEPRECATED
1835  void
1836  set_manifold(const types::manifold_id number);
1837 
1850  void
1851  reset_manifold(const types::manifold_id manifold_number);
1852 
1864  void
1865  reset_all_manifolds();
1866 
1875  void
1876  set_all_manifold_ids(const types::manifold_id number);
1877 
1886  void
1887  set_all_manifold_ids_on_boundary(const types::manifold_id number);
1888 
1898  void
1899  set_all_manifold_ids_on_boundary(const types::boundary_id b_id,
1900  const types::manifold_id number);
1901 
1911  const Manifold<dim, spacedim> &
1912  get_manifold(const types::manifold_id number) const;
1913 
1925  std::vector<types::boundary_id>
1926  get_boundary_ids() const;
1927 
1939  std::vector<types::manifold_id>
1940  get_manifold_ids() const;
1941 
1968  virtual void
1969  copy_triangulation(const Triangulation<dim, spacedim> &other_tria);
1970 
2016  virtual void
2017  create_triangulation(const std::vector<Point<spacedim>> &vertices,
2018  const std::vector<CellData<dim>> & cells,
2019  const SubCellData & subcelldata);
2020 
2029  virtual void
2030  create_triangulation_compatibility(
2031  const std::vector<Point<spacedim>> &vertices,
2032  const std::vector<CellData<dim>> & cells,
2033  const SubCellData & subcelldata);
2034 
2041  void
2042  flip_all_direction_flags();
2043 
2055  void
2056  set_all_refine_flags();
2057 
2072  void
2073  refine_global(const unsigned int times = 1);
2074 
2106  virtual void
2107  execute_coarsening_and_refinement();
2108 
2138  virtual bool
2139  prepare_coarsening_and_refinement();
2140 
2141  /*
2142  * @}
2143  */
2144 
2159  {
2176  CELL_INVALID
2177  };
2178 
2184  template <typename T>
2186  {
2187  using result_type = T;
2188 
2189  template <typename InputIterator>
2190  T
2191  operator()(InputIterator first, InputIterator last) const
2192  {
2193  return std::accumulate(first, last, T());
2194  }
2195  };
2196 
2206  struct Signals
2207  {
2215  boost::signals2::signal<void()> create;
2216 
2224  boost::signals2::signal<void()> pre_refinement;
2225 
2231  boost::signals2::signal<void()> post_refinement;
2232 
2239  boost::signals2::signal<void()> pre_partition;
2240 
2248  boost::signals2::signal<void()> mesh_movement;
2249 
2257  boost::signals2::signal<void(
2258  const typename Triangulation<dim, spacedim>::cell_iterator &cell)>
2260 
2267  boost::signals2::signal<void(
2268  const typename Triangulation<dim, spacedim>::cell_iterator &cell)>
2270 
2277  boost::signals2::signal<void(
2278  const Triangulation<dim, spacedim> &destination_tria)>
2280 
2294  boost::signals2::signal<void()> clear;
2295 
2305  boost::signals2::signal<void()> any_change;
2306 
2330  boost::signals2::signal<unsigned int(const cell_iterator &,
2331  const CellStatus),
2334 
2346  boost::signals2::signal<void()> pre_distributed_refinement;
2347 
2356  boost::signals2::signal<void()> post_distributed_refinement;
2357 
2368  boost::signals2::signal<void()> pre_distributed_repartition;
2369 
2375  boost::signals2::signal<void()> post_distributed_repartition;
2376 
2383  boost::signals2::signal<void()> pre_distributed_save;
2384 
2390  boost::signals2::signal<void()> post_distributed_save;
2391 
2398  boost::signals2::signal<void()> pre_distributed_load;
2399 
2405  boost::signals2::signal<void()> post_distributed_load;
2406  };
2407 
2411  mutable Signals signals;
2412 
2413  /*
2414  * @}
2415  */
2416 
2426  void
2427  save_refine_flags(std::ostream &out) const;
2428 
2432  void
2433  save_refine_flags(std::vector<bool> &v) const;
2434 
2438  void
2439  load_refine_flags(std::istream &in);
2440 
2444  void
2445  load_refine_flags(const std::vector<bool> &v);
2446 
2450  void
2451  save_coarsen_flags(std::ostream &out) const;
2452 
2456  void
2457  save_coarsen_flags(std::vector<bool> &v) const;
2458 
2462  void
2463  load_coarsen_flags(std::istream &out);
2464 
2468  void
2469  load_coarsen_flags(const std::vector<bool> &v);
2470 
2475  bool
2476  get_anisotropic_refinement_flag() const;
2477 
2478  /*
2479  * @}
2480  */
2481 
2491  void
2492  clear_user_flags();
2493 
2499  void
2500  save_user_flags(std::ostream &out) const;
2501 
2507  void
2508  save_user_flags(std::vector<bool> &v) const;
2509 
2514  void
2515  load_user_flags(std::istream &in);
2516 
2521  void
2522  load_user_flags(const std::vector<bool> &v);
2523 
2528  void
2529  clear_user_flags_line();
2530 
2535  void
2536  save_user_flags_line(std::ostream &out) const;
2537 
2543  void
2544  save_user_flags_line(std::vector<bool> &v) const;
2545 
2550  void
2551  load_user_flags_line(std::istream &in);
2552 
2557  void
2558  load_user_flags_line(const std::vector<bool> &v);
2559 
2564  void
2565  clear_user_flags_quad();
2566 
2571  void
2572  save_user_flags_quad(std::ostream &out) const;
2573 
2579  void
2580  save_user_flags_quad(std::vector<bool> &v) const;
2581 
2586  void
2587  load_user_flags_quad(std::istream &in);
2588 
2593  void
2594  load_user_flags_quad(const std::vector<bool> &v);
2595 
2596 
2601  void
2602  clear_user_flags_hex();
2603 
2608  void
2609  save_user_flags_hex(std::ostream &out) const;
2610 
2616  void
2617  save_user_flags_hex(std::vector<bool> &v) const;
2618 
2623  void
2624  load_user_flags_hex(std::istream &in);
2625 
2630  void
2631  load_user_flags_hex(const std::vector<bool> &v);
2632 
2638  void
2639  clear_user_data();
2640 
2646  void
2647  save_user_indices(std::vector<unsigned int> &v) const;
2648 
2653  void
2654  load_user_indices(const std::vector<unsigned int> &v);
2655 
2661  void
2662  save_user_pointers(std::vector<void *> &v) const;
2663 
2668  void
2669  load_user_pointers(const std::vector<void *> &v);
2670 
2676  void
2677  save_user_indices_line(std::vector<unsigned int> &v) const;
2678 
2683  void
2684  load_user_indices_line(const std::vector<unsigned int> &v);
2685 
2691  void
2692  save_user_indices_quad(std::vector<unsigned int> &v) const;
2693 
2698  void
2699  load_user_indices_quad(const std::vector<unsigned int> &v);
2700 
2706  void
2707  save_user_indices_hex(std::vector<unsigned int> &v) const;
2708 
2713  void
2714  load_user_indices_hex(const std::vector<unsigned int> &v);
2720  void
2721  save_user_pointers_line(std::vector<void *> &v) const;
2722 
2727  void
2728  load_user_pointers_line(const std::vector<void *> &v);
2729 
2735  void
2736  save_user_pointers_quad(std::vector<void *> &v) const;
2737 
2742  void
2743  load_user_pointers_quad(const std::vector<void *> &v);
2744 
2750  void
2751  save_user_pointers_hex(std::vector<void *> &v) const;
2752 
2757  void
2758  load_user_pointers_hex(const std::vector<void *> &v);
2759 
2760  /*
2761  * @}
2762  */
2763 
2787  begin(const unsigned int level = 0) const;
2788 
2820  begin_active(const unsigned int level = 0) const;
2821 
2827  end() const;
2828 
2848  end(const unsigned int level) const;
2849 
2870  end_active(const unsigned int level) const;
2871 
2872 
2877  last() const;
2878 
2883  last_active() const;
2884 
2900  cell_iterators() const;
2901 
2938  active_cell_iterators() const;
2939 
2956  cell_iterators_on_level(const unsigned int level) const;
2957 
2974  active_cell_iterators_on_level(const unsigned int level) const;
2975 
2976  /*
2977  * @}
2978  */
2979 
2980  /*-------------------------------------------------------------------------*/
2981 
2991  begin_face() const;
2992 
2997  begin_active_face() const;
2998 
3004  end_face() const;
3005 
3025  active_face_iterators() const;
3026 
3027  /*
3028  * @}
3029  */
3030 
3031  /*-------------------------------------------------------------------------*/
3032 
3043  begin_vertex() const;
3044 
3051  begin_active_vertex() const;
3052 
3059  end_vertex() const;
3060 
3061  /*
3062  * @}
3063  */
3064 
3082  unsigned int
3083  n_lines() const;
3084 
3088  unsigned int
3089  n_lines(const unsigned int level) const;
3090 
3094  unsigned int
3095  n_active_lines() const;
3096 
3100  unsigned int
3101  n_active_lines(const unsigned int level) const;
3102 
3106  unsigned int
3107  n_quads() const;
3108 
3112  unsigned int
3113  n_quads(const unsigned int level) const;
3114 
3118  unsigned int
3119  n_active_quads() const;
3120 
3124  unsigned int
3125  n_active_quads(const unsigned int level) const;
3126 
3130  unsigned int
3131  n_hexs() const;
3132 
3137  unsigned int
3138  n_hexs(const unsigned int level) const;
3139 
3143  unsigned int
3144  n_active_hexs() const;
3145 
3150  unsigned int
3151  n_active_hexs(const unsigned int level) const;
3152 
3157  unsigned int
3158  n_cells() const;
3159 
3164  unsigned int
3165  n_cells(const unsigned int level) const;
3166 
3171  unsigned int
3172  n_active_cells() const;
3173 
3181  virtual types::global_dof_index
3182  n_global_active_cells() const;
3183 
3184 
3189  unsigned int
3190  n_active_cells(const unsigned int level) const;
3191 
3197  unsigned int
3198  n_faces() const;
3199 
3205  unsigned int
3206  n_active_faces() const;
3207 
3225  unsigned int
3226  n_levels() const;
3227 
3234  virtual unsigned int
3235  n_global_levels() const;
3236 
3246  virtual bool
3247  has_hanging_nodes() const;
3248 
3256  unsigned int
3257  n_vertices() const;
3258 
3267  const std::vector<Point<spacedim>> &
3268  get_vertices() const;
3269 
3274  unsigned int
3275  n_used_vertices() const;
3276 
3280  bool
3281  vertex_used(const unsigned int index) const;
3282 
3287  const std::vector<bool> &
3288  get_used_vertices() const;
3289 
3301  unsigned int
3302  max_adjacent_cells() const;
3303 
3310  virtual types::subdomain_id
3311  locally_owned_subdomain() const;
3312 
3323  get_triangulation();
3324 
3330  get_triangulation() const;
3331 
3332 
3333  /*
3334  * @}
3335  */
3336 
3351  unsigned int
3352  n_raw_lines() const;
3353 
3363  unsigned int
3364  n_raw_lines(const unsigned int level) const;
3365 
3375  unsigned int
3376  n_raw_quads() const;
3377 
3387  unsigned int
3388  n_raw_quads(const unsigned int level) const;
3389 
3399  unsigned int
3400  n_raw_hexs(const unsigned int level) const;
3401 
3411  unsigned int
3412  n_raw_cells(const unsigned int level) const;
3413 
3425  unsigned int
3426  n_raw_faces() const;
3427 
3428  /*
3429  * @}
3430  */
3431 
3440  virtual std::size_t
3441  memory_consumption() const;
3442 
3451  template <class Archive>
3452  void
3453  save(Archive &ar, const unsigned int version) const;
3454 
3470  template <class Archive>
3471  void
3472  load(Archive &ar, const unsigned int version);
3473 
3474 
3490  virtual void
3491  add_periodicity(
3492  const std::vector<GridTools::PeriodicFacePair<cell_iterator>> &);
3493 
3497  const std::map<
3498  std::pair<cell_iterator, unsigned int>,
3499  std::pair<std::pair<cell_iterator, unsigned int>, std::bitset<3>>> &
3500  get_periodic_face_map() const;
3501 
3513  virtual unsigned int
3514  coarse_cell_id_to_coarse_cell_index(
3515  const types::coarse_cell_id coarse_cell_id) const;
3516 
3517 
3528  virtual types::coarse_cell_id
3529  coarse_cell_index_to_coarse_cell_id(
3530  const unsigned int coarse_cell_index) const;
3531 
3532 
3533  BOOST_SERIALIZATION_SPLIT_MEMBER()
3534 
3535 
3545  DeclException2(ExcInvalidLevel,
3546  int,
3547  int,
3548  << "You are requesting information from refinement level "
3549  << arg1
3550  << " of a triangulation, but this triangulation only has "
3551  << arg2 << " refinement levels. The given level " << arg1
3552  << " must be *less* than " << arg2 << ".");
3560  ExcTriangulationNotEmpty,
3561  int,
3562  int,
3563  << "You are trying to perform an operation on a triangulation "
3564  << "that is only allowed if the triangulation is currently empty. "
3565  << "However, it currently stores " << arg1 << " vertices and has "
3566  << "cells on " << arg2 << " levels.");
3572  DeclException0(ExcGridReadError);
3583  DeclException1(ExcEmptyLevel,
3584  int,
3585  << "You tried to do something on level " << arg1
3586  << ", but this level is empty.");
3592  DeclException0(ExcNonOrientableTriangulation);
3593 
3601  DeclException1(ExcBoundaryIdNotFound,
3602  types::boundary_id,
3603  << "The given boundary_id " << arg1
3604  << " is not defined in this Triangulation!");
3605 
3612  ExcInconsistentCoarseningFlags,
3613  "A cell is flagged for coarsening, but either not all of its siblings "
3614  "are active or flagged for coarsening as well. Please clean up all "
3615  "coarsen flags on your triangulation via "
3616  "Triangulation::prepare_coarsening_and_refinement() beforehand!");
3617 
3618  /*
3619  * @}
3620  */
3621 
3622 protected:
3627  MeshSmoothing smooth_grid;
3628 
3642  static void
3643  write_bool_vector(const unsigned int magic_number1,
3644  const std::vector<bool> &v,
3645  const unsigned int magic_number2,
3646  std::ostream & out);
3647 
3652  static void
3653  read_bool_vector(const unsigned int magic_number1,
3654  std::vector<bool> &v,
3655  const unsigned int magic_number2,
3656  std::istream & in);
3657 
3662  void
3663  update_periodic_face_map();
3664 
3665 
3666 private:
3673  std::vector<GridTools::PeriodicFacePair<cell_iterator>>
3674  periodic_face_pairs_level_0;
3675 
3680  std::map<std::pair<cell_iterator, unsigned int>,
3681  std::pair<std::pair<cell_iterator, unsigned int>, std::bitset<3>>>
3682  periodic_face_map;
3683 
3698  using raw_face_iterator =
3699  TriaRawIterator<TriaAccessor<dim - 1, dim, spacedim>>;
3700  using raw_vertex_iterator =
3701  TriaRawIterator<::TriaAccessor<0, dim, spacedim>>;
3702  using raw_line_iterator = typename IteratorSelector::raw_line_iterator;
3703  using raw_quad_iterator = typename IteratorSelector::raw_quad_iterator;
3704  using raw_hex_iterator = typename IteratorSelector::raw_hex_iterator;
3705 
3711  begin_raw(const unsigned int level = 0) const;
3712 
3718  end_raw(const unsigned int level) const;
3719 
3720  /*
3721  * @}
3722  */
3723 
3735  raw_line_iterator
3736  begin_raw_line(const unsigned int level = 0) const;
3737 
3756  begin_line(const unsigned int level = 0) const;
3757 
3776  begin_active_line(const unsigned int level = 0) const;
3777 
3783  end_line() const;
3784 
3785  /*
3786  * @}
3787  */
3788 
3813  raw_quad_iterator
3814  begin_raw_quad(const unsigned int level = 0) const;
3815 
3834  begin_quad(const unsigned int level = 0) const;
3835 
3854  begin_active_quad(const unsigned int level = 0) const;
3855 
3861  end_quad() const;
3862 
3863  /*
3864  * @}
3865  */
3866 
3890  raw_hex_iterator
3891  begin_raw_hex(const unsigned int level = 0) const;
3892 
3910  hex_iterator
3911  begin_hex(const unsigned int level = 0) const;
3912 
3931  begin_active_hex(const unsigned int level = 0) const;
3932 
3937  hex_iterator
3938  end_hex() const;
3939 
3940  /*
3941  * @}
3942  */
3943 
3944 
3958  void
3959  clear_despite_subscriptions();
3960 
3967  void
3968  reset_active_cell_indices();
3969 
3983  DistortedCellList
3984  execute_refinement();
3985 
3992  void
3993  execute_coarsening();
3994 
3999  void
4000  fix_coarsen_flags();
4001 
4006  std::vector<std::unique_ptr<
4007  ::internal::TriangulationImplementation::TriaLevel<dim>>>
4008  levels;
4009 
4015  std::unique_ptr<::internal::TriangulationImplementation::TriaFaces<dim>>
4016  faces;
4017 
4018 
4022  std::vector<Point<spacedim>> vertices;
4023 
4027  std::vector<bool> vertices_used;
4028 
4033  std::map<types::manifold_id, std::unique_ptr<const Manifold<dim, spacedim>>>
4034  manifold;
4035 
4039  bool anisotropic_refinement;
4040 
4041 
4046  const bool check_for_distorted_cells;
4047 
4057  ::internal::TriangulationImplementation::NumberCache<dim> number_cache;
4058 
4073  std::unique_ptr<std::map<unsigned int, types::boundary_id>>
4074  vertex_to_boundary_id_map_1d;
4075 
4076 
4096  std::unique_ptr<std::map<unsigned int, types::manifold_id>>
4097  vertex_to_manifold_id_map_1d;
4098 
4099  // make a couple of classes friends
4100  template <int, int, int>
4101  friend class TriaAccessorBase;
4102  template <int, int, int>
4103  friend class TriaAccessor;
4104  friend class TriaAccessor<0, 1, spacedim>;
4105 
4106  friend class CellAccessor<dim, spacedim>;
4107 
4108  friend struct ::internal::TriaAccessorImplementation::Implementation;
4109 
4110  friend class hp::DoFHandler<dim, spacedim>;
4111 
4112  friend struct ::internal::TriangulationImplementation::Implementation;
4113 
4114  template <typename>
4115  friend class ::internal::TriangulationImplementation::TriaObjects;
4116 
4117  // explicitly check for sensible template arguments, but not on windows
4118  // because MSVC creates bogus warnings during normal compilation
4119 #ifndef DEAL_II_MSVC
4120  static_assert(dim <= spacedim,
4121  "The dimension <dim> of a Triangulation must be less than or "
4122  "equal to the space dimension <spacedim> in which it lives.");
4123 #endif
4124 };
4125 
4126 
4127 #ifndef DOXYGEN
4128 
4129 
4130 template <int structdim>
4131 template <class Archive>
4132 void
4133 CellData<structdim>::serialize(Archive &ar, const unsigned int /*version*/)
4134 {
4135  ar &vertices;
4136  ar &material_id;
4137  ar &boundary_id;
4138  ar &manifold_id;
4139 }
4140 
4141 
4142 
4143 namespace internal
4144 {
4145  namespace TriangulationImplementation
4146  {
4147  template <class Archive>
4148  void
4149  NumberCache<1>::serialize(Archive &ar, const unsigned int)
4150  {
4151  ar &n_levels;
4152  ar &n_lines &n_lines_level;
4153  ar &n_active_lines &n_active_lines_level;
4154  }
4155 
4156 
4157  template <class Archive>
4158  void
4159  NumberCache<2>::serialize(Archive &ar, const unsigned int version)
4160  {
4161  this->NumberCache<1>::serialize(ar, version);
4162 
4163  ar &n_quads &n_quads_level;
4164  ar &n_active_quads &n_active_quads_level;
4165  }
4166 
4167 
4168  template <class Archive>
4169  void
4170  NumberCache<3>::serialize(Archive &ar, const unsigned int version)
4171  {
4172  this->NumberCache<2>::serialize(ar, version);
4173 
4174  ar &n_hexes &n_hexes_level;
4175  ar &n_active_hexes &n_active_hexes_level;
4176  }
4177 
4178  } // namespace TriangulationImplementation
4179 } // namespace internal
4180 
4181 
4182 template <int dim, int spacedim>
4183 inline bool
4184 Triangulation<dim, spacedim>::vertex_used(const unsigned int index) const
4185 {
4186  Assert(index < vertices_used.size(),
4187  ExcIndexRange(index, 0, vertices_used.size()));
4188  return vertices_used[index];
4189 }
4190 
4191 
4192 
4193 template <int dim, int spacedim>
4194 inline unsigned int
4196 {
4197  return number_cache.n_levels;
4198 }
4199 
4200 template <int dim, int spacedim>
4201 inline unsigned int
4203 {
4204  return number_cache.n_levels;
4205 }
4206 
4207 
4208 template <int dim, int spacedim>
4209 inline unsigned int
4211 {
4212  return vertices.size();
4213 }
4214 
4215 
4216 
4217 template <int dim, int spacedim>
4218 inline const std::vector<Point<spacedim>> &
4220 {
4221  return vertices;
4222 }
4223 
4224 
4225 template <int dim, int spacedim>
4226 template <class Archive>
4227 void
4228 Triangulation<dim, spacedim>::save(Archive &ar, const unsigned int) const
4229 {
4230  // as discussed in the documentation, do not store the signals as
4231  // well as boundary and manifold description but everything else
4232  ar &smooth_grid;
4233 
4234  unsigned int n_levels = levels.size();
4235  ar & n_levels;
4236  for (unsigned int i = 0; i < levels.size(); ++i)
4237  ar &levels[i];
4238 
4239  // boost dereferences a nullptr when serializing a nullptr
4240  // at least up to 1.65.1. This causes problems with clang-5.
4241  // Therefore, work around it.
4242  bool faces_is_nullptr = (faces.get() == nullptr);
4243  ar & faces_is_nullptr;
4244  if (!faces_is_nullptr)
4245  ar &faces;
4246 
4247  ar &vertices;
4248  ar &vertices_used;
4249 
4250  ar &anisotropic_refinement;
4251  ar &number_cache;
4252 
4253  ar &check_for_distorted_cells;
4254 
4255  if (dim == 1)
4256  {
4257  ar &vertex_to_boundary_id_map_1d;
4258  ar &vertex_to_manifold_id_map_1d;
4259  }
4260 }
4261 
4262 
4263 
4264 template <int dim, int spacedim>
4265 template <class Archive>
4266 void
4267 Triangulation<dim, spacedim>::load(Archive &ar, const unsigned int)
4268 {
4269  // clear previous content. this also calls the respective signal
4270  clear();
4271 
4272  // as discussed in the documentation, do not store the signals as
4273  // well as boundary and manifold description but everything else
4274  ar &smooth_grid;
4275 
4276  unsigned int size;
4277  ar & size;
4278  levels.resize(size);
4279  for (unsigned int i = 0; i < levels.size(); ++i)
4280  {
4281  std::unique_ptr<internal::TriangulationImplementation::TriaLevel<dim>>
4282  level;
4283  ar &level;
4284  levels[i] = std::move(level);
4285  }
4286 
4287  // Workaround for nullptr, see in save().
4288  bool faces_is_nullptr = true;
4289  ar & faces_is_nullptr;
4290  if (!faces_is_nullptr)
4291  ar &faces;
4292 
4293  ar &vertices;
4294  ar &vertices_used;
4295 
4296  ar &anisotropic_refinement;
4297  ar &number_cache;
4298 
4299  // the levels do not serialize the active_cell_indices because
4300  // they are easy enough to rebuild upon re-loading data. do
4301  // this here. don't forget to first resize the fields appropriately
4302  {
4303  for (unsigned int l = 0; l < levels.size(); ++l)
4304  levels[l]->active_cell_indices.resize(levels[l]->refine_flags.size());
4305  reset_active_cell_indices();
4306  }
4307 
4308 
4309  bool my_check_for_distorted_cells;
4310  ar & my_check_for_distorted_cells;
4311 
4312  Assert(my_check_for_distorted_cells == check_for_distorted_cells,
4313  ExcMessage("The triangulation loaded into here must have the "
4314  "same setting with regard to reporting distorted "
4315  "cell as the one previously stored."));
4316 
4317  if (dim == 1)
4318  {
4319  ar &vertex_to_boundary_id_map_1d;
4320  ar &vertex_to_manifold_id_map_1d;
4321  }
4322 
4323  // trigger the create signal to indicate
4324  // that new content has been imported into
4325  // the triangulation
4326  signals.create();
4327 }
4328 
4329 
4330 
4331 template <int dim, int spacedim>
4332 inline unsigned int
4334  const types::coarse_cell_id coarse_cell_id) const
4335 {
4336  return coarse_cell_id;
4337 }
4338 
4339 
4340 
4341 template <int dim, int spacedim>
4342 inline types::coarse_cell_id
4344  const unsigned int coarse_cell_index) const
4345 {
4346  return coarse_cell_index;
4347 }
4348 
4349 
4350 
4351 /* -------------- declaration of explicit specializations ------------- */
4352 
4353 template <>
4354 unsigned int
4355 Triangulation<1, 1>::n_raw_lines(const unsigned int level) const;
4356 template <>
4357 unsigned int
4359 template <>
4360 unsigned int
4361 Triangulation<1, 1>::n_quads(const unsigned int level) const;
4362 template <>
4363 unsigned int
4364 Triangulation<1, 1>::n_raw_quads(const unsigned int level) const;
4365 template <>
4366 unsigned int
4367 Triangulation<2, 2>::n_raw_quads(const unsigned int level) const;
4368 template <>
4369 unsigned int
4370 Triangulation<1, 1>::n_raw_hexs(const unsigned int level) const;
4371 template <>
4372 unsigned int
4373 Triangulation<1, 1>::n_active_quads(const unsigned int level) const;
4374 template <>
4375 unsigned int
4377 template <>
4378 unsigned int
4380 
4381 
4382 // -------------------------------------------------------------------
4383 // -- Explicit specializations for codimension one grids
4384 
4385 
4386 template <>
4387 unsigned int
4388 Triangulation<1, 2>::n_raw_lines(const unsigned int level) const;
4389 template <>
4390 unsigned int
4392 template <>
4393 unsigned int
4394 Triangulation<1, 2>::n_quads(const unsigned int level) const;
4395 template <>
4396 unsigned int
4397 Triangulation<1, 2>::n_raw_quads(const unsigned int level) const;
4398 template <>
4399 unsigned int
4400 Triangulation<2, 3>::n_raw_quads(const unsigned int level) const;
4401 template <>
4402 unsigned int
4403 Triangulation<1, 2>::n_raw_hexs(const unsigned int level) const;
4404 template <>
4405 unsigned int
4406 Triangulation<1, 2>::n_active_quads(const unsigned int level) const;
4407 template <>
4408 unsigned int
4410 template <>
4411 unsigned int
4413 
4414 // -------------------------------------------------------------------
4415 // -- Explicit specializations for codimension two grids
4416 
4417 
4418 template <>
4419 unsigned int
4420 Triangulation<1, 3>::n_raw_lines(const unsigned int level) const;
4421 template <>
4422 unsigned int
4424 template <>
4425 unsigned int
4426 Triangulation<1, 3>::n_quads(const unsigned int level) const;
4427 template <>
4428 unsigned int
4429 Triangulation<1, 3>::n_raw_quads(const unsigned int level) const;
4430 template <>
4431 unsigned int
4432 Triangulation<2, 3>::n_raw_quads(const unsigned int level) const;
4433 template <>
4434 unsigned int
4435 Triangulation<1, 3>::n_raw_hexs(const unsigned int level) const;
4436 template <>
4437 unsigned int
4438 Triangulation<1, 3>::n_active_quads(const unsigned int level) const;
4439 template <>
4440 unsigned int
4442 template <>
4443 unsigned int
4445 
4446 
4447 #endif // DOXYGEN
4448 
4449 DEAL_II_NAMESPACE_CLOSE
4450 
4451 // Include tria_accessor.h here, so that it is possible for an end
4452 // user to use the iterators of Triangulation<dim> directly without
4453 // the need to include tria_accessor.h separately. (Otherwise the
4454 // iterators are an 'opaque' or 'incomplete' type.)
4455 #include <deal.II/grid/tria_accessor.h>
4456 
4457 #endif
std::vector< CellData< 1 > > boundary_lines
Definition: tria.h:276
boost::signals2::signal< void(const Triangulation< dim, spacedim > &destination_tria)> copy
Definition: tria.h:2279
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:2231
unsigned int n_vertices() const
unsigned int manifold_id
Definition: types.h:137
boost::signals2::signal< void(const typename Triangulation< dim, spacedim >::cell_iterator &cell)> pre_coarsening_on_cell
Definition: tria.h:2259
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:541
Definition: tria.h:136
boost::signals2::signal< void()> post_distributed_repartition
Definition: tria.h:2375
typename IteratorSelector::line_iterator line_iterator
Definition: tria.h:1610
typename IteratorSelector::hex_iterator hex_iterator
Definition: tria.h:1658
typename IteratorSelector::active_line_iterator active_line_iterator
Definition: tria.h:1625
boost::signals2::signal< void()> pre_distributed_refinement
Definition: tria.h:2346
unsigned int material_id
Definition: types.h:148
STL namespace.
types::boundary_id boundary_id
Definition: tria.h:173
unsigned int coarse_cell_id
Definition: types.h:109
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
static ::ExceptionBase & ExcFacesHaveNoLevel()
unsigned int n_levels() const
static ::ExceptionBase & ExcMessage(std::string arg1)
boost::signals2::signal< void()> pre_refinement
Definition: tria.h:2224
unsigned int subdomain_id
Definition: types.h:43
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
Definition: types.h:31
unsigned int n_raw_quads() const
Definition: tria.cc:13052
#define Assert(cond, exc)
Definition: exceptions.h:1407
Signals signals
Definition: tria.h:2411
boost::signals2::signal< void()> clear
Definition: tria.h:2294
::internal::TriangulationImplementation::Iterators< dim, spacedim > IteratorSelector
Definition: tria.h:1303
unsigned int n_raw_hexs(const unsigned int level) const
Definition: tria.cc:13109
unsigned int max_adjacent_cells() const
Definition: tria.cc:13228
void load(Archive &ar, const unsigned int version)
std::list< typename Triangulation< dim, spacedim >::cell_iterator > distorted_cells
Definition: tria.h:1704
unsigned int n_quads() const
Definition: tria.cc:13003
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
#define DeclException0(Exception0)
Definition: exceptions.h:473
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:2269
typename IteratorSelector::active_quad_iterator active_quad_iterator
Definition: tria.h:1649
types::material_id material_id
Definition: tria.h:162
const std::vector< Point< spacedim > > & get_vertices() const
boost::signals2::signal< void()> mesh_movement
Definition: tria.h:2248
unsigned int n_raw_lines() const
Definition: tria.cc:12818
boost::signals2::signal< void()> pre_distributed_save
Definition: tria.h:2383
void save(Archive &ar, const unsigned int version) const
typename IteratorSelector::quad_iterator quad_iterator
Definition: tria.h:1634
boost::signals2::signal< void()> create
Definition: tria.h:2215
Definition: hp.h:117
unsigned int n_active_quads() const
Definition: tria.cc:13071
types::manifold_id manifold_id
Definition: tria.h:184
boost::signals2::signal< void()> pre_partition
Definition: tria.h:2239
typename IteratorSelector::active_hex_iterator active_hex_iterator
Definition: tria.h:1669
unsigned int global_dof_index
Definition: types.h:89
boost::signals2::signal< void()> post_distributed_save
Definition: tria.h:2390
std::vector< unsigned int > n_active_hexes_level
Definition: tria.h:485
virtual unsigned int n_global_levels() const
std::vector< CellData< 2 > > boundary_quads
Definition: tria.h:284
boost::signals2::signal< unsigned int(const cell_iterator &, const CellStatus), CellWeightSum< unsigned int > > cell_weight
Definition: tria.h:2333
boost::signals2::signal< void()> pre_distributed_load
Definition: tria.h:2398
std::vector< unsigned int > n_active_lines_level
Definition: tria.h:366
std::vector< unsigned int > n_active_quads_level
Definition: tria.h:425
boost::signals2::signal< void()> post_distributed_load
Definition: tria.h:2405
unsigned int boundary_id
Definition: types.h:125
boost::signals2::signal< void()> pre_distributed_repartition
Definition: tria.h:2368
bool vertex_used(const unsigned int index) const
boost::signals2::signal< void()> any_change
Definition: tria.h:2305
void serialize(Archive &ar, const unsigned int version)
boost::signals2::signal< void()> post_distributed_refinement
Definition: tria.h:2356