Reference documentation for deal.II version Git a5ed68a04a 2019-09-22 06:50:58 -0600
\(\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 
209 
210 
266 {
273  std::vector<CellData<1>> boundary_lines;
274 
281  std::vector<CellData<2>> boundary_quads;
282 
291  bool
292  check_consistency(const unsigned int dim) const;
293 };
294 
295 
296 /*------------------------------------------------------------------------*/
297 
298 
299 namespace internal
300 {
305  namespace TriangulationImplementation
306  {
320  template <int dim>
321  struct NumberCache
322  {};
323 
337  template <>
338  struct NumberCache<1>
339  {
343  unsigned int n_levels;
344 
348  unsigned int n_lines;
349 
353  std::vector<unsigned int> n_lines_level;
354 
358  unsigned int n_active_lines;
359 
363  std::vector<unsigned int> n_active_lines_level;
364 
368  NumberCache();
369 
374  std::size_t
375  memory_consumption() const;
376 
381  template <class Archive>
382  void
383  serialize(Archive &ar, const unsigned int version);
384  };
385 
386 
401  template <>
402  struct NumberCache<2> : public NumberCache<1>
403  {
407  unsigned int n_quads;
408 
412  std::vector<unsigned int> n_quads_level;
413 
417  unsigned int n_active_quads;
418 
422  std::vector<unsigned int> n_active_quads_level;
423 
427  NumberCache();
428 
433  std::size_t
434  memory_consumption() const;
435 
440  template <class Archive>
441  void
442  serialize(Archive &ar, const unsigned int version);
443  };
444 
445 
461  template <>
462  struct NumberCache<3> : public NumberCache<2>
463  {
467  unsigned int n_hexes;
468 
472  std::vector<unsigned int> n_hexes_level;
473 
477  unsigned int n_active_hexes;
478 
482  std::vector<unsigned int> n_active_hexes_level;
483 
487  NumberCache();
488 
493  std::size_t
494  memory_consumption() const;
495 
500  template <class Archive>
501  void
502  serialize(Archive &ar, const unsigned int version);
503  };
504  } // namespace TriangulationImplementation
505 } // namespace internal
506 
507 
508 /*------------------------------------------------------------------------*/
509 
510 
1291 template <int dim, int spacedim = dim>
1293 {
1294 private:
1299  using IteratorSelector =
1300  ::internal::TriangulationImplementation::Iterators<dim, spacedim>;
1301 
1302 public:
1308  {
1313  none = 0x0,
1357  limit_level_difference_at_vertices = 0x1,
1378  eliminate_unrefined_islands = 0x2,
1394  patch_level_1 = 0x4,
1418  coarsest_level_1 = 0x8,
1443  allow_anisotropic_smoothing = 0x10,
1476  eliminate_refined_inner_islands = 0x100,
1481  eliminate_refined_boundary_islands = 0x200,
1487  do_not_produce_unrefined_islands = 0x400,
1488 
1493  smoothing_on_refinement =
1494  (limit_level_difference_at_vertices | eliminate_unrefined_islands),
1499  smoothing_on_coarsening =
1500  (eliminate_refined_inner_islands | eliminate_refined_boundary_islands |
1501  do_not_produce_unrefined_islands),
1502 
1508  maximum_smoothing = 0xffff ^ allow_anisotropic_smoothing
1509  };
1510 
1527 
1545 
1559  using face_iterator = TriaIterator<TriaAccessor<dim - 1, dim, spacedim>>;
1560 
1572  using active_face_iterator =
1573  TriaActiveIterator<TriaAccessor<dim - 1, dim, spacedim>>;
1574 
1584 
1597  using active_vertex_iterator =
1599 
1607  using line_iterator = typename IteratorSelector::line_iterator;
1608 
1622  using active_line_iterator = typename IteratorSelector::active_line_iterator;
1623 
1631  using quad_iterator = typename IteratorSelector::quad_iterator;
1632 
1646  using active_quad_iterator = typename IteratorSelector::active_quad_iterator;
1647 
1655  using hex_iterator = typename IteratorSelector::hex_iterator;
1656 
1666  using active_hex_iterator = typename IteratorSelector::active_hex_iterator;
1667 
1687  {
1694  virtual ~DistortedCellList() noexcept override = default;
1695 
1700  std::list<typename Triangulation<dim, spacedim>::cell_iterator>
1702  };
1703 
1707  static const unsigned int dimension = dim;
1708 
1712  static const unsigned int space_dimension = spacedim;
1713 
1728  Triangulation(const MeshSmoothing smooth_grid = none,
1729  const bool check_for_distorted_cells = false);
1730 
1746  Triangulation(const Triangulation<dim, spacedim> &) = delete;
1747 
1755 
1759  Triangulation &
1760  operator=(Triangulation<dim, spacedim> &&tria) noexcept;
1761 
1765  virtual ~Triangulation() override;
1766 
1773  virtual void
1774  clear();
1775 
1781  virtual void
1782  set_mesh_smoothing(const MeshSmoothing mesh_smoothing);
1783 
1787  virtual const MeshSmoothing &
1788  get_mesh_smoothing() const;
1789 
1812  void
1813  set_manifold(const types::manifold_id number,
1814  const Manifold<dim, spacedim> &manifold_object);
1815 
1816 
1831  DEAL_II_DEPRECATED
1832  void
1833  set_manifold(const types::manifold_id number);
1834 
1847  void
1848  reset_manifold(const types::manifold_id manifold_number);
1849 
1861  void
1862  reset_all_manifolds();
1863 
1872  void
1873  set_all_manifold_ids(const types::manifold_id number);
1874 
1883  void
1884  set_all_manifold_ids_on_boundary(const types::manifold_id number);
1885 
1895  void
1896  set_all_manifold_ids_on_boundary(const types::boundary_id b_id,
1897  const types::manifold_id number);
1898 
1908  const Manifold<dim, spacedim> &
1909  get_manifold(const types::manifold_id number) const;
1910 
1922  std::vector<types::boundary_id>
1923  get_boundary_ids() const;
1924 
1936  std::vector<types::manifold_id>
1937  get_manifold_ids() const;
1938 
1965  virtual void
1966  copy_triangulation(const Triangulation<dim, spacedim> &other_tria);
1967 
2013  virtual void
2014  create_triangulation(const std::vector<Point<spacedim>> &vertices,
2015  const std::vector<CellData<dim>> & cells,
2016  const SubCellData & subcelldata);
2017 
2026  virtual void
2027  create_triangulation_compatibility(
2028  const std::vector<Point<spacedim>> &vertices,
2029  const std::vector<CellData<dim>> & cells,
2030  const SubCellData & subcelldata);
2031 
2038  void
2039  flip_all_direction_flags();
2040 
2052  void
2053  set_all_refine_flags();
2054 
2069  void
2070  refine_global(const unsigned int times = 1);
2071 
2103  virtual void
2104  execute_coarsening_and_refinement();
2105 
2135  virtual bool
2136  prepare_coarsening_and_refinement();
2137 
2138  /*
2139  * @}
2140  */
2141 
2156  {
2173  CELL_INVALID
2174  };
2175 
2181  template <typename T>
2183  {
2184  using result_type = T;
2185 
2186  template <typename InputIterator>
2187  T
2188  operator()(InputIterator first, InputIterator last) const
2189  {
2190  return std::accumulate(first, last, T());
2191  }
2192  };
2193 
2203  struct Signals
2204  {
2212  boost::signals2::signal<void()> create;
2213 
2221  boost::signals2::signal<void()> pre_refinement;
2222 
2228  boost::signals2::signal<void()> post_refinement;
2229 
2236  boost::signals2::signal<void()> pre_partition;
2237 
2245  boost::signals2::signal<void()> mesh_movement;
2246 
2254  boost::signals2::signal<void(
2255  const typename Triangulation<dim, spacedim>::cell_iterator &cell)>
2257 
2264  boost::signals2::signal<void(
2265  const typename Triangulation<dim, spacedim>::cell_iterator &cell)>
2267 
2274  boost::signals2::signal<void(
2275  const Triangulation<dim, spacedim> &destination_tria)>
2277 
2291  boost::signals2::signal<void()> clear;
2292 
2302  boost::signals2::signal<void()> any_change;
2303 
2327  boost::signals2::signal<unsigned int(const cell_iterator &,
2328  const CellStatus),
2331 
2343  boost::signals2::signal<void()> pre_distributed_refinement;
2344 
2353  boost::signals2::signal<void()> post_distributed_refinement;
2354 
2365  boost::signals2::signal<void()> pre_distributed_repartition;
2366 
2372  boost::signals2::signal<void()> post_distributed_repartition;
2373 
2380  boost::signals2::signal<void()> pre_distributed_save;
2381 
2387  boost::signals2::signal<void()> post_distributed_save;
2388 
2395  boost::signals2::signal<void()> pre_distributed_load;
2396 
2402  boost::signals2::signal<void()> post_distributed_load;
2403  };
2404 
2408  mutable Signals signals;
2409 
2410  /*
2411  * @}
2412  */
2413 
2423  void
2424  save_refine_flags(std::ostream &out) const;
2425 
2429  void
2430  save_refine_flags(std::vector<bool> &v) const;
2431 
2435  void
2436  load_refine_flags(std::istream &in);
2437 
2441  void
2442  load_refine_flags(const std::vector<bool> &v);
2443 
2447  void
2448  save_coarsen_flags(std::ostream &out) const;
2449 
2453  void
2454  save_coarsen_flags(std::vector<bool> &v) const;
2455 
2459  void
2460  load_coarsen_flags(std::istream &out);
2461 
2465  void
2466  load_coarsen_flags(const std::vector<bool> &v);
2467 
2472  bool
2473  get_anisotropic_refinement_flag() const;
2474 
2475  /*
2476  * @}
2477  */
2478 
2488  void
2489  clear_user_flags();
2490 
2496  void
2497  save_user_flags(std::ostream &out) const;
2498 
2504  void
2505  save_user_flags(std::vector<bool> &v) const;
2506 
2511  void
2512  load_user_flags(std::istream &in);
2513 
2518  void
2519  load_user_flags(const std::vector<bool> &v);
2520 
2525  void
2526  clear_user_flags_line();
2527 
2532  void
2533  save_user_flags_line(std::ostream &out) const;
2534 
2540  void
2541  save_user_flags_line(std::vector<bool> &v) const;
2542 
2547  void
2548  load_user_flags_line(std::istream &in);
2549 
2554  void
2555  load_user_flags_line(const std::vector<bool> &v);
2556 
2561  void
2562  clear_user_flags_quad();
2563 
2568  void
2569  save_user_flags_quad(std::ostream &out) const;
2570 
2576  void
2577  save_user_flags_quad(std::vector<bool> &v) const;
2578 
2583  void
2584  load_user_flags_quad(std::istream &in);
2585 
2590  void
2591  load_user_flags_quad(const std::vector<bool> &v);
2592 
2593 
2598  void
2599  clear_user_flags_hex();
2600 
2605  void
2606  save_user_flags_hex(std::ostream &out) const;
2607 
2613  void
2614  save_user_flags_hex(std::vector<bool> &v) const;
2615 
2620  void
2621  load_user_flags_hex(std::istream &in);
2622 
2627  void
2628  load_user_flags_hex(const std::vector<bool> &v);
2629 
2635  void
2636  clear_user_data();
2637 
2643  void
2644  save_user_indices(std::vector<unsigned int> &v) const;
2645 
2650  void
2651  load_user_indices(const std::vector<unsigned int> &v);
2652 
2658  void
2659  save_user_pointers(std::vector<void *> &v) const;
2660 
2665  void
2666  load_user_pointers(const std::vector<void *> &v);
2667 
2673  void
2674  save_user_indices_line(std::vector<unsigned int> &v) const;
2675 
2680  void
2681  load_user_indices_line(const std::vector<unsigned int> &v);
2682 
2688  void
2689  save_user_indices_quad(std::vector<unsigned int> &v) const;
2690 
2695  void
2696  load_user_indices_quad(const std::vector<unsigned int> &v);
2697 
2703  void
2704  save_user_indices_hex(std::vector<unsigned int> &v) const;
2705 
2710  void
2711  load_user_indices_hex(const std::vector<unsigned int> &v);
2717  void
2718  save_user_pointers_line(std::vector<void *> &v) const;
2719 
2724  void
2725  load_user_pointers_line(const std::vector<void *> &v);
2726 
2732  void
2733  save_user_pointers_quad(std::vector<void *> &v) const;
2734 
2739  void
2740  load_user_pointers_quad(const std::vector<void *> &v);
2741 
2747  void
2748  save_user_pointers_hex(std::vector<void *> &v) const;
2749 
2754  void
2755  load_user_pointers_hex(const std::vector<void *> &v);
2756 
2757  /*
2758  * @}
2759  */
2760 
2784  begin(const unsigned int level = 0) const;
2785 
2817  begin_active(const unsigned int level = 0) const;
2818 
2824  end() const;
2825 
2845  end(const unsigned int level) const;
2846 
2867  end_active(const unsigned int level) const;
2868 
2869 
2874  last() const;
2875 
2880  last_active() const;
2881 
2897  cell_iterators() const;
2898 
2935  active_cell_iterators() const;
2936 
2953  cell_iterators_on_level(const unsigned int level) const;
2954 
2971  active_cell_iterators_on_level(const unsigned int level) const;
2972 
2973  /*
2974  * @}
2975  */
2976 
2977  /*-------------------------------------------------------------------------*/
2978 
2988  begin_face() const;
2989 
2994  begin_active_face() const;
2995 
3001  end_face() const;
3002 
3022  active_face_iterators() const;
3023 
3024  /*
3025  * @}
3026  */
3027 
3028  /*-------------------------------------------------------------------------*/
3029 
3040  begin_vertex() const;
3041 
3048  begin_active_vertex() const;
3049 
3056  end_vertex() const;
3057 
3058  /*
3059  * @}
3060  */
3061 
3079  unsigned int
3080  n_lines() const;
3081 
3085  unsigned int
3086  n_lines(const unsigned int level) const;
3087 
3091  unsigned int
3092  n_active_lines() const;
3093 
3097  unsigned int
3098  n_active_lines(const unsigned int level) const;
3099 
3103  unsigned int
3104  n_quads() const;
3105 
3109  unsigned int
3110  n_quads(const unsigned int level) const;
3111 
3115  unsigned int
3116  n_active_quads() const;
3117 
3121  unsigned int
3122  n_active_quads(const unsigned int level) const;
3123 
3127  unsigned int
3128  n_hexs() const;
3129 
3134  unsigned int
3135  n_hexs(const unsigned int level) const;
3136 
3140  unsigned int
3141  n_active_hexs() const;
3142 
3147  unsigned int
3148  n_active_hexs(const unsigned int level) const;
3149 
3154  unsigned int
3155  n_cells() const;
3156 
3161  unsigned int
3162  n_cells(const unsigned int level) const;
3163 
3168  unsigned int
3169  n_active_cells() const;
3170 
3178  virtual types::global_dof_index
3179  n_global_active_cells() const;
3180 
3181 
3186  unsigned int
3187  n_active_cells(const unsigned int level) const;
3188 
3194  unsigned int
3195  n_faces() const;
3196 
3202  unsigned int
3203  n_active_faces() const;
3204 
3222  unsigned int
3223  n_levels() const;
3224 
3231  virtual unsigned int
3232  n_global_levels() const;
3233 
3243  virtual bool
3244  has_hanging_nodes() const;
3245 
3253  unsigned int
3254  n_vertices() const;
3255 
3264  const std::vector<Point<spacedim>> &
3265  get_vertices() const;
3266 
3271  unsigned int
3272  n_used_vertices() const;
3273 
3277  bool
3278  vertex_used(const unsigned int index) const;
3279 
3284  const std::vector<bool> &
3285  get_used_vertices() const;
3286 
3298  unsigned int
3299  max_adjacent_cells() const;
3300 
3307  virtual types::subdomain_id
3308  locally_owned_subdomain() const;
3309 
3320  get_triangulation();
3321 
3327  get_triangulation() const;
3328 
3329 
3330  /*
3331  * @}
3332  */
3333 
3348  unsigned int
3349  n_raw_lines() const;
3350 
3360  unsigned int
3361  n_raw_lines(const unsigned int level) const;
3362 
3372  unsigned int
3373  n_raw_quads() const;
3374 
3384  unsigned int
3385  n_raw_quads(const unsigned int level) const;
3386 
3396  unsigned int
3397  n_raw_hexs(const unsigned int level) const;
3398 
3408  unsigned int
3409  n_raw_cells(const unsigned int level) const;
3410 
3422  unsigned int
3423  n_raw_faces() const;
3424 
3425  /*
3426  * @}
3427  */
3428 
3437  virtual std::size_t
3438  memory_consumption() const;
3439 
3448  template <class Archive>
3449  void
3450  save(Archive &ar, const unsigned int version) const;
3451 
3467  template <class Archive>
3468  void
3469  load(Archive &ar, const unsigned int version);
3470 
3471 
3487  virtual void
3488  add_periodicity(
3489  const std::vector<GridTools::PeriodicFacePair<cell_iterator>> &);
3490 
3494  const std::map<
3495  std::pair<cell_iterator, unsigned int>,
3496  std::pair<std::pair<cell_iterator, unsigned int>, std::bitset<3>>> &
3497  get_periodic_face_map() const;
3498 
3510  virtual unsigned int
3511  coarse_cell_id_to_coarse_cell_index(
3512  const types::coarse_cell_id coarse_cell_id) const;
3513 
3514 
3525  virtual types::coarse_cell_id
3526  coarse_cell_index_to_coarse_cell_id(
3527  const unsigned int coarse_cell_index) const;
3528 
3529 
3530  BOOST_SERIALIZATION_SPLIT_MEMBER()
3531 
3532 
3542  DeclException2(ExcInvalidLevel,
3543  int,
3544  int,
3545  << "You are requesting information from refinement level "
3546  << arg1
3547  << " of a triangulation, but this triangulation only has "
3548  << arg2 << " refinement levels. The given level " << arg1
3549  << " must be *less* than " << arg2 << ".");
3557  ExcTriangulationNotEmpty,
3558  int,
3559  int,
3560  << "You are trying to perform an operation on a triangulation "
3561  << "that is only allowed if the triangulation is currently empty. "
3562  << "However, it currently stores " << arg1 << " vertices and has "
3563  << "cells on " << arg2 << " levels.");
3569  DeclException0(ExcGridReadError);
3580  DeclException1(ExcEmptyLevel,
3581  int,
3582  << "You tried to do something on level " << arg1
3583  << ", but this level is empty.");
3589  DeclException0(ExcNonOrientableTriangulation);
3590 
3598  DeclException1(ExcBoundaryIdNotFound,
3599  types::boundary_id,
3600  << "The given boundary_id " << arg1
3601  << " is not defined in this Triangulation!");
3602 
3609  ExcInconsistentCoarseningFlags,
3610  "A cell is flagged for coarsening, but either not all of its siblings "
3611  "are active or flagged for coarsening as well. Please clean up all "
3612  "coarsen flags on your triangulation via "
3613  "Triangulation::prepare_coarsening_and_refinement() beforehand!");
3614 
3615  /*
3616  * @}
3617  */
3618 
3619 protected:
3624  MeshSmoothing smooth_grid;
3625 
3639  static void
3640  write_bool_vector(const unsigned int magic_number1,
3641  const std::vector<bool> &v,
3642  const unsigned int magic_number2,
3643  std::ostream & out);
3644 
3649  static void
3650  read_bool_vector(const unsigned int magic_number1,
3651  std::vector<bool> &v,
3652  const unsigned int magic_number2,
3653  std::istream & in);
3654 
3659  void
3660  update_periodic_face_map();
3661 
3662 
3663 private:
3670  std::vector<GridTools::PeriodicFacePair<cell_iterator>>
3671  periodic_face_pairs_level_0;
3672 
3677  std::map<std::pair<cell_iterator, unsigned int>,
3678  std::pair<std::pair<cell_iterator, unsigned int>, std::bitset<3>>>
3679  periodic_face_map;
3680 
3695  using raw_face_iterator =
3696  TriaRawIterator<TriaAccessor<dim - 1, dim, spacedim>>;
3697  using raw_vertex_iterator =
3698  TriaRawIterator<::TriaAccessor<0, dim, spacedim>>;
3699  using raw_line_iterator = typename IteratorSelector::raw_line_iterator;
3700  using raw_quad_iterator = typename IteratorSelector::raw_quad_iterator;
3701  using raw_hex_iterator = typename IteratorSelector::raw_hex_iterator;
3702 
3708  begin_raw(const unsigned int level = 0) const;
3709 
3715  end_raw(const unsigned int level) const;
3716 
3717  /*
3718  * @}
3719  */
3720 
3732  raw_line_iterator
3733  begin_raw_line(const unsigned int level = 0) const;
3734 
3753  begin_line(const unsigned int level = 0) const;
3754 
3773  begin_active_line(const unsigned int level = 0) const;
3774 
3780  end_line() const;
3781 
3782  /*
3783  * @}
3784  */
3785 
3810  raw_quad_iterator
3811  begin_raw_quad(const unsigned int level = 0) const;
3812 
3831  begin_quad(const unsigned int level = 0) const;
3832 
3851  begin_active_quad(const unsigned int level = 0) const;
3852 
3858  end_quad() const;
3859 
3860  /*
3861  * @}
3862  */
3863 
3887  raw_hex_iterator
3888  begin_raw_hex(const unsigned int level = 0) const;
3889 
3907  hex_iterator
3908  begin_hex(const unsigned int level = 0) const;
3909 
3928  begin_active_hex(const unsigned int level = 0) const;
3929 
3934  hex_iterator
3935  end_hex() const;
3936 
3937  /*
3938  * @}
3939  */
3940 
3941 
3955  void
3956  clear_despite_subscriptions();
3957 
3964  void
3965  reset_active_cell_indices();
3966 
3980  DistortedCellList
3981  execute_refinement();
3982 
3989  void
3990  execute_coarsening();
3991 
3996  void
3997  fix_coarsen_flags();
3998 
4003  std::vector<std::unique_ptr<
4004  ::internal::TriangulationImplementation::TriaLevel<dim>>>
4005  levels;
4006 
4012  std::unique_ptr<::internal::TriangulationImplementation::TriaFaces<dim>>
4013  faces;
4014 
4015 
4019  std::vector<Point<spacedim>> vertices;
4020 
4024  std::vector<bool> vertices_used;
4025 
4030  std::map<types::manifold_id, std::unique_ptr<const Manifold<dim, spacedim>>>
4031  manifold;
4032 
4036  bool anisotropic_refinement;
4037 
4038 
4043  const bool check_for_distorted_cells;
4044 
4054  ::internal::TriangulationImplementation::NumberCache<dim> number_cache;
4055 
4070  std::unique_ptr<std::map<unsigned int, types::boundary_id>>
4071  vertex_to_boundary_id_map_1d;
4072 
4073 
4093  std::unique_ptr<std::map<unsigned int, types::manifold_id>>
4094  vertex_to_manifold_id_map_1d;
4095 
4096  // make a couple of classes friends
4097  template <int, int, int>
4098  friend class TriaAccessorBase;
4099  template <int, int, int>
4100  friend class TriaAccessor;
4101  friend class TriaAccessor<0, 1, spacedim>;
4102 
4103  friend class CellAccessor<dim, spacedim>;
4104 
4105  friend struct ::internal::TriaAccessorImplementation::Implementation;
4106 
4107  friend class hp::DoFHandler<dim, spacedim>;
4108 
4109  friend struct ::internal::TriangulationImplementation::Implementation;
4110 
4111  template <typename>
4112  friend class ::internal::TriangulationImplementation::TriaObjects;
4113 
4114  // explicitly check for sensible template arguments, but not on windows
4115  // because MSVC creates bogus warnings during normal compilation
4116 #ifndef DEAL_II_MSVC
4117  static_assert(dim <= spacedim,
4118  "The dimension <dim> of a Triangulation must be less than or "
4119  "equal to the space dimension <spacedim> in which it lives.");
4120 #endif
4121 };
4122 
4123 
4124 #ifndef DOXYGEN
4125 
4126 
4127 template <int structdim>
4128 template <class Archive>
4129 void
4130 CellData<structdim>::serialize(Archive &ar, const unsigned int /*version*/)
4131 {
4132  ar &vertices;
4133  ar &material_id;
4134  ar &boundary_id;
4135  ar &manifold_id;
4136 }
4137 
4138 
4139 
4140 namespace internal
4141 {
4142  namespace TriangulationImplementation
4143  {
4144  template <class Archive>
4145  void
4146  NumberCache<1>::serialize(Archive &ar, const unsigned int)
4147  {
4148  ar &n_levels;
4149  ar &n_lines &n_lines_level;
4150  ar &n_active_lines &n_active_lines_level;
4151  }
4152 
4153 
4154  template <class Archive>
4155  void
4156  NumberCache<2>::serialize(Archive &ar, const unsigned int version)
4157  {
4158  this->NumberCache<1>::serialize(ar, version);
4159 
4160  ar &n_quads &n_quads_level;
4161  ar &n_active_quads &n_active_quads_level;
4162  }
4163 
4164 
4165  template <class Archive>
4166  void
4167  NumberCache<3>::serialize(Archive &ar, const unsigned int version)
4168  {
4169  this->NumberCache<2>::serialize(ar, version);
4170 
4171  ar &n_hexes &n_hexes_level;
4172  ar &n_active_hexes &n_active_hexes_level;
4173  }
4174 
4175  } // namespace TriangulationImplementation
4176 } // namespace internal
4177 
4178 
4179 template <int dim, int spacedim>
4180 inline bool
4181 Triangulation<dim, spacedim>::vertex_used(const unsigned int index) const
4182 {
4183  Assert(index < vertices_used.size(),
4184  ExcIndexRange(index, 0, vertices_used.size()));
4185  return vertices_used[index];
4186 }
4187 
4188 
4189 
4190 template <int dim, int spacedim>
4191 inline unsigned int
4193 {
4194  return number_cache.n_levels;
4195 }
4196 
4197 template <int dim, int spacedim>
4198 inline unsigned int
4200 {
4201  return number_cache.n_levels;
4202 }
4203 
4204 
4205 template <int dim, int spacedim>
4206 inline unsigned int
4208 {
4209  return vertices.size();
4210 }
4211 
4212 
4213 
4214 template <int dim, int spacedim>
4215 inline const std::vector<Point<spacedim>> &
4217 {
4218  return vertices;
4219 }
4220 
4221 
4222 template <int dim, int spacedim>
4223 template <class Archive>
4224 void
4225 Triangulation<dim, spacedim>::save(Archive &ar, const unsigned int) const
4226 {
4227  // as discussed in the documentation, do not store the signals as
4228  // well as boundary and manifold description but everything else
4229  ar &smooth_grid;
4230 
4231  unsigned int n_levels = levels.size();
4232  ar & n_levels;
4233  for (unsigned int i = 0; i < levels.size(); ++i)
4234  ar &levels[i];
4235 
4236  // boost dereferences a nullptr when serializing a nullptr
4237  // at least up to 1.65.1. This causes problems with clang-5.
4238  // Therefore, work around it.
4239  bool faces_is_nullptr = (faces.get() == nullptr);
4240  ar & faces_is_nullptr;
4241  if (!faces_is_nullptr)
4242  ar &faces;
4243 
4244  ar &vertices;
4245  ar &vertices_used;
4246 
4247  ar &anisotropic_refinement;
4248  ar &number_cache;
4249 
4250  ar &check_for_distorted_cells;
4251 
4252  if (dim == 1)
4253  {
4254  ar &vertex_to_boundary_id_map_1d;
4255  ar &vertex_to_manifold_id_map_1d;
4256  }
4257 }
4258 
4259 
4260 
4261 template <int dim, int spacedim>
4262 template <class Archive>
4263 void
4264 Triangulation<dim, spacedim>::load(Archive &ar, const unsigned int)
4265 {
4266  // clear previous content. this also calls the respective signal
4267  clear();
4268 
4269  // as discussed in the documentation, do not store the signals as
4270  // well as boundary and manifold description but everything else
4271  ar &smooth_grid;
4272 
4273  unsigned int size;
4274  ar & size;
4275  levels.resize(size);
4276  for (unsigned int i = 0; i < levels.size(); ++i)
4277  {
4278  std::unique_ptr<internal::TriangulationImplementation::TriaLevel<dim>>
4279  level;
4280  ar &level;
4281  levels[i] = std::move(level);
4282  }
4283 
4284  // Workaround for nullptr, see in save().
4285  bool faces_is_nullptr = true;
4286  ar & faces_is_nullptr;
4287  if (!faces_is_nullptr)
4288  ar &faces;
4289 
4290  ar &vertices;
4291  ar &vertices_used;
4292 
4293  ar &anisotropic_refinement;
4294  ar &number_cache;
4295 
4296  // the levels do not serialize the active_cell_indices because
4297  // they are easy enough to rebuild upon re-loading data. do
4298  // this here. don't forget to first resize the fields appropriately
4299  {
4300  for (unsigned int l = 0; l < levels.size(); ++l)
4301  levels[l]->active_cell_indices.resize(levels[l]->refine_flags.size());
4302  reset_active_cell_indices();
4303  }
4304 
4305 
4306  bool my_check_for_distorted_cells;
4307  ar & my_check_for_distorted_cells;
4308 
4309  Assert(my_check_for_distorted_cells == check_for_distorted_cells,
4310  ExcMessage("The triangulation loaded into here must have the "
4311  "same setting with regard to reporting distorted "
4312  "cell as the one previously stored."));
4313 
4314  if (dim == 1)
4315  {
4316  ar &vertex_to_boundary_id_map_1d;
4317  ar &vertex_to_manifold_id_map_1d;
4318  }
4319 
4320  // trigger the create signal to indicate
4321  // that new content has been imported into
4322  // the triangulation
4323  signals.create();
4324 }
4325 
4326 
4327 
4328 template <int dim, int spacedim>
4329 inline unsigned int
4331  const types::coarse_cell_id coarse_cell_id) const
4332 {
4333  return coarse_cell_id;
4334 }
4335 
4336 
4337 
4338 template <int dim, int spacedim>
4339 inline types::coarse_cell_id
4341  const unsigned int coarse_cell_index) const
4342 {
4343  return coarse_cell_index;
4344 }
4345 
4346 
4347 
4348 /* -------------- declaration of explicit specializations ------------- */
4349 
4350 template <>
4351 unsigned int
4352 Triangulation<1, 1>::n_raw_lines(const unsigned int level) const;
4353 template <>
4354 unsigned int
4356 template <>
4357 unsigned int
4358 Triangulation<1, 1>::n_quads(const unsigned int level) const;
4359 template <>
4360 unsigned int
4361 Triangulation<1, 1>::n_raw_quads(const unsigned int level) const;
4362 template <>
4363 unsigned int
4364 Triangulation<2, 2>::n_raw_quads(const unsigned int level) const;
4365 template <>
4366 unsigned int
4367 Triangulation<1, 1>::n_raw_hexs(const unsigned int level) const;
4368 template <>
4369 unsigned int
4370 Triangulation<1, 1>::n_active_quads(const unsigned int level) const;
4371 template <>
4372 unsigned int
4374 template <>
4375 unsigned int
4377 
4378 
4379 // -------------------------------------------------------------------
4380 // -- Explicit specializations for codimension one grids
4381 
4382 
4383 template <>
4384 unsigned int
4385 Triangulation<1, 2>::n_raw_lines(const unsigned int level) const;
4386 template <>
4387 unsigned int
4389 template <>
4390 unsigned int
4391 Triangulation<1, 2>::n_quads(const unsigned int level) const;
4392 template <>
4393 unsigned int
4394 Triangulation<1, 2>::n_raw_quads(const unsigned int level) const;
4395 template <>
4396 unsigned int
4397 Triangulation<2, 3>::n_raw_quads(const unsigned int level) const;
4398 template <>
4399 unsigned int
4400 Triangulation<1, 2>::n_raw_hexs(const unsigned int level) const;
4401 template <>
4402 unsigned int
4403 Triangulation<1, 2>::n_active_quads(const unsigned int level) const;
4404 template <>
4405 unsigned int
4407 template <>
4408 unsigned int
4410 
4411 // -------------------------------------------------------------------
4412 // -- Explicit specializations for codimension two grids
4413 
4414 
4415 template <>
4416 unsigned int
4417 Triangulation<1, 3>::n_raw_lines(const unsigned int level) const;
4418 template <>
4419 unsigned int
4421 template <>
4422 unsigned int
4423 Triangulation<1, 3>::n_quads(const unsigned int level) const;
4424 template <>
4425 unsigned int
4426 Triangulation<1, 3>::n_raw_quads(const unsigned int level) const;
4427 template <>
4428 unsigned int
4429 Triangulation<2, 3>::n_raw_quads(const unsigned int level) const;
4430 template <>
4431 unsigned int
4432 Triangulation<1, 3>::n_raw_hexs(const unsigned int level) const;
4433 template <>
4434 unsigned int
4435 Triangulation<1, 3>::n_active_quads(const unsigned int level) const;
4436 template <>
4437 unsigned int
4439 template <>
4440 unsigned int
4442 
4443 
4444 #endif // DOXYGEN
4445 
4446 DEAL_II_NAMESPACE_CLOSE
4447 
4448 // Include tria_accessor.h here, so that it is possible for an end
4449 // user to use the iterators of Triangulation<dim> directly without
4450 // the need to include tria_accessor.h separately. (Otherwise the
4451 // iterators are an 'opaque' or 'incomplete' type.)
4452 #include <deal.II/grid/tria_accessor.h>
4453 
4454 #endif
std::vector< CellData< 1 > > boundary_lines
Definition: tria.h:273
boost::signals2::signal< void(const Triangulation< dim, spacedim > &destination_tria)> copy
Definition: tria.h:2276
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:2228
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:2256
#define DeclException2(Exception2, type1, type2, outsequence)
Definition: exceptions.h:541
Definition: tria.h:136
boost::signals2::signal< void()> post_distributed_repartition
Definition: tria.h:2372
typename IteratorSelector::line_iterator line_iterator
Definition: tria.h:1607
typename IteratorSelector::hex_iterator hex_iterator
Definition: tria.h:1655
typename IteratorSelector::active_line_iterator active_line_iterator
Definition: tria.h:1622
boost::signals2::signal< void()> pre_distributed_refinement
Definition: tria.h:2343
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:2221
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:2408
boost::signals2::signal< void()> clear
Definition: tria.h:2291
::internal::TriangulationImplementation::Iterators< dim, spacedim > IteratorSelector
Definition: tria.h:1300
unsigned int n_raw_hexs(const unsigned int level) const
Definition: tria.cc:13109
unsigned int max_adjacent_cells() const
Definition: tria.cc:13232
void load(Archive &ar, const unsigned int version)
std::list< typename Triangulation< dim, spacedim >::cell_iterator > distorted_cells
Definition: tria.h:1701
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:2266
typename IteratorSelector::active_quad_iterator active_quad_iterator
Definition: tria.h:1646
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:2245
unsigned int n_raw_lines() const
Definition: tria.cc:12818
boost::signals2::signal< void()> pre_distributed_save
Definition: tria.h:2380
void save(Archive &ar, const unsigned int version) const
typename IteratorSelector::quad_iterator quad_iterator
Definition: tria.h:1631
boost::signals2::signal< void()> create
Definition: tria.h:2212
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:2236
typename IteratorSelector::active_hex_iterator active_hex_iterator
Definition: tria.h:1666
unsigned int global_dof_index
Definition: types.h:89
boost::signals2::signal< void()> post_distributed_save
Definition: tria.h:2387
std::vector< unsigned int > n_active_hexes_level
Definition: tria.h:482
virtual unsigned int n_global_levels() const
std::vector< CellData< 2 > > boundary_quads
Definition: tria.h:281
boost::signals2::signal< unsigned int(const cell_iterator &, const CellStatus), CellWeightSum< unsigned int > > cell_weight
Definition: tria.h:2330
boost::signals2::signal< void()> pre_distributed_load
Definition: tria.h:2395
std::vector< unsigned int > n_active_lines_level
Definition: tria.h:363
std::vector< unsigned int > n_active_quads_level
Definition: tria.h:422
boost::signals2::signal< void()> post_distributed_load
Definition: tria.h:2402
unsigned int boundary_id
Definition: types.h:125
boost::signals2::signal< void()> pre_distributed_repartition
Definition: tria.h:2365
bool vertex_used(const unsigned int index) const
boost::signals2::signal< void()> any_change
Definition: tria.h:2302
void serialize(Archive &ar, const unsigned int version)
boost::signals2::signal< void()> post_distributed_refinement
Definition: tria.h:2353