Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Modules Pages
tria.h
1 // ---------------------------------------------------------------------
2 // @f$Id: tria.h 31641 2013-11-13 16:40:04Z bangerth @f$
3 //
4 // Copyright (C) 1998 - 2013 by the deal.II authors
5 //
6 // This file is part of the deal.II library.
7 //
8 // The deal.II library is free software; you can use it, redistribute
9 // it, and/or modify it under the terms of the GNU Lesser General
10 // Public License as published by the Free Software Foundation; either
11 // version 2.1 of the License, or (at your option) any later version.
12 // The full text of the license can be found in the file LICENSE at
13 // the top level of the deal.II distribution.
14 //
15 // ---------------------------------------------------------------------
16 
17 #ifndef __deal2__tria_h
18 #define __deal2__tria_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/base/point.h>
23 #include <deal.II/base/subscriptor.h>
24 #include <deal.II/base/smartpointer.h>
25 #include <deal.II/base/geometry_info.h>
26 #include <deal.II/base/std_cxx1x/function.h>
27 #include <deal.II/grid/tria_iterator_selector.h>
28 #include <deal.II/grid/tria_faces.h>
29 #include <deal.II/grid/tria_levels.h>
30 
31 #include <boost/signals2.hpp>
32 #include <boost/serialization/vector.hpp>
33 #include <boost/serialization/map.hpp>
34 #include <boost/serialization/split_member.hpp>
35 
36 #include <vector>
37 #include <list>
38 #include <map>
39 
40 
42 
43 template <int dim, int spacedim> class Boundary;
44 template <int dim, int spacedim> class StraightBoundary;
45 
46 template <int, int, int> class TriaAccessor;
47 template <int spacedim> class TriaAccessor<0,1,spacedim>;
48 
49 namespace internal
50 {
51  namespace Triangulation
52  {
53  template <int dim> class TriaLevel;
54  template <int dim> class TriaFaces;
55 
56  template <typename> class TriaObjects;
57 
65  struct Implementation;
66  }
67 
68  namespace TriaAccessor
69  {
70  struct Implementation;
71  }
72 }
73 
74 template <int dim, int spacedim> class DoFHandler;
75 namespace hp
76 {
77  template <int dim, int spacedim> class DoFHandler;
78 }
79 template <int dim, int spacedim> class MGDoFHandler;
80 
81 
82 /*------------------------------------------------------------------------*/
83 
97 template <int structdim>
98 struct CellData
99 {
105 
118  union
119  {
120  types::boundary_id boundary_id;
121  types::material_id material_id;
122  };
123 
128  CellData ();
129 };
130 
131 
132 
159 {
165  std::vector<CellData<1> > boundary_lines;
166 
172  std::vector<CellData<2> > boundary_quads;
173 
187  bool check_consistency (const unsigned int dim) const;
188 };
189 
190 
191 /*------------------------------------------------------------------------*/
192 
193 
194 namespace internal
195 {
200  namespace Triangulation
201  {
202 
217  template <int dim>
218  struct NumberCache
219  {
220  };
221 
235  template <>
236  struct NumberCache<1>
237  {
243  unsigned int n_levels;
244 
249  unsigned int n_lines;
250 
255  std::vector<unsigned int> n_lines_level;
256 
261  unsigned int n_active_lines;
262 
267  std::vector<unsigned int> n_active_lines_level;
268 
273  NumberCache ();
274 
280  std::size_t memory_consumption () const;
281 
286  template <class Archive>
287  void serialize (Archive &ar,
288  const unsigned int version);
289  };
290 
291 
307  template <>
308  struct NumberCache<2> : public NumberCache<1>
309  {
314  unsigned int n_quads;
315 
320  std::vector<unsigned int> n_quads_level;
321 
326  unsigned int n_active_quads;
327 
332  std::vector<unsigned int> n_active_quads_level;
333 
338  NumberCache ();
339 
345  std::size_t memory_consumption () const;
346 
351  template <class Archive>
352  void serialize (Archive &ar,
353  const unsigned int version);
354  };
355 
356 
372  template <>
373  struct NumberCache<3> : public NumberCache<2>
374  {
379  unsigned int n_hexes;
380 
385  std::vector<unsigned int> n_hexes_level;
386 
391  unsigned int n_active_hexes;
392 
397  std::vector<unsigned int> n_active_hexes_level;
398 
403  NumberCache ();
404 
410  std::size_t memory_consumption () const;
411 
416  template <class Archive>
417  void serialize (Archive &ar,
418  const unsigned int version);
419  };
420  }
421 }
422 
423 
424 /*------------------------------------------------------------------------*/
425 
426 
1215 template <int dim, int spacedim=dim>
1216 class Triangulation : public Subscriptor
1217 {
1218 private:
1219 
1225  typedef ::internal::Triangulation::Iterators<dim, spacedim> IteratorSelector;
1226 
1227 public:
1235 
1244  {
1248  none = 0x0,
1422 
1436 
1443  };
1444 
1445  typedef TriaIterator <CellAccessor<dim,spacedim> > cell_iterator;
1446  typedef TriaActiveIterator<CellAccessor<dim,spacedim> > active_cell_iterator;
1447 
1448  typedef TriaIterator <TriaAccessor<dim-1, dim, spacedim> > face_iterator;
1449  typedef TriaActiveIterator<TriaAccessor<dim-1, dim, spacedim> > active_face_iterator;
1450 
1451  typedef typename IteratorSelector::line_iterator line_iterator;
1452  typedef typename IteratorSelector::active_line_iterator active_line_iterator;
1453 
1454  typedef typename IteratorSelector::quad_iterator quad_iterator;
1455  typedef typename IteratorSelector::active_quad_iterator active_quad_iterator;
1456 
1457  typedef typename IteratorSelector::hex_iterator hex_iterator;
1458  typedef typename IteratorSelector::active_hex_iterator active_hex_iterator;
1459 
1476  {
1477  public:
1490  virtual ~RefinementListener ();
1491 
1506  virtual
1507  void
1509 
1524  virtual
1525  void
1527 
1549  virtual
1550  void
1552  const Triangulation<dim, spacedim> &new_tria);
1553 
1572  virtual
1573  void
1575  };
1576 
1612  {
1624  virtual ~DistortedCellList () throw();
1625 
1633  std::list<typename Triangulation<dim,spacedim>::cell_iterator>
1635  };
1636 
1637 
1642  static const unsigned int dimension = dim;
1643 
1648  static const unsigned int space_dimension = spacedim;
1649 
1677  const bool check_for_distorted_cells = false);
1678 
1712 
1717  virtual ~Triangulation ();
1718 
1728  virtual void clear ();
1729 
1738  virtual void set_mesh_smoothing (const MeshSmoothing mesh_smoothing);
1739 
1799  void set_boundary (const types::boundary_id number,
1800  const Boundary<dim,spacedim> &boundary_object);
1801 
1815  void set_boundary (const types::boundary_id number);
1816 
1828  const Boundary<dim,spacedim> &get_boundary (const types::boundary_id number) const;
1829 
1847  std::vector<types::boundary_id> get_boundary_indicators() const;
1848 
1898  virtual void copy_triangulation (const Triangulation<dim, spacedim> &old_tria);
1899 
1985  virtual void create_triangulation (const std::vector<Point<spacedim> > &vertices,
1986  const std::vector<CellData<dim> > &cells,
1987  const SubCellData &subcelldata);
1988 
2005  virtual void create_triangulation_compatibility (
2006  const std::vector<Point<spacedim> > &vertices,
2007  const std::vector<CellData<dim> > &cells,
2008  const SubCellData &subcelldata);
2009 
2020  void flip_all_direction_flags();
2021 
2042  void distort_random (const double factor,
2043  const bool keep_boundary=true) DEAL_II_DEPRECATED;
2044 
2045 
2062  void set_all_refine_flags ();
2063 
2087  void refine_global (const unsigned int times = 1);
2088 
2146  virtual void execute_coarsening_and_refinement ();
2147 
2223  void add_refinement_listener (RefinementListener &listener) const DEAL_II_DEPRECATED;
2224 
2241  void remove_refinement_listener (RefinementListener &listener) const DEAL_II_DEPRECATED;
2242 
2251  struct Signals
2252  {
2253  boost::signals2::signal<void ()> create;
2254  boost::signals2::signal<void ()> pre_refinement;
2255  boost::signals2::signal<void ()> post_refinement;
2256  boost::signals2::signal<void (const Triangulation<dim, spacedim> &original_tria)> copy;
2257  boost::signals2::signal<void ()> clear;
2258  boost::signals2::signal<void ()> any_change;
2259  };
2260 
2264  mutable Signals signals;
2265 
2279  void save_refine_flags (std::ostream &out) const;
2280 
2285  void save_refine_flags (std::vector<bool> &v) const;
2286 
2291  void load_refine_flags (std::istream &in);
2292 
2297  void load_refine_flags (const std::vector<bool> &v);
2298 
2302  void save_coarsen_flags (std::ostream &out) const;
2303 
2308  void save_coarsen_flags (std::vector<bool> &v) const;
2309 
2313  void load_coarsen_flags (std::istream &out);
2314 
2318  void load_coarsen_flags (const std::vector<bool> &v);
2319 
2325  bool get_anisotropic_refinement_flag() const;
2337  void clear_user_flags ();
2338 
2347  void save_user_flags (std::ostream &out) const;
2348 
2356  void save_user_flags (std::vector<bool> &v) const;
2357 
2363  void load_user_flags (std::istream &in);
2364 
2370  void load_user_flags (const std::vector<bool> &v);
2371 
2376  void clear_user_flags_line ();
2377 
2382  void save_user_flags_line (std::ostream &out) const;
2383 
2391  void save_user_flags_line (std::vector<bool> &v) const;
2392 
2397  void load_user_flags_line (std::istream &in);
2398 
2403  void load_user_flags_line (const std::vector<bool> &v);
2404 
2409  void clear_user_flags_quad ();
2410 
2415  void save_user_flags_quad (std::ostream &out) const;
2416 
2424  void save_user_flags_quad (std::vector<bool> &v) const;
2425 
2430  void load_user_flags_quad (std::istream &in);
2431 
2436  void load_user_flags_quad (const std::vector<bool> &v);
2437 
2438 
2443  void clear_user_flags_hex ();
2444 
2449  void save_user_flags_hex (std::ostream &out) const;
2450 
2458  void save_user_flags_hex (std::vector<bool> &v) const;
2459 
2464  void load_user_flags_hex (std::istream &in);
2465 
2470  void load_user_flags_hex (const std::vector<bool> &v);
2471 
2478  void clear_user_data ();
2479 
2488 
2495  void save_user_indices (std::vector<unsigned int> &v) const;
2496 
2502  void load_user_indices (const std::vector<unsigned int> &v);
2503 
2510  void save_user_pointers (std::vector<void *> &v) const;
2511 
2517  void load_user_pointers (const std::vector<void *> &v);
2518 
2525  void save_user_indices_line (std::vector<unsigned int> &v) const;
2526 
2532  void load_user_indices_line (const std::vector<unsigned int> &v);
2533 
2540  void save_user_indices_quad (std::vector<unsigned int> &v) const;
2541 
2547  void load_user_indices_quad (const std::vector<unsigned int> &v);
2548 
2555  void save_user_indices_hex (std::vector<unsigned int> &v) const;
2556 
2562  void load_user_indices_hex (const std::vector<unsigned int> &v);
2569  void save_user_pointers_line (std::vector<void *> &v) const;
2570 
2576  void load_user_pointers_line (const std::vector<void *> &v);
2577 
2584  void save_user_pointers_quad (std::vector<void *> &v) const;
2585 
2591  void load_user_pointers_quad (const std::vector<void *> &v);
2592 
2599  void save_user_pointers_hex (std::vector<void *> &v) const;
2600 
2606  void load_user_pointers_hex (const std::vector<void *> &v);
2618  cell_iterator begin (const unsigned int level = 0) const;
2619 
2632  active_cell_iterator begin_active(const unsigned int level = 0) const;
2633 
2640  cell_iterator end () const;
2641 
2648  cell_iterator end (const unsigned int level) const;
2649 
2655  active_cell_iterator end_active (const unsigned int level) const;
2656 
2657 
2662  cell_iterator last () const;
2663 
2668  active_cell_iterator last_active () const;
2671  /*---------------------------------------*/
2672  /*---------------------------------------*/
2673 
2681  face_iterator begin_face () const;
2682 
2687  active_face_iterator begin_active_face() const;
2688 
2695  face_iterator end_face () const;
2699  /*---------------------------------------*/
2700 
2705 
2727  unsigned int n_lines () const;
2728 
2733  unsigned int n_lines (const unsigned int level) const;
2734 
2738  unsigned int n_active_lines () const;
2739 
2744  unsigned int n_active_lines (const unsigned int level) const;
2745 
2750  unsigned int n_quads () const;
2751 
2756  unsigned int n_quads (const unsigned int level) const;
2757 
2762  unsigned int n_active_quads () const;
2763 
2768  unsigned int n_active_quads (const unsigned int level) const;
2769 
2774  unsigned int n_hexs() const;
2775 
2781  unsigned int n_hexs(const unsigned int level) const;
2782 
2787  unsigned int n_active_hexs() const;
2788 
2794  unsigned int n_active_hexs(const unsigned int level) const;
2795 
2802  unsigned int n_cells () const;
2803 
2810  unsigned int n_cells (const unsigned int level) const;
2811 
2817  unsigned int n_active_cells () const;
2818 
2825  unsigned int n_active_cells (const unsigned int level) const;
2826 
2835  unsigned int n_faces () const;
2836 
2845  unsigned int n_active_faces () const;
2846 
2857  unsigned int n_levels () const;
2858 
2865  virtual unsigned int n_global_levels () const;
2866 
2881  unsigned int n_vertices () const;
2882 
2900  const std::vector<Point<spacedim> > &
2901  get_vertices () const;
2902 
2909  unsigned int n_used_vertices () const;
2910 
2915  bool vertex_used (const unsigned int index) const;
2916 
2924  const std::vector<bool> &
2925  get_used_vertices () const;
2926 
2950  unsigned int max_adjacent_cells () const;
2951 
2963  virtual types::subdomain_id locally_owned_subdomain () const;
2970 
2988  unsigned int n_raw_lines () const;
2989 
3007  unsigned int n_raw_lines (const unsigned int level) const;
3008 
3026  unsigned int n_raw_quads () const;
3027 
3045  unsigned int n_raw_quads (const unsigned int level) const;
3046 
3064  unsigned int n_raw_hexs () const;
3065 
3083  unsigned int n_raw_hexs (const unsigned int level) const;
3084 
3102  unsigned int n_raw_cells (const unsigned int level) const;
3103 
3123  unsigned int n_raw_faces () const;
3124 
3139  virtual std::size_t memory_consumption () const;
3140 
3153  template <class Archive>
3154  void save (Archive &ar,
3155  const unsigned int version) const;
3156 
3181  template <class Archive>
3182  void load (Archive &ar,
3183  const unsigned int version);
3184 
3185  BOOST_SERIALIZATION_SPLIT_MEMBER()
3186 
3195  DeclException1 (ExcInvalidLevel,
3196  int,
3197  << "The given level " << arg1
3198  << " is not in the valid range!");
3208  DeclException0 (ExcTriangulationNotEmpty);
3214  DeclException0 (ExcGridReadError);
3219  DeclException0 (ExcFacesHaveNoLevel);
3226  DeclException1 (ExcEmptyLevel,
3227  int,
3228  << "You tried to do something on level " << arg1
3229  << ", but this level is empty.");
3234  DeclException0 (ExcNonOrientableTriangulation);
3236 protected:
3244 
3271  static void write_bool_vector (const unsigned int magic_number1,
3272  const std::vector<bool> &v,
3273  const unsigned int magic_number2,
3274  std::ostream &out);
3275 
3281  static void read_bool_vector (const unsigned int magic_number1,
3282  std::vector<bool> &v,
3283  const unsigned int magic_number2,
3284  std::istream &in);
3285 
3286 private:
3305  typedef TriaRawIterator <TriaAccessor<dim-1, dim, spacedim> > raw_face_iterator;
3306  typedef typename IteratorSelector::raw_line_iterator raw_line_iterator;
3307  typedef typename IteratorSelector::raw_quad_iterator raw_quad_iterator;
3308  typedef typename IteratorSelector::raw_hex_iterator raw_hex_iterator;
3309 
3316  raw_cell_iterator begin_raw (const unsigned int level = 0) const;
3317 
3324  raw_cell_iterator end_raw (const unsigned int level) const;
3327  /*---------------------------------------*/
3328 
3333 
3344  raw_line_iterator
3345  begin_raw_line (const unsigned int level = 0) const;
3346 
3351  line_iterator
3352  begin_line (const unsigned int level = 0) const;
3353 
3358  active_line_iterator
3359  begin_active_line(const unsigned int level = 0) const;
3360 
3367  line_iterator end_line () const;
3370  /*---------------------------------------*/
3371 
3375 
3386  raw_quad_iterator
3387  begin_raw_quad (const unsigned int level = 0) const;
3388 
3393  quad_iterator
3394  begin_quad (const unsigned int level = 0) const;
3395 
3400  active_quad_iterator
3401  begin_active_quad (const unsigned int level = 0) const;
3402 
3409  quad_iterator
3410  end_quad () const;
3413  /*---------------------------------------*/
3414 
3418 
3426  raw_hex_iterator
3427  begin_raw_hex (const unsigned int level = 0) const;
3428 
3433  hex_iterator
3434  begin_hex (const unsigned int level = 0) const;
3435 
3440  active_hex_iterator
3441  begin_active_hex (const unsigned int level = 0) const;
3442 
3449  hex_iterator
3450  end_hex () const;
3478 
3499  DistortedCellList execute_refinement ();
3500 
3509  void execute_coarsening ();
3510 
3516  void fix_coarsen_flags ();
3517 
3523  std::vector<::internal::Triangulation::TriaLevel<dim>*> levels;
3524 
3532  ::internal::Triangulation::TriaFaces<dim> *faces;
3533 
3534 
3539  std::vector<Point<spacedim> > vertices;
3540 
3545  std::vector<bool> vertices_used;
3546 
3553  std::map<types::boundary_id, SmartPointer<const Boundary<dim, spacedim> , Triangulation<dim, spacedim> > > boundary;
3554 
3561 
3562 
3570 
3587  ::internal::Triangulation::NumberCache<dim> number_cache;
3588 
3612  std::map<unsigned int, types::boundary_id> *vertex_to_boundary_id_map_1d;
3613 
3625  mutable
3626  std::multimap<const RefinementListener *, std::vector<boost::signals2::connection> >
3628 
3629  // make a couple of classes
3630  // friends
3631  template <int,int,int> friend class TriaAccessorBase;
3632  template <int,int,int> friend class TriaAccessor;
3633  friend class TriaAccessor<0, 1, spacedim>;
3634 
3635  friend class CellAccessor<dim, spacedim>;
3636 
3637  friend struct ::internal::TriaAccessor::Implementation;
3638 
3639  friend class hp::DoFHandler<dim,spacedim>;
3640 
3641  friend struct ::internal::Triangulation::Implementation;
3642 
3643  template <typename>
3644  friend class ::internal::Triangulation::TriaObjects;
3645 };
3646 
3647 
3648 #ifndef DOXYGEN
3649 
3650 
3651 template <int structdim>
3652 inline
3654 {
3655  for (unsigned int i=0; i<GeometryInfo<structdim>::vertices_per_cell; ++i)
3656  vertices[i] = numbers::invalid_unsigned_int;
3657 
3658  material_id = 0;
3659 }
3660 
3661 
3662 
3663 namespace internal
3664 {
3665  namespace Triangulation
3666  {
3667  template <class Archive>
3668  void NumberCache<1>::serialize (Archive &ar,
3669  const unsigned int)
3670  {
3671  ar &n_levels;
3672  ar &n_lines &n_lines_level;
3673  ar &n_active_lines &n_active_lines_level;
3674  }
3675 
3676 
3677  template <class Archive>
3678  void NumberCache<2>::serialize (Archive &ar,
3679  const unsigned int version)
3680  {
3681  this->NumberCache<1>::serialize (ar, version);
3682 
3683  ar &n_quads &n_quads_level;
3684  ar &n_active_quads &n_active_quads_level;
3685  }
3686 
3687 
3688  template <class Archive>
3689  void NumberCache<3>::serialize (Archive &ar,
3690  const unsigned int version)
3691  {
3692  this->NumberCache<2>::serialize (ar, version);
3693 
3694  ar &n_hexes &n_hexes_level;
3695  ar &n_active_hexes &n_active_hexes_level;
3696  }
3697 
3698  }
3699 }
3700 
3701 
3702 template <int dim, int spacedim>
3703 inline
3704 bool
3705 Triangulation<dim,spacedim>::vertex_used(const unsigned int index) const
3706 {
3707  Assert (index < vertices_used.size(),
3708  ExcIndexRange(index, 0, vertices_used.size()));
3709  return vertices_used[index];
3710 }
3711 
3712 
3713 
3714 template <int dim, int spacedim>
3715 inline
3716 unsigned int Triangulation<dim, spacedim>::n_levels () const
3717 {
3718  return number_cache.n_levels;
3719 }
3720 
3721 template <int dim, int spacedim>
3722 inline
3724 {
3725  return number_cache.n_levels;
3726 }
3727 
3728 
3729 template <int dim, int spacedim>
3730 inline
3731 unsigned int
3733 {
3734  return vertices.size();
3735 }
3736 
3737 
3738 
3739 template <int dim, int spacedim>
3740 inline
3741 const std::vector<Point<spacedim> > &
3743 {
3744  return vertices;
3745 }
3746 
3747 
3748 template <int dim, int spacedim>
3749 template <class Archive>
3750 void
3752  const unsigned int) const
3753 {
3754  // as discussed in the documentation, do
3755  // not store the signals as well as
3756  // boundary and manifold descrption
3757  // but everything else
3758  ar &smooth_grid;
3759  ar &levels;
3760  ar &faces;
3761  ar &vertices;
3762  ar &vertices_used;
3763 
3765  ar &number_cache;
3766 
3768 
3769  if (dim == 1)
3771 }
3772 
3773 
3774 
3775 template <int dim, int spacedim>
3776 template <class Archive>
3777 void
3779  const unsigned int)
3780 {
3781  // clear previous content. this also calls
3782  // the respective signal
3783  clear ();
3784 
3785  // as discussed in the documentation, do
3786  // not store the signals as well as
3787  // boundary and manifold description
3788  // but everything else
3789  ar &smooth_grid;
3790  ar &levels;
3791  ar &faces;
3792  ar &vertices;
3793  ar &vertices_used;
3794 
3796  ar &number_cache;
3797 
3798  bool my_check_for_distorted_cells;
3799  ar &my_check_for_distorted_cells;
3800 
3801  Assert (my_check_for_distorted_cells == check_for_distorted_cells,
3802  ExcMessage ("The triangulation loaded into here must have the "
3803  "same setting with regard to reporting distorted "
3804  "cell as the one previously stored."));
3805 
3806  if (dim == 1)
3808 
3809  // trigger the create signal to indicate
3810  // that new content has been imported into
3811  // the triangulation
3812  signals.create();
3813 }
3814 
3815 
3816 /* -------------- declaration of explicit specializations ------------- */
3817 
3818 template <> unsigned int Triangulation<1,1>::n_raw_lines (const unsigned int level) const;
3819 template <> unsigned int Triangulation<1,1>::n_quads () const;
3820 template <> unsigned int Triangulation<1,1>::n_quads (const unsigned int level) const;
3821 template <> unsigned int Triangulation<1,1>::n_raw_quads (const unsigned int level) const;
3822 template <> unsigned int Triangulation<2,2>::n_raw_quads (const unsigned int level) const;
3823 template <> unsigned int Triangulation<1,1>::n_raw_hexs (const unsigned int level) const;
3824 template <> unsigned int Triangulation<1,1>::n_active_quads (const unsigned int level) const;
3825 template <> unsigned int Triangulation<1,1>::n_active_quads () const;
3826 template <> unsigned int Triangulation<1,1>::max_adjacent_cells () const;
3827 
3828 
3829 // -------------------------------------------------------------------
3830 // -- Explicit specializations for codimension one grids
3831 
3832 
3833 template <> unsigned int Triangulation<1,2>::n_raw_lines (const unsigned int level) const;
3834 template <> unsigned int Triangulation<1,2>::n_quads () const;
3835 template <> unsigned int Triangulation<1,2>::n_quads (const unsigned int level) const;
3836 template <> unsigned int Triangulation<1,2>::n_raw_quads (const unsigned int level) const;
3837 template <> unsigned int Triangulation<2,3>::n_raw_quads (const unsigned int level) const;
3838 template <> unsigned int Triangulation<1,2>::n_raw_hexs (const unsigned int level) const;
3839 template <> unsigned int Triangulation<1,2>::n_active_quads (const unsigned int level) const;
3840 template <> unsigned int Triangulation<1,2>::n_active_quads () const;
3841 template <> unsigned int Triangulation<1,2>::max_adjacent_cells () const;
3842 
3843 // -------------------------------------------------------------------
3844 // -- Explicit specializations for codimension two grids
3845 
3846 
3847 template <> unsigned int Triangulation<1,3>::n_raw_lines (const unsigned int level) const;
3848 template <> unsigned int Triangulation<1,3>::n_quads () const;
3849 template <> unsigned int Triangulation<1,3>::n_quads (const unsigned int level) const;
3850 template <> unsigned int Triangulation<1,3>::n_raw_quads (const unsigned int level) const;
3851 template <> unsigned int Triangulation<2,3>::n_raw_quads (const unsigned int level) const;
3852 template <> unsigned int Triangulation<1,3>::n_raw_hexs (const unsigned int level) const;
3853 template <> unsigned int Triangulation<1,3>::n_active_quads (const unsigned int level) const;
3854 template <> unsigned int Triangulation<1,3>::n_active_quads () const;
3855 template <> unsigned int Triangulation<1,3>::max_adjacent_cells () const;
3856 
3857 
3858 // -------------------------------------------------------------------
3859 
3860 
3861 
3862 
3863 #endif // DOXYGEN
3864 
3865 DEAL_II_NAMESPACE_CLOSE
3866 
3867 // Include tria_accessor.h here, so that it is possible for an end user to
3868 // use the iterators of Triangulation<dim> directly without the need to
3869 // include tria_accessor.h separately. (Otherwise the iterators are an
3870 // 'opaque' or 'incomplete' type.)
3871 #include <deal.II/grid/tria_accessor.h>
3872 
3873 #endif
std::vector< CellData< 1 > > boundary_lines
Definition: tria.h:165
::internal::Triangulation::Iterators< dim, spacedim > IteratorSelector
Definition: tria.h:1225
static const unsigned int dimension
Definition: tria.h:1642
void load_user_pointers_line(const std::vector< void * > &v)
static const unsigned int invalid_unsigned_int
Definition: types.h:191
unsigned int n_used_vertices() const
void save_user_pointers(std::vector< void * > &v) const
virtual types::subdomain_id locally_owned_subdomain() const
virtual void pre_refinement_notification(const Triangulation< dim, spacedim > &tria)
void fix_coarsen_flags()
Iterator points to a valid object.
std::map< unsigned int, types::boundary_id > * vertex_to_boundary_id_map_1d
Definition: tria.h:3612
cell_iterator begin(const unsigned int level=0) const
Definition: tria.h:98
bool get_anisotropic_refinement_flag() const
unsigned int n_active_cells() const
quad_iterator end_quad() const
void refine_global(const unsigned int times=1)
unsigned int n_raw_quads() const
hex_iterator begin_hex(const unsigned int level=0) const
void load_user_flags_quad(std::istream &in)
std::list< typename Triangulation< dim, spacedim >::cell_iterator > distorted_cells
Definition: tria.h:1634
unsigned char material_id
Definition: types.h:142
bool check_consistency(const unsigned int dim) const
cell_iterator last() const
face_iterator begin_face() const
void save_user_flags_line(std::ostream &out) const
::ExceptionBase & ExcMessage(std::string arg1)
line_iterator begin_line(const unsigned int level=0) const
std::vector< unsigned int > n_active_lines_level
Definition: tria.h:267
std::vector< unsigned int > n_hexes_level
Definition: tria.h:385
bool vertex_used(const unsigned int index) const
MeshSmoothing smooth_grid
Definition: tria.h:3243
std::vector< bool > vertices_used
Definition: tria.h:3545
virtual void create_triangulation(const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata)
raw_cell_iterator end_raw(const unsigned int level) const
std::vector< unsigned int > n_quads_level
Definition: tria.h:320
void load(Archive &ar, const unsigned int version)
void clear_user_flags()
STL namespace.
unsigned int max_adjacent_cells() const
void save_coarsen_flags(std::ostream &out) const
void save(Archive &ar, const unsigned int version) const
Signals signals
Definition: tria.h:2264
active_cell_iterator begin_active(const unsigned int level=0) const
unsigned int n_cells() const
std::multimap< const RefinementListener *, std::vector< boost::signals2::connection > > refinement_listener_map
Definition: tria.h:3627
std::vector< unsigned int > n_active_quads_level
Definition: tria.h:332
virtual void create_notification(const Triangulation< dim, spacedim > &tria)
void clear_user_flags_hex()
static void write_bool_vector(const unsigned int magic_number1, const std::vector< bool > &v, const unsigned int magic_number2, std::ostream &out)
bool anisotropic_refinement
Definition: tria.h:3560
void save_user_flags_quad(std::ostream &out) const
virtual void execute_coarsening_and_refinement()
virtual void post_refinement_notification(const Triangulation< dim, spacedim > &tria)
void set_all_refine_flags()
void save_user_flags(std::ostream &out) const
const Boundary< dim, spacedim > & get_boundary(const types::boundary_id number) const
void save_user_indices_line(std::vector< unsigned int > &v) const
active_line_iterator begin_active_line(const unsigned int level=0) const
void load_user_flags_line(std::istream &in)
void save_user_indices_quad(std::vector< unsigned int > &v) const
::internal::Triangulation::NumberCache< dim > number_cache
Definition: tria.h:3587
void clear_user_flags_quad()
unsigned int n_raw_faces() const
unsigned int n_lines() const
DeclException1(ExcInvalidLevel, int,<< "The given level "<< arg1<< " is not in the valid range!")
DistortedCellList execute_refinement()
void save_user_pointers_line(std::vector< void * > &v) const
void save_user_pointers_quad(std::vector< void * > &v) const
raw_quad_iterator begin_raw_quad(const unsigned int level=0) const
raw_cell_iterator begin_raw(const unsigned int level=0) const
void load_user_indices_line(const std::vector< unsigned int > &v)
active_quad_iterator begin_active_quad(const unsigned int level=0) const
void load_user_indices(const std::vector< unsigned int > &v)
Definition: types.h:29
#define Assert(cond, exc)
Definition: exceptions.h:299
const bool check_for_distorted_cells
Definition: tria.h:3569
unsigned int n_faces() const
unsigned int n_quads() const
raw_line_iterator begin_raw_line(const unsigned int level=0) const
void save_refine_flags(std::ostream &out) const
unsigned int n_raw_cells(const unsigned int level) const
void set_boundary(const types::boundary_id number, const Boundary< dim, spacedim > &boundary_object)
active_face_iterator begin_active_face() const
::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
raw_hex_iterator begin_raw_hex(const unsigned int level=0) const
face_iterator end_face() const
unsigned int n_active_faces() const
void load_user_indices_hex(const std::vector< unsigned int > &v)
std::vector<::internal::Triangulation::TriaLevel< dim > * > levels
Definition: tria.h:3523
bool prepare_coarsening_and_refinement()
virtual void create_triangulation_compatibility(const std::vector< Point< spacedim > > &vertices, const std::vector< CellData< dim > > &cells, const SubCellData &subcelldata)
BlockCompressedSparsityPattern CompressedBlockSparsityPattern DEAL_II_DEPRECATED
void add_refinement_listener(RefinementListener &listener) const DEAL_II_DEPRECATED
quad_iterator begin_quad(const unsigned int level=0) const
unsigned int n_vertices() const
static const unsigned int space_dimension
Definition: tria.h:1648
Definition: hp.h:103
void distort_random(const double factor, const bool keep_boundary=true) DEAL_II_DEPRECATED
virtual void copy_triangulation(const Triangulation< dim, spacedim > &old_tria)
std::map< types::boundary_id, SmartPointer< const Boundary< dim, spacedim >, Triangulation< dim, spacedim > > > boundary
Definition: tria.h:3553
unsigned int n_hexs() const
void save_user_flags_hex(std::ostream &out) const
void load_coarsen_flags(std::istream &out)
void clear_user_pointers() DEAL_II_DEPRECATED
const std::vector< bool > & get_used_vertices() const
void load_user_pointers(const std::vector< void * > &v)
void load_user_flags_hex(std::istream &in)
void clear_user_flags_line()
void execute_coarsening()
unsigned int n_active_quads() const
std::vector< unsigned int > n_active_hexes_level
Definition: tria.h:397
void remove_refinement_listener(RefinementListener &listener) const DEAL_II_DEPRECATED
void load_user_indices_quad(const std::vector< unsigned int > &v)
unsigned int n_levels() const
void clear_user_data()
virtual void set_mesh_smoothing(const MeshSmoothing mesh_smoothing)
active_hex_iterator begin_active_hex(const unsigned int level=0) const
void save_user_indices(std::vector< unsigned int > &v) const
::internal::Triangulation::TriaFaces< dim > * faces
Definition: tria.h:3532
void load_user_flags(std::istream &in)
std::vector< CellData< 2 > > boundary_quads
Definition: tria.h:172
std::vector< types::boundary_id > get_boundary_indicators() const
void load_refine_flags(std::istream &in)
active_cell_iterator end_active(const unsigned int level) const
void save_user_pointers_hex(std::vector< void * > &v) const
Triangulation(const MeshSmoothing smooth_grid=none, const bool check_for_distorted_cells=false)
void save_user_indices_hex(std::vector< unsigned int > &v) const
std::vector< Point< spacedim > > vertices
Definition: tria.h:3539
unsigned int n_raw_hexs() const
hex_iterator end_hex() const
static const StraightBoundary< dim, spacedim > straight_boundary
Definition: tria.h:1234
cell_iterator end() const
active_cell_iterator last_active() const
void load_user_pointers_quad(const std::vector< void * > &v)
unsigned char boundary_id
Definition: types.h:128
DeclException0(ExcTriangulationNotEmpty)
unsigned int n_active_lines() const
virtual ~Triangulation()
unsigned int n_active_hexs() const
virtual unsigned int n_global_levels() const
const std::vector< Point< spacedim > > & get_vertices() const
virtual std::size_t memory_consumption() const
line_iterator end_line() const
void load_user_pointers_hex(const std::vector< void * > &v)
void flip_all_direction_flags()
unsigned int n_raw_lines() const
virtual void clear()
virtual void copy_notification(const Triangulation< dim, spacedim > &old_tria, const Triangulation< dim, spacedim > &new_tria)
void clear_despite_subscriptions()
static void read_bool_vector(const unsigned int magic_number1, std::vector< bool > &v, const unsigned int magic_number2, std::istream &in)
unsigned int vertices[GeometryInfo< structdim >::vertices_per_cell]
Definition: tria.h:104
std::vector< unsigned int > n_lines_level
Definition: tria.h:255