Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
grid_generator.h
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1999 - 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_grid_generator_h
17 #define dealii_grid_generator_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/exceptions.h>
23 #include <deal.II/base/function.h>
24 #include <deal.II/base/point.h>
25 #include <deal.II/base/table.h>
26 
27 #include <deal.II/grid/tria.h>
28 
29 #include <array>
30 #include <map>
31 
32 DEAL_II_NAMESPACE_OPEN
33 
56 namespace GridGenerator
57 {
61 
90  template <int dim, int spacedim>
91  void
93  const double left = 0.,
94  const double right = 1.,
95  const bool colorize = false);
96 
121  template <int dim>
122  void
124  const std::vector<Point<dim>> &vertices);
125 
139  template <int dim, int spacedim>
140  void
142  const unsigned int repetitions,
143  const double left = 0.,
144  const double right = 1.);
145 
180  template <int dim, int spacedim>
181  void
183  const Point<dim> & p1,
184  const Point<dim> & p2,
185  const bool colorize = false);
186 
237  template <int dim, int spacedim>
238  void
240  const std::vector<unsigned int> &repetitions,
241  const Point<dim> & p1,
242  const Point<dim> & p2,
243  const bool colorize = false);
244 
260  template <int dim>
261  void
263  const std::vector<std::vector<double>> &step_sizes,
264  const Point<dim> & p_1,
265  const Point<dim> & p_2,
266  const bool colorize = false);
267 
278  template <int dim>
279  void
281  const std::vector<std::vector<double>> &spacing,
282  const Point<dim> & p,
283  const Table<dim, types::material_id> &material_id,
284  const bool colorize = false);
285 
312  template <int dim, int spacedim>
313  void
315  const std::vector<unsigned int> &holes);
316 
360  template <int dim>
361  void
363  const double inner_radius = 0.4,
364  const double outer_radius = 1.,
365  const double pad_bottom = 2.,
366  const double pad_top = 2.,
367  const double pad_left = 1.,
368  const double pad_right = 1.,
369  const Point<dim> center = Point<dim>(),
370  const types::manifold_id polar_manifold_id = 0,
371  const types::manifold_id tfi_manifold_id = 1,
372  const double L = 1.,
373  const unsigned int n_slices = 2,
374  const bool colorize = false);
375 
462  template <int dim>
463  void
465  const double shell_region_width = 0.03,
466  const unsigned int n_shells = 2,
467  const double skewness = 2.0,
468  const bool colorize = false);
469 
492  template <int dim, int spacedim>
493  void
495  const std::vector<Point<spacedim>> &vertices,
496  const bool colorize = false);
497 
508  template <int dim>
509  void
511  const Point<dim> (&corners)[dim],
512  const bool colorize = false);
513 
528  template <int dim>
529  void
531  const Point<dim> (&corners)[dim],
532  const bool colorize = false);
533 
544  template <int dim>
545  void
547  const unsigned int n_subdivisions,
548  const Point<dim> (&corners)[dim],
549  const bool colorize = false);
550 
558  template <int dim>
559  void
561 #ifndef _MSC_VER
562  const unsigned int (&n_subdivisions)[dim],
563 #else
564  const unsigned int *n_subdivisions,
565 #endif
566  const Point<dim> (&corners)[dim],
567  const bool colorize = false);
568 
592  template <int dim, int spacedim>
593  void
595  const Point<spacedim> & origin,
596  const std::array<Tensor<1, spacedim>, dim> &edges,
597  const std::vector<unsigned int> &subdivisions = {},
598  const bool colorize = false);
599 
615  template <int dim>
616  void
618  const double left = 0.,
619  const double right = 1.,
620  const double thickness = 1.,
621  const bool colorize = false);
622 
658  template <int dim>
659  void
661  const Point<dim> & center = Point<dim>(),
662  const double radius = 1.,
663  const bool attach_spherical_manifold_on_boundary_cells = false);
664 
690  template <int spacedim>
692  const Point<spacedim> &center = Point<spacedim>(),
693  const double radius = 1.);
694 
706  template <int dim>
707  void
709  const Point<dim> & center = Point<dim>(),
710  const double radius = 1.);
711 
723  template <int dim>
724  void
726  const Point<dim> & center = Point<dim>(),
727  const double radius = 1.);
728 
756  template <int dim>
757  void
759  const double radius = 1.,
760  const double half_length = 1.);
761 
793  template <int dim>
794  void
796  const double radius_0 = 1.0,
797  const double radius_1 = 0.5,
798  const double half_length = 1.0);
799 
829  template <int dim, int spacedim>
830  void
832  const std::vector<unsigned int> &sizes,
833  const bool colorize_cells = false);
834 
869  template <int dim>
870  void
872  const double left = -1.,
873  const double right = 1.,
874  const bool colorize = false);
875 
894  template <int dim>
895  void
897  const double left = 0.,
898  const double right = 1.,
899  const bool colorize = false);
900 
939  template <int dim>
940  void
942  const Point<dim> & center,
943  const double inner_radius,
944  const double outer_radius,
945  const unsigned int n_cells = 0,
946  bool colorize = false);
947 
973  template <int dim>
974  void
976  const Point<dim> & center,
977  const double inner_radius,
978  const double outer_radius,
979  const unsigned int n_cells = 0,
980  const bool colorize = false);
981 
982 
1007  template <int dim>
1008  void
1010  const Point<dim> & center,
1011  const double inner_radius,
1012  const double outer_radius,
1013  const unsigned int n_cells = 0,
1014  const bool colorize = false);
1015 
1035  template <int dim>
1036  void
1038  const double length,
1039  const double inner_radius,
1040  const double outer_radius,
1041  const unsigned int n_radial_cells = 0,
1042  const unsigned int n_axial_cells = 0);
1043 
1044 
1045 
1067  template <int dim, int spacedim>
1068  void
1069  torus(Triangulation<dim, spacedim> &tria, const double R, const double r);
1070 
1100  template <int dim>
1101  void
1103  const double inner_radius = .25,
1104  const double outer_radius = .5,
1105  const double L = .5,
1106  const unsigned int repetitions = 1,
1107  const bool colorize = false);
1108 
1180  template <int dim>
1181  void
1183  const Point<dim> & center,
1184  const double inner_radius = 0.125,
1185  const double outer_radius = 0.25,
1186  const unsigned int n_shells = 1,
1187  const double skewness = 0.1,
1188  const unsigned int n_cells_per_shell = 0,
1189  const bool colorize = false);
1190 
1204  void moebius(Triangulation<3, 3> &tria,
1205  const unsigned int n_cells,
1206  const unsigned int n_rotations,
1207  const double R,
1208  const double r);
1209 
1211 
1215 
1276  template <int dim, int spacedim>
1277  void
1278  merge_triangulations(const Triangulation<dim, spacedim> &triangulation_1,
1279  const Triangulation<dim, spacedim> &triangulation_2,
1281  const double duplicated_vertex_tolerance = 1.0e-12,
1282  const bool copy_manifold_ids = false);
1283 
1298  template <int dim, int spacedim>
1299  void
1301  const std::initializer_list<const Triangulation<dim, spacedim> *const>
1302  & triangulations,
1304  const double duplicated_vertex_tolerance = 1.0e-12,
1305  const bool copy_manifold_ids = false);
1306 
1343  template <int dim, int spacedim>
1344  void
1346  const Triangulation<dim, spacedim> &triangulation_1,
1347  const Triangulation<dim, spacedim> &triangulation_2,
1348  Triangulation<dim, spacedim> & result);
1349 
1385  template <int dim, int spacedim>
1386  void
1388  const Triangulation<dim, spacedim> &input_triangulation,
1390  & cells_to_remove,
1392 
1447  void
1449  const Triangulation<2, 2> & input,
1450  const unsigned int n_slices,
1451  const double height,
1452  Triangulation<3, 3> & result,
1453  const bool copy_manifold_ids = false,
1454  const std::vector<types::manifold_id> &manifold_priorities = {});
1455 
1476  void
1478  const Triangulation<2, 2> & input,
1479  const std::vector<double> & slice_coordinates,
1480  Triangulation<3, 3> & result,
1481  const bool copy_manifold_ids = false,
1482  const std::vector<types::manifold_id> &manifold_priorities = {});
1483 
1484 
1515  template <int dim, int spacedim1, int spacedim2>
1516  void
1518  Triangulation<dim, spacedim2> & out_tria);
1519 
1521 
1527 
1529 #ifdef _MSC_VER
1530  // Microsoft's VC++ has a bug where it doesn't want to recognize that
1531  // an implementation (definition) of the extract_boundary_mesh function
1532  // matches a declaration. This can apparently only be avoided by
1533  // doing some contortion with the return type using the following
1534  // intermediate type. This is only used when using MS VC++ and uses
1535  // the direct way of doing it otherwise
1536  template <template <int, int> class MeshType, int dim, int spacedim>
1537  struct ExtractBoundaryMesh
1538  {
1539  using return_type =
1540  std::map<typename MeshType<dim - 1, spacedim>::cell_iterator,
1541  typename MeshType<dim, spacedim>::face_iterator>;
1542  };
1543 #endif
1544 
1620  template <template <int, int> class MeshType, int dim, int spacedim>
1621 #ifndef _MSC_VER
1622  std::map<typename MeshType<dim - 1, spacedim>::cell_iterator,
1623  typename MeshType<dim, spacedim>::face_iterator>
1624 #else
1625  typename ExtractBoundaryMesh<MeshType, dim, spacedim>::return_type
1626 #endif
1627  extract_boundary_mesh(const MeshType<dim, spacedim> & volume_mesh,
1628  MeshType<dim - 1, spacedim> & surface_mesh,
1629  const std::set<types::boundary_id> &boundary_ids =
1630  std::set<types::boundary_id>());
1631 
1633 
1634 
1638 
1640 
1649  int,
1650  << "The number of repetitions " << arg1 << " must be >=1.");
1655  int,
1656  << "The vector of repetitions must have " << arg1
1657  << " elements.");
1658 
1663  "The input to this function is oriented in a way that will"
1664  " cause all cells to have negative measure.");
1666 
1667 #ifndef DOXYGEN
1668  // These functions are only implemented with specializations; declare them
1669  // here
1670  template <>
1672  const double,
1673  const double,
1674  const double,
1675  const unsigned int,
1676  const bool);
1677 
1678  template <>
1680  const double,
1681  const double,
1682  const double,
1683  const unsigned int,
1684  const bool);
1685 
1686  template <>
1688  const double,
1689  const double,
1690  const double,
1691  const unsigned int,
1692  const bool);
1693 
1694  template <>
1696  const double,
1697  const unsigned int,
1698  const double,
1699  const bool);
1700 
1701  template <>
1703  const double,
1704  const unsigned int,
1705  const double,
1706  const bool);
1707 
1708  template <>
1710  const double,
1711  const unsigned int,
1712  const double,
1713  const bool);
1714 #endif
1715 } // namespace GridGenerator
1716 
1717 
1718 
1719 DEAL_II_NAMESPACE_CLOSE
1720 
1721 #endif
void create_union_triangulation(const Triangulation< dim, spacedim > &triangulation_1, const Triangulation< dim, spacedim > &triangulation_2, Triangulation< dim, spacedim > &result)
unsigned int manifold_id
Definition: types.h:123
void torus(Triangulation< dim, spacedim > &tria, const double R, const double r)
void simplex(Triangulation< dim, dim > &tria, const std::vector< Point< dim >> &vertices)
Triangulation of a d-simplex with (d+1) vertices and mesh cells.
void channel_with_cylinder(Triangulation< dim > &tria, const double shell_region_width=0.03, const unsigned int n_shells=2, const double skewness=2.0, const bool colorize=false)
void concentric_hyper_shells(Triangulation< dim > &triangulation, const Point< dim > &center, const double inner_radius=0.125, const double outer_radius=0.25, const unsigned int n_shells=1, const double skewness=0.1, const unsigned int n_cells_per_shell=0, const bool colorize=false)
void general_cell(Triangulation< dim, spacedim > &tria, const std::vector< Point< spacedim >> &vertices, const bool colorize=false)
void truncated_cone(Triangulation< dim > &tria, const double radius_0=1.0, const double radius_1=0.5, const double half_length=1.0)
void hyper_rectangle(Triangulation< dim, spacedim > &tria, const Point< dim > &p1, const Point< dim > &p2, const bool colorize=false)
void moebius(Triangulation< 3, 3 > &tria, const unsigned int n_cells, const unsigned int n_rotations, const double R, const double r)
void flatten_triangulation(const Triangulation< dim, spacedim1 > &in_tria, Triangulation< dim, spacedim2 > &out_tria)
void quarter_hyper_shell(Triangulation< dim > &tria, const Point< dim > &center, const double inner_radius, const double outer_radius, const unsigned int n_cells=0, const bool colorize=false)
void enclosed_hyper_cube(Triangulation< dim > &tria, const double left=0., const double right=1., const double thickness=1., const bool colorize=false)
void subdivided_parallelepiped(Triangulation< dim > &tria, const unsigned int n_subdivisions, const Point< dim >(&corners)[dim], const bool colorize=false)
static ::ExceptionBase & ExcInvalidRepetitionsDimension(int arg1)
void parallelogram(Triangulation< dim > &tria, const Point< dim >(&corners)[dim], const bool colorize=false)
void hyper_shell(Triangulation< dim > &tria, const Point< dim > &center, const double inner_radius, const double outer_radius, const unsigned int n_cells=0, bool colorize=false)
static ::ExceptionBase & ExcInvalidInputOrientation()
void cylinder_shell(Triangulation< dim > &tria, const double length, const double inner_radius, const double outer_radius, const unsigned int n_radial_cells=0, const unsigned int n_axial_cells=0)
void hyper_cross(Triangulation< dim, spacedim > &tria, const std::vector< unsigned int > &sizes, const bool colorize_cells=false)
A center cell with stacks of cell protruding from each surface.
void hyper_L(Triangulation< dim > &tria, const double left=-1., const double right=1., const bool colorize=false)
void half_hyper_ball(Triangulation< dim > &tria, const Point< dim > &center=Point< dim >(), const double radius=1.)
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
void hyper_sphere(Triangulation< spacedim - 1, spacedim > &tria, const Point< spacedim > &center=Point< spacedim >(), const double radius=1.)
void cylinder(Triangulation< dim > &tria, const double radius=1., const double half_length=1.)
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:496
#define DeclException0(Exception0)
Definition: exceptions.h:473
void quarter_hyper_ball(Triangulation< dim > &tria, const Point< dim > &center=Point< dim >(), const double radius=1.)
void subdivided_hyper_cube(Triangulation< dim, spacedim > &tria, const unsigned int repetitions, const double left=0., const double right=1.)
void create_triangulation_with_removed_cells(const Triangulation< dim, spacedim > &input_triangulation, const std::set< typename Triangulation< dim, spacedim >::active_cell_iterator > &cells_to_remove, Triangulation< dim, spacedim > &result)
void extrude_triangulation(const Triangulation< 2, 2 > &input, const unsigned int n_slices, const double height, Triangulation< 3, 3 > &result, const bool copy_manifold_ids=false, const std::vector< types::manifold_id > &manifold_priorities={})
void hyper_cube_slit(Triangulation< dim > &tria, const double left=0., const double right=1., const bool colorize=false)
void subdivided_hyper_rectangle(Triangulation< dim, spacedim > &tria, const std::vector< unsigned int > &repetitions, const Point< dim > &p1, const Point< dim > &p2, const bool colorize=false)
void merge_triangulations(const Triangulation< dim, spacedim > &triangulation_1, const Triangulation< dim, spacedim > &triangulation_2, Triangulation< dim, spacedim > &result, const double duplicated_vertex_tolerance=1.0e-12, const bool copy_manifold_ids=false)
static ::ExceptionBase & ExcInvalidRadii()
void cheese(Triangulation< dim, spacedim > &tria, const std::vector< unsigned int > &holes)
Rectangular domain with rectangular pattern of holes.
__global__ void set(Number *val, const Number s, const size_type N)
void parallelepiped(Triangulation< dim > &tria, const Point< dim >(&corners)[dim], const bool colorize=false)
void half_hyper_shell(Triangulation< dim > &tria, const Point< dim > &center, const double inner_radius, const double outer_radius, const unsigned int n_cells=0, const bool colorize=false)
void hyper_cube(Triangulation< dim, spacedim > &tria, const double left=0., const double right=1., const bool colorize=false)
std::map< typename MeshType< dim - 1, spacedim >::cell_iterator, typename MeshType< dim, spacedim >::face_iterator > extract_boundary_mesh(const MeshType< dim, spacedim > &volume_mesh, MeshType< dim - 1, spacedim > &surface_mesh, const std::set< types::boundary_id > &boundary_ids=std::set< types::boundary_id >())
void hyper_ball(Triangulation< dim > &tria, const Point< dim > &center=Point< dim >(), const double radius=1., const bool attach_spherical_manifold_on_boundary_cells=false)
Definition: table.h:37
void plate_with_a_hole(Triangulation< dim > &tria, const double inner_radius=0.4, const double outer_radius=1., const double pad_bottom=2., const double pad_top=2., const double pad_left=1., const double pad_right=1., const Point< dim > center=Point< dim >(), const types::manifold_id polar_manifold_id=0, const types::manifold_id tfi_manifold_id=1, const double L=1., const unsigned int n_slices=2, const bool colorize=false)
Rectangular plate with an (offset) cylindrical hole.
void hyper_cube_with_cylindrical_hole(Triangulation< dim > &triangulation, const double inner_radius=.25, const double outer_radius=.5, const double L=.5, const unsigned int repetitions=1, const bool colorize=false)
static ::ExceptionBase & ExcInvalidRepetitions(int arg1)