Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
grid_reordering_internal.h
1 // ---------------------------------------------------------------------
2 // @f$Id: grid_reordering_internal.h 30036 2013-07-18 16:55:32Z maier @f$
3 //
4 // Copyright (C) 2003 - 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__grid_reordering_internal_h
18 #define __deal2__grid_reordering_internal_h
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/grid/tria.h>
23 
24 #include <map>
25 #include <vector>
26 
27 DEAL_II_NAMESPACE_OPEN
28 
29 
30 namespace internal
31 {
38  namespace GridReordering2d
39  {
40 
59  bool
60  is_consistent (const std::vector<CellData<2> > &cells);
61 
62 
84  {
85  public:
91  static const int EdgeToNode[4][2];
92 
99  static const int NodeToEdge[4][2];
100 
106  static const int DefaultOrientation[4][2];
107  };
108 
109 
116  class MQuad
117  {
118  public:
125  MQuad (const unsigned int v0,
126  const unsigned int v1,
127  const unsigned int v2,
128  const unsigned int v3,
129  const unsigned int s0,
130  const unsigned int s1,
131  const unsigned int s2,
132  const unsigned int s3,
133  const CellData<2> &cd);
134 
138  unsigned int v[4];
142  unsigned int side[4];
143 
150  };
151 
159  struct MSide
160  {
164  MSide (const unsigned int initv0,
165  const unsigned int initv1);
166 
172  bool operator==(const MSide &s2) const;
173 
177  bool operator!=(const MSide &s2) const;
178 
179  unsigned int v0;
180  unsigned int v1;
181  unsigned int Q0;
182  unsigned int Q1;
183 
187  unsigned int lsn0, lsn1;
188  bool Oriented;
189 
193  struct SideRectify;
194 
201  struct SideSortLess;
202  };
203 
204 
205 
213  {
214  public:
215 
220  void reorient(std::vector<CellData<2> > &quads);
221  private:
222 
231  void build_graph (const std::vector<CellData<2> > &inquads);
232 
266  void orient();
267 
273  void get_quads(std::vector<CellData<2> > &outquads) const;
274 
286  void orient_side (const unsigned int quadnum,
287  const unsigned int localsidenum);
288 
294  bool is_fully_oriented_quad (const unsigned int quadnum) const;
295 
301  bool is_oriented_side (const unsigned int quadnum,
302  const unsigned int lsn) const;
303 
309  bool is_side_default_oriented (const unsigned int qnum,
310  const unsigned int lsn) const;
311 
320  bool get_unoriented_quad (unsigned int &UnOrQLoc) const;
321 
329  bool get_unoriented_side (const unsigned int quadnum,
330  unsigned int &sidenum) const;
331 
350  bool side_hop (unsigned int &qnum,
351  unsigned int &lsn) const;
352 
357  std::vector<MSide> sides;
362  std::vector<MQuad> mquads;
363  };
364  } // namespace GridReordering2d
365 
366 
373  namespace GridReordering3d
374  {
427  {
432 
436  bool operator == (const EdgeOrientation &edge_orientation) const;
437 
441  bool operator != (const EdgeOrientation &edge_orientation) const;
442  };
443 
460  struct CheapEdge
461  {
465  const unsigned int node0;
466 
470  const unsigned int node1;
471 
477  CheapEdge (const unsigned int n0,
478  const unsigned int n1);
479 
484  bool operator< (const CheapEdge &e2) const;
485  };
486 
487 
488 
493  struct Edge
494  {
498  Edge (const unsigned int n0,
499  const unsigned int n1);
500 
504  unsigned int nodes[2];
505 
515 
530  unsigned int group;
531 
535  std::vector<unsigned int> neighboring_cubes;
536  };
537 
559  struct Cell
560  {
564  Cell ();
565 
570 
575 
585 
594  };
595 
596 
606  class Mesh
607  {
608  public:
612  Mesh (const std::vector<CellData<3> > &incubes);
613 
621  void
622  export_to_deal_format (std::vector<CellData<3> > &outcubes) const;
623 
624  private:
628  std::vector<Edge> edge_list;
629 
633  std::vector<Cell> cell_list;
634 
639  void sanity_check() const;
640 
648  void build_connectivity ();
649 
654  Mesh (const Mesh &);
655 
661  Mesh &operator=(const Mesh &);
662 
668  void sanity_check_node (const Cell &cell,
669  const unsigned int local_node_num) const;
670 
675  friend class Orienter;
676  };
677 
678 
687  class Orienter
688  {
689  public:
706  static
707  bool
708  orient_mesh (std::vector<CellData<3> > &incubes);
709 
710  private:
718 
723  unsigned int cur_posn;
724 
729  unsigned int marker_cube;
730 
736  unsigned int cur_edge_group;
737 
746  std::vector<int> sheet_to_process;
747 
748 
757 
771  Orienter (const std::vector<CellData<3> > &incubes);
772 
780  bool orient_edges ();
781 
788  void orient_cubes ();
789 
790  bool get_next_unoriented_cube ();
791 
798  bool is_oriented (const unsigned int cell_num) const;
799 
800  bool orient_edges_in_current_cube ();
801  bool orient_edge_set_in_current_cube (const unsigned int edge_set);
802  bool orient_next_unoriented_edge ();
803 
815  bool cell_is_consistent (const unsigned int cell_num) const;
816 
817 
818  void get_adjacent_cubes ();
819  bool get_next_active_cube ();
820  };
821  } // namespace GridReordering3d
822 } // namespace internal
823 
824 
825 DEAL_II_NAMESPACE_CLOSE
826 
827 #endif
bool get_unoriented_quad(unsigned int &UnOrQLoc) const
bool get_unoriented_side(const unsigned int quadnum, unsigned int &sidenum) const
bool cell_is_consistent(const unsigned int cell_num) const
Mesh & operator=(const Mesh &)
bool operator<(const CheapEdge &e2) const
bool is_oriented(const unsigned int cell_num) const
MQuad(const unsigned int v0, const unsigned int v1, const unsigned int v2, const unsigned int v3, const unsigned int s0, const unsigned int s1, const unsigned int s2, const unsigned int s3, const CellData< 2 > &cd)
bool operator==(const EdgeOrientation &edge_orientation) const
void reorient(std::vector< CellData< 2 > > &quads)
bool side_hop(unsigned int &qnum, unsigned int &lsn) const
unsigned int nodes[GeometryInfo< 3 >::vertices_per_cell]
bool is_side_default_oriented(const unsigned int qnum, const unsigned int lsn) const
void build_graph(const std::vector< CellData< 2 > > &inquads)
Orienter(const std::vector< CellData< 3 > > &incubes)
std::vector< unsigned int > neighboring_cubes
MSide(const unsigned int initv0, const unsigned int initv1)
bool operator!=(const MSide &s2) const
EdgeOrientation local_orientation_flags[GeometryInfo< 3 >::lines_per_cell]
CheapEdge(const unsigned int n0, const unsigned int n1)
static bool orient_mesh(std::vector< CellData< 3 > > &incubes)
unsigned int edges[GeometryInfo< 3 >::lines_per_cell]
bool operator!=(const EdgeOrientation &edge_orientation) const
Mesh(const std::vector< CellData< 3 > > &incubes)
bool operator==(const MSide &s2) const
bool is_consistent(const std::vector< CellData< 2 > > &cells)
bool is_fully_oriented_quad(const unsigned int quadnum) const
void sanity_check_node(const Cell &cell, const unsigned int local_node_num) const
void export_to_deal_format(std::vector< CellData< 3 > > &outcubes) const
Edge(const unsigned int n0, const unsigned int n1)
void orient_side(const unsigned int quadnum, const unsigned int localsidenum)
bool is_oriented_side(const unsigned int quadnum, const unsigned int lsn) const
void get_quads(std::vector< CellData< 2 > > &outquads) const