Reference documentation for deal.II version 8.1.0
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
grid_tools.h
1 // ---------------------------------------------------------------------
2 // @f$Id: grid_tools.h 31940 2013-12-08 15:49:12Z heister @f$
3 //
4 // Copyright (C) 2001 - 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_tools_H
18 #define __deal2__grid_tools_H
19 
20 
21 #include <deal.II/base/config.h>
22 #include <deal.II/fe/mapping.h>
23 #include <deal.II/grid/tria.h>
24 #include <deal.II/grid/tria_accessor.h>
25 #include <deal.II/grid/tria_iterator.h>
26 #include <deal.II/fe/mapping_q1.h>
27 
28 #include <bitset>
29 #include <list>
30 
31 DEAL_II_NAMESPACE_OPEN
32 
33 
34 template <int, int> class DoFHandler;
35 template <int, int> class Mapping;
36 namespace hp
37 {
38  template <int, int> class DoFHandler;
39  template <int, int> class MappingCollection;
40 }
41 
42 class SparsityPattern;
43 
44 
53 namespace GridTools
54 {
67  template <int dim, int spacedim>
68  double diameter (const Triangulation<dim, spacedim> &tria);
69 
80  template <int dim, int spacedim>
81  double volume (const Triangulation<dim,spacedim> &tria,
83 
100  template <int dim>
101  double cell_measure (const std::vector<Point<dim> > &all_vertices,
102  const unsigned int (&vertex_indices)[GeometryInfo<dim>::vertices_per_cell]);
103 
136  template <int dim, int spacedim>
137  void delete_unused_vertices (std::vector<Point<spacedim> > &vertices,
138  std::vector<CellData<dim> > &cells,
139  SubCellData &subcelldata);
140 
164  template <int dim, int spacedim>
165  void delete_duplicated_vertices (std::vector<Point<spacedim> > &all_vertices,
166  std::vector<CellData<dim> > &cells,
167  SubCellData &subcelldata,
168  std::vector<unsigned int> &considered_vertices,
169  const double tol=1e-12);
170 
201  template <int dim, typename Transformation, int spacedim>
202  void transform (const Transformation &transformation,
203  Triangulation<dim,spacedim> &triangulation);
204 
216  template <int dim, int spacedim>
217  void shift (const Point<spacedim> &shift_vector,
218  Triangulation<dim,spacedim> &triangulation);
219 
220 
236  void rotate (const double angle,
237  Triangulation<2> &triangulation);
238 
267  template <int dim>
268  void laplace_transform (const std::map<unsigned int,Point<dim> > &new_points,
269  Triangulation<dim> &tria);
270 
285  template <int dim, int spacedim>
286  void scale (const double scaling_factor,
287  Triangulation<dim, spacedim> &triangulation);
288 
307  template <int dim, int spacedim>
308  void distort_random (const double factor,
309  Triangulation<dim, spacedim> &triangulation,
310  const bool keep_boundary=true);
311 
324  template <int dim, template <int, int> class Container, int spacedim>
325  unsigned int
326  find_closest_vertex (const Container<dim, spacedim> &container,
327  const Point<spacedim> &p);
328 
355  template<int dim, template <int, int> class Container, int spacedim>
356  std::vector<typename Container<dim,spacedim>::active_cell_iterator>
357  find_cells_adjacent_to_vertex (const Container<dim,spacedim> &container,
358  const unsigned int vertex_index);
359 
360 
388  template <int dim, template <int,int> class Container, int spacedim>
389  typename Container<dim,spacedim>::active_cell_iterator
390  find_active_cell_around_point (const Container<dim,spacedim> &container,
391  const Point<spacedim> &p);
392 
436  template <int dim, template<int, int> class Container, int spacedim>
437  std::pair<typename Container<dim,spacedim>::active_cell_iterator, Point<dim> >
439  const Container<dim,spacedim> &container,
440  const Point<spacedim> &p);
441 
462  template <int dim, int spacedim>
463  std::pair<typename hp::DoFHandler<dim,spacedim>::active_cell_iterator, Point<dim> >
465  const hp::DoFHandler<dim,spacedim> &container,
466  const Point<spacedim> &p);
467 
491  template <class Container>
492  std::vector<typename Container::active_cell_iterator>
493  get_active_child_cells (const typename Container::cell_iterator &cell);
494 
500  template <class Container>
501  void
502  get_active_neighbors (const typename Container::active_cell_iterator &cell,
503  std::vector<typename Container::active_cell_iterator> &active_neighbors);
504 
516  template <int dim, int spacedim>
517  void
519  SparsityPattern &connectivity);
520 
541  template <int dim, int spacedim>
542  void
543  partition_triangulation (const unsigned int n_partitions,
544  Triangulation<dim, spacedim> &triangulation);
545 
622  template <int dim, int spacedim>
623  void
624  partition_triangulation (const unsigned int n_partitions,
625  const SparsityPattern &cell_connection_graph,
626  Triangulation<dim,spacedim> &triangulation);
627 
643  template <int dim, int spacedim>
644  void
646  std::vector<types::subdomain_id> &subdomain);
647 
670  template <int dim, int spacedim>
671  unsigned int
673  const types::subdomain_id subdomain);
674 
731  template <typename Container>
732  std::list<std::pair<typename Container::cell_iterator,
733  typename Container::cell_iterator> >
734  get_finest_common_cells (const Container &mesh_1,
735  const Container &mesh_2);
736 
751  template <int dim, int spacedim>
752  bool
754  const Triangulation<dim, spacedim> &mesh_2);
755 
769  template <typename Container>
770  bool
771  have_same_coarse_mesh (const Container &mesh_1,
772  const Container &mesh_2);
773 
780  template <int dim, int spacedim>
781  double
783 
788  template <int dim, int spacedim>
789  double
791 
841  template <int dim, int spacedim>
842  void
844  const Triangulation<dim, spacedim> &triangulation_2,
846 
872  template <int dim, int spacedim>
875  Triangulation<dim,spacedim> &triangulation);
876 
963  template <template <int,int> class Container, int dim, int spacedim>
964  std::map<typename Container<dim-1,spacedim>::cell_iterator,
965  typename Container<dim,spacedim>::face_iterator>
966  extract_boundary_mesh (const Container<dim,spacedim> &volume_mesh,
967  Container<dim-1,spacedim> &surface_mesh,
968  const std::set<types::boundary_id> &boundary_ids
969  = std::set<types::boundary_id>());
970 
976  template<typename CellIterator>
977  struct PeriodicFacePair
978  {
979  CellIterator cell[2];
980  unsigned int face_idx[2];
981  std::bitset<3> orientation;
982  };
983 
1047  template<typename FaceIterator>
1048  bool
1049  orthogonal_equality (std::bitset<3> &orientation,
1050  const FaceIterator &face1,
1051  const FaceIterator &face2,
1052  const int direction,
1053  const ::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset);
1054 
1055 
1059  template<typename FaceIterator>
1060  bool
1061  orthogonal_equality (const FaceIterator &face1,
1062  const FaceIterator &face2,
1063  const int direction,
1064  const ::Tensor<1,FaceIterator::AccessorType::space_dimension> &offset);
1065 
1066 
1105  template<typename CONTAINER>
1106  void
1108  (const CONTAINER &container,
1109  const types::boundary_id b_id1,
1110  const types::boundary_id b_id2,
1111  const int direction,
1112  std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> > &matched_pairs,
1113  const ::Tensor<1,CONTAINER::space_dimension> &offset = ::Tensor<1,CONTAINER::space_dimension>());
1114 
1115 
1136  template<typename CONTAINER>
1137  void
1139  (const CONTAINER &container,
1140  const types::boundary_id b_id,
1141  const int direction,
1142  std::vector<PeriodicFacePair<typename CONTAINER::cell_iterator> > &matched_pairs,
1143  const ::Tensor<1,CONTAINER::space_dimension> &offset = ::Tensor<1,CONTAINER::space_dimension>());
1144 
1145 
1149  DeclException1 (ExcInvalidNumberOfPartitions,
1150  int,
1151  << "The number of partitions you gave is " << arg1
1152  << ", but must be greater than zero.");
1156  DeclException1 (ExcNonExistentSubdomain,
1157  int,
1158  << "The subdomain id " << arg1
1159  << " has no cells associated with it.");
1163  DeclException0 (ExcTriangulationHasBeenRefined);
1167  DeclException1 (ExcScalingFactorNotPositive,
1168  double,
1169  << "The scaling factor must be positive, but is " << arg1);
1173  template <int N>
1174  DeclException1 (ExcPointNotFoundInCoarseGrid,
1175  Point<N>,
1176  << "The point <" << arg1
1177  << "> could not be found inside any of the "
1178  << "coarse grid cells.");
1182  template <int N>
1183  DeclException1 (ExcPointNotFound,
1184  Point<N>,
1185  << "The point <" << arg1
1186  << "> could not be found inside any of the "
1187  << "subcells of a coarse grid cell.");
1188 
1189  DeclException1 (ExcVertexNotUsed,
1190  unsigned int,
1191  << "The given vertex " << arg1
1192  << " is not used in the given triangulation");
1193 
1194 
1195 } /*namespace GridTools*/
1196 
1197 
1198 
1199 /* ----------------- Template function --------------- */
1200 
1201 namespace GridTools
1202 {
1203  template <int dim, typename Predicate, int spacedim>
1204  void transform (const Predicate &predicate,
1205  Triangulation<dim, spacedim> &triangulation)
1206  {
1207  // ensure that all the cells of the
1208  // triangulation are on the same level
1209  Assert (triangulation.n_levels() ==
1210  static_cast<unsigned int>(triangulation.begin_active()->level()+1),
1211  ExcMessage ("Not all cells of this triangulation are at the same "
1212  "refinement level, as is required for this function."));
1213 
1214  std::vector<bool> treated_vertices (triangulation.n_vertices(),
1215  false);
1216 
1217  // loop over all active cells, and
1218  // transform those vertices that
1219  // have not yet been touched. note
1220  // that we get to all vertices in
1221  // the triangulation by only
1222  // visiting the active cells.
1224  cell = triangulation.begin_active (),
1225  endc = triangulation.end ();
1226  for (; cell!=endc; ++cell)
1227  for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
1228  if (treated_vertices[cell->vertex_index(v)] == false)
1229  {
1230  // transform this vertex
1231  cell->vertex(v) = predicate(cell->vertex(v));
1232  // and mark it as treated
1233  treated_vertices[cell->vertex_index(v)] = true;
1234  };
1235  }
1236 
1237 
1238 
1239  template <class DH>
1240  std::vector<typename DH::active_cell_iterator>
1241  get_active_child_cells (const typename DH::cell_iterator &cell)
1242  {
1243  std::vector<typename DH::active_cell_iterator> child_cells;
1244 
1245  if (cell->has_children())
1246  {
1247  for (unsigned int child=0;
1248  child<cell->n_children(); ++child)
1249  if (cell->child (child)->has_children())
1250  {
1251  const std::vector<typename DH::active_cell_iterator>
1252  children = get_active_child_cells<DH> (cell->child(child));
1253  child_cells.insert (child_cells.end(),
1254  children.begin(), children.end());
1255  }
1256  else
1257  child_cells.push_back (cell->child(child));
1258  }
1259 
1260  return child_cells;
1261  }
1262 
1263 
1264 
1265  template <class Container>
1266  void
1267  get_active_neighbors(const typename Container::active_cell_iterator &cell,
1268  std::vector<typename Container::active_cell_iterator> &active_neighbors)
1269  {
1270  active_neighbors.clear ();
1271  for (unsigned int n=0; n<GeometryInfo<Container::dimension>::faces_per_cell; ++n)
1272  if (! cell->at_boundary(n))
1273  {
1274  if (Container::dimension == 1)
1275  {
1276  // check children of neighbor. note
1277  // that in 1d children of the neighbor
1278  // may be further refined. In 1d the
1279  // case is simple since we know what
1280  // children bound to the present cell
1281  typename Container::cell_iterator
1282  neighbor_child = cell->neighbor(n);
1283  if (!neighbor_child->active())
1284  {
1285  while (neighbor_child->has_children())
1286  neighbor_child = neighbor_child->child (n==0 ? 1 : 0);
1287 
1288  Assert (neighbor_child->neighbor(n==0 ? 1 : 0)==cell,
1289  ExcInternalError());
1290  }
1291  active_neighbors.push_back (neighbor_child);
1292  }
1293  else
1294  {
1295  if (cell->face(n)->has_children())
1296  // this neighbor has children. find
1297  // out which border to the present
1298  // cell
1299  for (unsigned int c=0; c<cell->face(n)->number_of_children(); ++c)
1300  active_neighbors.push_back (cell->neighbor_child_on_subface(n,c));
1301  else
1302  {
1303  // the neighbor must be active
1304  // himself
1305  Assert(cell->neighbor(n)->active(), ExcInternalError());
1306  active_neighbors.push_back(cell->neighbor(n));
1307  }
1308  }
1309  }
1310  }
1311 
1312 
1313 
1314 
1315 // declaration of explicit specializations
1316 
1317  template <>
1318  double
1319  cell_measure<3>(const std::vector<Point<3> > &all_vertices,
1320  const unsigned int (&vertex_indices) [GeometryInfo<3>::vertices_per_cell]);
1321 
1322  template <>
1323  double
1324  cell_measure<2>(const std::vector<Point<2> > &all_vertices,
1325  const unsigned int (&vertex_indices) [GeometryInfo<2>::vertices_per_cell]);
1326 }
1327 
1328 
1329 
1330 DEAL_II_NAMESPACE_CLOSE
1331 
1332 /*---------------------------- grid_tools.h ---------------------------*/
1333 /* end of #ifndef __deal2__grid_tools_H */
1334 #endif
1335 /*---------------------------- grid_tools.h ---------------------------*/
void collect_periodic_faces(const CONTAINER &container, const types::boundary_id b_id1, const types::boundary_id b_id2, const int direction, std::vector< PeriodicFacePair< typename CONTAINER::cell_iterator > > &matched_pairs, const ::Tensor< 1, CONTAINER::space_dimension > &offset=::Tensor< 1, CONTAINER::space_dimension >())
std::map< typename Container< dim-1, spacedim >::cell_iterator, typename Container< dim, spacedim >::face_iterator > extract_boundary_mesh(const Container< dim, spacedim > &volume_mesh, Container< dim-1, spacedim > &surface_mesh, const std::set< types::boundary_id > &boundary_ids=std::set< types::boundary_id >())
active_cell_iterator begin_active(const unsigned int level=0) const
double diameter(const Triangulation< dim, spacedim > &tria)
void distort_random(const double factor, Triangulation< dim, spacedim > &triangulation, const bool keep_boundary=true)
void delete_duplicated_vertices(std::vector< Point< spacedim > > &all_vertices, std::vector< CellData< dim > > &cells, SubCellData &subcelldata, std::vector< unsigned int > &considered_vertices, const double tol=1e-12)
double maximal_cell_diameter(const Triangulation< dim, spacedim > &triangulation)
bool have_same_coarse_mesh(const Triangulation< dim, spacedim > &mesh_1, const Triangulation< dim, spacedim > &mesh_2)
void scale(const double scaling_factor, Triangulation< dim, spacedim > &triangulation)
double volume(const Triangulation< dim, spacedim > &tria, const Mapping< dim, spacedim > &mapping=(StaticMappingQ1< dim, spacedim >::mapping))
::ExceptionBase & ExcMessage(std::string arg1)
void get_active_neighbors(const typename Container::active_cell_iterator &cell, std::vector< typename Container::active_cell_iterator > &active_neighbors)
Definition: grid_tools.h:1267
std::vector< typename Container< dim, spacedim >::active_cell_iterator > find_cells_adjacent_to_vertex(const Container< dim, spacedim > &container, const unsigned int vertex_index)
unsigned int find_closest_vertex(const Container< dim, spacedim > &container, const Point< spacedim > &p)
unsigned int n_levels() const
cell_iterator end() const
void create_union_triangulation(const Triangulation< dim, spacedim > &triangulation_1, const Triangulation< dim, spacedim > &triangulation_2, Triangulation< dim, spacedim > &result)
Triangulation< dim, spacedim >::DistortedCellList fix_up_distorted_child_cells(const typename Triangulation< dim, spacedim >::DistortedCellList &distorted_cells, Triangulation< dim, spacedim > &triangulation)
double cell_measure(const std::vector< Point< dim > > &all_vertices, const unsigned int(&vertex_indices)[GeometryInfo< dim >::vertices_per_cell])
bool orthogonal_equality(std::bitset< 3 > &orientation, const FaceIterator &face1, const FaceIterator &face2, const int direction, const ::Tensor< 1, FaceIterator::AccessorType::space_dimension > &offset)
#define Assert(cond, exc)
Definition: exceptions.h:299
std::list< std::pair< typename Container::cell_iterator, typename Container::cell_iterator > > get_finest_common_cells(const Container &mesh_1, const Container &mesh_2)
DeclException0(ExcTriangulationHasBeenRefined)
void shift(const Point< spacedim > &shift_vector, Triangulation< dim, spacedim > &triangulation)
void laplace_transform(const std::map< unsigned int, Point< dim > > &new_points, Triangulation< dim > &tria)
void rotate(const double angle, Triangulation< 2 > &triangulation)
void partition_triangulation(const unsigned int n_partitions, Triangulation< dim, spacedim > &triangulation)
unsigned int subdomain_id
Definition: types.h:43
double minimal_cell_diameter(const Triangulation< dim, spacedim > &triangulation)
std::vector< typename Container::active_cell_iterator > get_active_child_cells(const typename Container::cell_iterator &cell)
void transform(const Transformation &transformation, Triangulation< dim, spacedim > &triangulation)
void get_subdomain_association(const Triangulation< dim, spacedim > &triangulation, std::vector< types::subdomain_id > &subdomain)
void delete_unused_vertices(std::vector< Point< spacedim > > &vertices, std::vector< CellData< dim > > &cells, SubCellData &subcelldata)
unsigned char boundary_id
Definition: types.h:128
void get_face_connectivity_of_cells(const Triangulation< dim, spacedim > &triangulation, SparsityPattern &connectivity)
DeclException1(ExcInvalidNumberOfPartitions, int,<< "The number of partitions you gave is "<< arg1<< ", but must be greater than zero.")
Container< dim, spacedim >::active_cell_iterator find_active_cell_around_point(const Container< dim, spacedim > &container, const Point< spacedim > &p)
::ExceptionBase & ExcInternalError()
unsigned int n_vertices() const
unsigned int count_cells_with_subdomain_association(const Triangulation< dim, spacedim > &triangulation, const types::subdomain_id subdomain)