deal.II version GIT relicensing-3132-g78fd2c863f 2025-04-24 02:10:00+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
grid_generator.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1999 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_grid_generator_h
16#define dealii_grid_generator_h
17
18
19#include <deal.II/base/config.h>
20
23#include <deal.II/base/point.h>
24#include <deal.II/base/table.h>
25
27#include <deal.II/grid/tria.h>
28
30
31#include <array>
32#include <limits>
33#include <map>
34#include <set>
35
36
38
39#ifndef DOXYGEN
41#endif
42
66{
104 template <int dim, int spacedim>
105 void
107 const double left = 0.,
108 const double right = 1.,
109 const bool colorize = false);
110
145 template <int dim>
146 void
148 const std::vector<Point<dim>> &vertices);
149
155 template <int dim, int spacedim>
156 void
159
160
190 template <int dim, int spacedim>
191 void
193 const unsigned int repetitions,
194 const double left = 0.,
195 const double right = 1.,
196 const bool colorize = false);
197
236 template <int dim, int spacedim>
237 void
239 const Point<dim> &p1,
240 const Point<dim> &p2,
241 const bool colorize = false);
242
296 template <int dim, int spacedim>
297 void
299 const std::vector<unsigned int> &repetitions,
300 const Point<dim> &p1,
301 const Point<dim> &p2,
302 const bool colorize = false);
303
323 template <int dim>
324 void
326 const std::vector<std::vector<double>> &step_sizes,
327 const Point<dim> &p_1,
328 const Point<dim> &p_2,
329 const bool colorize = false);
330
345 template <int dim>
346 void
348 const std::vector<std::vector<double>> &spacing,
349 const Point<dim> &p,
350 const Table<dim, types::material_id> &material_id,
351 const bool colorize = false);
352
380 template <int dim, int spacedim>
381 void
383 const std::vector<unsigned int> &holes);
384
434 template <int dim>
435 void
437 const double inner_radius = 0.4,
438 const double outer_radius = 1.,
439 const double pad_bottom = 2.,
440 const double pad_top = 2.,
441 const double pad_left = 1.,
442 const double pad_right = 1.,
443 const Point<dim> &center = Point<dim>(),
444 const types::manifold_id polar_manifold_id = 0,
445 const types::manifold_id tfi_manifold_id = 1,
446 const double L = 1.,
447 const unsigned int n_slices = 2,
448 const bool colorize = false);
449
542 template <int dim>
543 void
545 const double shell_region_width = 0.03,
546 const unsigned int n_shells = 2,
547 const double skewness = 2.0,
548 const bool colorize = false);
549
573 template <int dim, int spacedim>
574 void
576 const std::vector<Point<spacedim>> &vertices,
577 const bool colorize = false);
578
598 template <int dim>
599 void
601 const Point<dim> (&corners)[dim],
602 const bool colorize = false);
603
622 template <int dim>
623 void
625 const Point<dim> (&corners)[dim],
626 const bool colorize = false);
627
643 template <int dim>
644 void
646 const unsigned int n_subdivisions,
647 const Point<dim> (&corners)[dim],
648 const bool colorize = false);
649
662 template <int dim>
663 void
665#ifndef _MSC_VER
666 const unsigned int (&n_subdivisions)[dim],
667#else
668 const unsigned int *n_subdivisions,
669#endif
670 const Point<dim> (&corners)[dim],
671 const bool colorize = false);
672
696 template <int dim, int spacedim>
697 void
699 const Point<spacedim> &origin,
700 const std::array<Tensor<1, spacedim>, dim> &edges,
701 const std::vector<unsigned int> &subdivisions = {},
702 const bool colorize = false);
703
724 template <int dim>
725 void
727 const double left = 0.,
728 const double right = 1.,
729 const double thickness = 1.,
730 const bool colorize = false);
731
814 template <int dim, int spacedim>
815 void
817 const Point<spacedim> &center = {},
818 const double radius = 1.,
819 const bool attach_spherical_manifold_on_boundary_cells = false);
820
856 template <int dim>
857 void
859 const Point<dim> &center = Point<dim>(),
860 const double radius = 1.);
861
880 void
882 const unsigned int n_rotate_middle_square);
883
905 void
907 const bool face_orientation,
908 const bool face_flip,
909 const bool face_rotation,
910 const bool manipulate_left_cube);
911
912
941 template <int spacedim>
942 void
944 const Point<spacedim> &center = Point<spacedim>(),
945 const double radius = 1.);
946
983 template <int dim>
984 void
986 const Point<dim> &center = Point<dim>(),
987 const double radius = 1.);
988
1004 template <int dim>
1005 void
1007 const Point<dim> &center = Point<dim>(),
1008 const double radius = 1.);
1009
1036 template <int dim>
1037 void
1039 const double radius = 1.,
1040 const double half_length = 1.);
1041
1042
1080 template <int dim>
1081 void
1083 const unsigned int x_subdivisions,
1084 const double radius = 1.,
1085 const double half_length = 1.);
1086
1087
1122 template <int dim>
1123 void
1125 const double radius_0 = 1.0,
1126 const double radius_1 = 0.5,
1127 const double half_length = 1.0);
1128
1252 template <int dim, int spacedim>
1253 void
1255 const std::vector<std::pair<Point<spacedim>, double>> &openings,
1256 const std::pair<Point<spacedim>, double> &bifurcation,
1257 const double aspect_ratio = 0.5);
1258
1285 template <int dim, int spacedim>
1286 void
1288 const std::vector<unsigned int> &sizes,
1289 const bool colorize_cells = false);
1290
1326 template <int dim>
1327 void
1329 const double left = -1.,
1330 const double right = 1.,
1331 const bool colorize = false);
1332
1364 template <int dim, int spacedim>
1365 void
1367 const std::vector<unsigned int> &repetitions,
1368 const Point<dim> &bottom_left,
1369 const Point<dim> &top_right,
1370 const std::vector<int> &n_cells_to_remove);
1371
1391 template <int dim>
1392 void
1394 const double left = 0.,
1395 const double right = 1.,
1396 const bool colorize = false);
1397
1459 template <int dim, int spacedim>
1460 void
1462 const Point<spacedim> &center,
1463 const double inner_radius,
1464 const double outer_radius,
1465 const unsigned int n_cells = 0,
1466 bool colorize = false);
1467
1501 template <int dim>
1502 void
1504 const Point<dim> &inner_center,
1505 const Point<dim> &outer_center,
1506 const double inner_radius,
1507 const double outer_radius,
1508 const unsigned int n_cells);
1509
1539 template <int dim>
1540 void
1542 const Point<dim> &center,
1543 const double inner_radius,
1544 const double outer_radius,
1545 const unsigned int n_cells = 0,
1546 const bool colorize = false);
1547
1548
1574 template <int dim>
1575 void
1577 const Point<dim> &center,
1578 const double inner_radius,
1579 const double outer_radius,
1580 const unsigned int n_cells = 0,
1581 const bool colorize = false);
1582
1614 template <int dim>
1615 void
1617 const double length,
1618 const double inner_radius,
1619 const double outer_radius,
1620 const unsigned int n_radial_cells = 0,
1621 const unsigned int n_axial_cells = 0,
1622 const bool colorize = false);
1623
1675 template <int dim, int spacedim>
1676 void
1678 const double centerline_radius,
1679 const double inner_radius,
1680 const unsigned int n_cells_toroidal = 6,
1681 const double phi = 2.0 * numbers::PI);
1682
1716 template <int dim, int spacedim>
1717 void
1719 const double inner_radius = .25,
1720 const double outer_radius = .5,
1721 const double L = .5,
1722 const unsigned int repetitions = 1,
1723 const bool colorize = false);
1724
1795 template <int dim>
1796 void
1798 const Point<dim> &center,
1799 const double inner_radius = 0.125,
1800 const double outer_radius = 0.25,
1801 const unsigned int n_shells = 1,
1802 const double skewness = 0.1,
1803 const unsigned int n_cells_per_shell = 0,
1804 const bool colorize = false);
1805
1819 void
1821 const unsigned int n_cells,
1822 const unsigned int n_rotations,
1823 const double R,
1824 const double r);
1825
1861 template <int dim, int spacedim>
1862 void
1865 const std::string &grid_generator_function_name,
1866 const std::string &grid_generator_function_arguments);
1867
1918 template <int dim>
1919 void
1924 const Point<3> &interior_point = Point<3>(),
1925 const double &outer_ball_radius = 1.0);
1926
1941 void
1943 Triangulation<3> &vol_tria,
2014 template <int dim, int spacedim>
2015 void
2017 const Triangulation<dim, spacedim> &triangulation_2,
2019 const double duplicated_vertex_tolerance = 1.0e-12,
2020 const bool copy_manifold_ids = false,
2021 const bool copy_boundary_ids = false);
2022
2042 template <int dim, int spacedim>
2043 void
2045 const std::vector<const Triangulation<dim, spacedim> *> &triangulations,
2047 const double duplicated_vertex_tolerance = 1.0e-12,
2048 const bool copy_manifold_ids = false,
2049 const bool copy_boundary_ids = false);
2050
2101 template <int dim, int spacedim = dim>
2102 void
2104 const std::vector<unsigned int> &extents,
2106
2144 template <int dim, int spacedim>
2145 void
2147 const Triangulation<dim, spacedim> &triangulation_1,
2148 const Triangulation<dim, spacedim> &triangulation_2,
2150
2187 template <int dim, int spacedim>
2188 void
2190 const Triangulation<dim, spacedim> &input_triangulation,
2192 &cells_to_remove,
2194
2265 void
2267 const Triangulation<2, 2> &input,
2268 const unsigned int n_slices,
2269 const double height,
2270 Triangulation<3, 3> &result,
2271 const bool copy_manifold_ids = false,
2272 const std::vector<types::manifold_id> &manifold_priorities = {});
2273
2280 void
2282 const Triangulation<2, 2> &input,
2283 const unsigned int n_slices,
2284 const double height,
2285 Triangulation<2, 2> &result,
2286 const bool copy_manifold_ids = false,
2287 const std::vector<types::manifold_id> &manifold_priorities = {});
2288
2289
2308 void
2310 const Triangulation<2, 2> &input,
2311 const std::vector<double> &slice_coordinates,
2312 Triangulation<3, 3> &result,
2313 const bool copy_manifold_ids = false,
2314 const std::vector<types::manifold_id> &manifold_priorities = {});
2315
2322 void
2324 const Triangulation<2, 2> &input,
2325 const std::vector<double> &slice_coordinates,
2326 Triangulation<2, 2> &result,
2327 const bool copy_manifold_ids = false,
2328 const std::vector<types::manifold_id> &manifold_priorities = {});
2329
2330
2331
2371 template <int dim, int spacedim1, int spacedim2>
2372 void
2375
2415 template <int dim, int spacedim>
2416 void
2419
2420 // Doxygen will not show the function above if we include the
2421 // specialization.
2422#ifndef DOXYGEN
2426 template <int spacedim>
2427 void
2429 Triangulation<1, spacedim> &out_tria);
2430#endif
2431
2432
2455 template <int dim, int spacedim>
2456 void
2459
2464 namespace Airfoil
2465 {
2471 {
2476 std::string airfoil_type;
2477
2485 std::string naca_id;
2486
2493
2499
2504 double height;
2505
2510
2527
2535
2539 unsigned int refinements;
2540
2544 unsigned int n_subdivision_x_0;
2545
2549 unsigned int n_subdivision_x_1;
2550
2554 unsigned int n_subdivision_x_2;
2555
2559 unsigned int n_subdivision_y;
2560
2573
2578
2584 void
2586 };
2587
2617 template <int dim>
2618 void
2621 const AdditionalData &additional_data = AdditionalData());
2622
2623
2624
2637 template <int dim>
2638 void
2641 std::vector<GridTools::PeriodicFacePair<
2642 typename Triangulation<dim, dim>::cell_iterator>> &periodic_faces,
2643 const AdditionalData &additional_data = AdditionalData());
2644
2645 } // namespace Airfoil
2646
2668 template <int dim, int spacedim>
2669 void
2672 const std::vector<unsigned int> &repetitions,
2673 const Point<dim> &p1,
2674 const Point<dim> &p2,
2675 const bool colorize = false,
2676 const bool periodic = false);
2677
2700 template <int dim, int spacedim>
2701 void
2703 const unsigned int repetitions,
2704 const double p1 = 0.0,
2705 const double p2 = 1.0,
2706 const bool colorize = false,
2707 const bool periodic = false);
2708
2718#ifdef _MSC_VER
2719 // Microsoft's VC++ has a bug where it doesn't want to recognize that
2720 // an implementation (definition) of the extract_boundary_mesh function
2721 // matches a declaration. This can apparently only be avoided by
2722 // doing some contortion with the return type using the following
2723 // intermediate type. This is only used when using MS VC++ and uses
2724 // the direct way of doing it otherwise
2725 template <template <int, int> class MeshType, int dim, int spacedim>
2727 (concepts::is_triangulation_or_dof_handler<MeshType<dim, spacedim>>))
2728 struct ExtractBoundaryMesh
2729 {
2730 using return_type =
2731 std::map<typename MeshType<dim - 1, spacedim>::cell_iterator,
2732 typename MeshType<dim, spacedim>::face_iterator>;
2733 };
2734#endif
2735
2833 template <template <int, int> class MeshType, int dim, int spacedim>
2835 (concepts::is_triangulation_or_dof_handler<MeshType<dim, spacedim>>))
2836#ifdef DOXYGEN
2837 return_type
2838#else
2839# ifndef _MSC_VER
2840 std::map<typename MeshType<dim - 1, spacedim>::cell_iterator,
2841 typename MeshType<dim, spacedim>::face_iterator>
2842# else
2843 typename ExtractBoundaryMesh<MeshType, dim, spacedim>::return_type
2844# endif
2845#endif
2846 extract_boundary_mesh(const MeshType<dim, spacedim> &volume_mesh,
2847 MeshType<dim - 1, spacedim> &surface_mesh,
2848 const std::set<types::boundary_id> &boundary_ids =
2849 std::set<types::boundary_id>());
2850
2868 int,
2869 << "The number of repetitions " << arg1 << " must be >=1.");
2874 int,
2875 << "The vector of repetitions must have " << arg1
2876 << " elements.");
2877
2882 "The input to this function is oriented in a way that will"
2883 " cause all cells to have negative measure.");
2886#ifndef DOXYGEN
2887 // These functions are only implemented with specializations; declare them
2888 // here
2889 template <>
2890 void
2892 const double,
2893 const double,
2894 const double,
2895 const unsigned int,
2896 const bool);
2897
2898 template <>
2899 void
2901 const double,
2902 const double,
2903 const double,
2904 const unsigned int,
2905 const bool);
2906
2907 template <>
2908 void
2910 const double,
2911 const double,
2912 const double,
2913 const unsigned int,
2914 const bool);
2915
2916 template <>
2917 void
2919 const double,
2920 const unsigned int,
2921 const double,
2922 const bool);
2923
2924 template <>
2925 void
2927 const double,
2928 const unsigned int,
2929 const double,
2930 const bool);
2931
2932 template <>
2933 void
2935 const double,
2936 const unsigned int,
2937 const double,
2938 const bool);
2939
2940
2941
2942#endif
2943} // namespace GridGenerator
2944
2945
2946
2948
2949#endif
Definition point.h:113
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:35
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:243
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:36
bool colorize
Definition grid_out.cc:4634
#define DeclException0(Exception0)
static ::ExceptionBase & ExcInvalidInputOrientation()
static ::ExceptionBase & ExcInvalidRepetitions(int arg1)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcInvalidRadii()
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcInvalidRepetitionsDimension(int arg1)
std::vector< index_type > data
Definition mpi.cc:746
void create_triangulation(Triangulation< dim, dim > &tria, const AdditionalData &additional_data=AdditionalData())
void parallelepiped(Triangulation< dim > &tria, const Point< dim >(&corners)[dim], const bool colorize=false)
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_cube_with_cylindrical_hole(Triangulation< dim, spacedim > &triangulation, const double inner_radius=.25, const double outer_radius=.5, const double L=.5, const unsigned int repetitions=1, const bool colorize=false)
void hyper_ball_balanced(Triangulation< dim > &tria, const Point< dim > &center=Point< dim >(), const double radius=1.)
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 generate_from_name_and_arguments(Triangulation< dim, spacedim > &tria, const std::string &grid_generator_function_name, const std::string &grid_generator_function_arguments)
void enclosed_hyper_cube(Triangulation< dim > &tria, const double left=0., const double right=1., const double thickness=1., const bool colorize=false)
void replicate_triangulation(const Triangulation< dim, spacedim > &input, const std::vector< unsigned int > &extents, Triangulation< dim, spacedim > &result)
Replicate a given triangulation in multiple coordinate axes.
void parallelogram(Triangulation< dim > &tria, const Point< dim >(&corners)[dim], const bool colorize=false)
void general_cell(Triangulation< dim, spacedim > &tria, const std::vector< Point< spacedim > > &vertices, const bool colorize=false)
void subdivided_hyper_cube(Triangulation< dim, spacedim > &tria, const unsigned int repetitions, const double left=0., const double right=1., const bool colorize=false)
void hyper_shell(Triangulation< dim, spacedim > &tria, const Point< spacedim > &center, const double inner_radius, const double outer_radius, const unsigned int n_cells=0, bool colorize=false)
void hyper_L(Triangulation< dim > &tria, const double left=-1., const double right=1., const bool colorize=false)
void hyper_cube_slit(Triangulation< dim > &tria, const double left=0., const double right=1., const bool colorize=false)
void eccentric_hyper_shell(Triangulation< dim > &triangulation, const Point< dim > &inner_center, const Point< dim > &outer_center, const double inner_radius, const double outer_radius, const unsigned int n_cells)
void hyper_rectangle(Triangulation< dim, spacedim > &tria, const Point< dim > &p1, const Point< dim > &p2, const bool colorize=false)
void surface_mesh_to_volumetric_mesh(const Triangulation< 2, 3 > &surface_tria, Triangulation< 3 > &vol_tria, const CGALWrappers::AdditionalData< 3 > &data=CGALWrappers::AdditionalData< 3 >{})
void subdivided_hyper_cube_with_simplices(Triangulation< dim, spacedim > &tria, const unsigned int repetitions, const double p1=0.0, const double p2=1.0, const bool colorize=false, const bool periodic=false)
void cylinder(Triangulation< dim > &tria, const double radius=1., const double half_length=1.)
void moebius(Triangulation< 3, 3 > &tria, const unsigned int n_cells, const unsigned int n_rotations, const double R, const double r)
void pipe_junction(Triangulation< dim, spacedim > &tria, const std::vector< std::pair< Point< spacedim >, double > > &openings, const std::pair< Point< spacedim >, double > &bifurcation, const double aspect_ratio=0.5)
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 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 quarter_hyper_ball(Triangulation< dim > &tria, const Point< dim > &center=Point< dim >(), const double radius=1.)
void cheese(Triangulation< dim, spacedim > &tria, const std::vector< unsigned int > &holes)
Rectangular domain with rectangular pattern of holes.
void implicit_function(Triangulation< dim, 3 > &tria, const Function< 3 > &implicit_function, const CGALWrappers::AdditionalData< dim > &data=CGALWrappers::AdditionalData< dim >{}, const Point< 3 > &interior_point=Point< 3 >(), const double &outer_ball_radius=1.0)
void create_union_triangulation(const Triangulation< dim, spacedim > &triangulation_1, const Triangulation< dim, spacedim > &triangulation_2, Triangulation< dim, spacedim > &result)
void hyper_ball(Triangulation< dim, spacedim > &tria, const Point< spacedim > &center={}, const double radius=1., const bool attach_spherical_manifold_on_boundary_cells=false)
void non_standard_orientation_mesh(Triangulation< 2 > &tria, const unsigned int n_rotate_middle_square)
return_type 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 alfeld_split_of_simplex_mesh(const Triangulation< dim, spacedim > &in_tria, Triangulation< dim, spacedim > &out_tria)
void subdivided_parallelepiped(Triangulation< dim > &tria, const unsigned int n_subdivisions, const Point< dim >(&corners)[dim], const bool colorize=false)
void subdivided_cylinder(Triangulation< dim > &tria, const unsigned int x_subdivisions, const double radius=1., const double half_length=1.)
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 subdivided_hyper_L(Triangulation< dim, spacedim > &tria, const std::vector< unsigned int > &repetitions, const Point< dim > &bottom_left, const Point< dim > &top_right, const std::vector< int > &n_cells_to_remove)
void subdivided_hyper_rectangle_with_simplices(Triangulation< dim, spacedim > &tria, const std::vector< unsigned int > &repetitions, const Point< dim > &p1, const Point< dim > &p2, const bool colorize=false, const bool periodic=false)
void hyper_sphere(Triangulation< spacedim - 1, spacedim > &tria, const Point< spacedim > &center=Point< spacedim >(), const double radius=1.)
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 convert_hypercube_to_simplex_mesh(const Triangulation< dim, spacedim > &in_tria, Triangulation< dim, spacedim > &out_tria)
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 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 hyper_cube(Triangulation< dim, spacedim > &tria, const double left=0., const double right=1., const bool colorize=false)
void torus(Triangulation< dim, spacedim > &tria, const double centerline_radius, const double inner_radius, const unsigned int n_cells_toroidal=6, const double phi=2.0 *numbers::PI)
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 simplex(Triangulation< dim, dim > &tria, const std::vector< Point< dim > > &vertices)
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 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, const bool copy_boundary_ids=false)
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
void half_hyper_ball(Triangulation< dim > &tria, const Point< dim > &center=Point< dim >(), const double radius=1.)
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, const bool colorize=false)
void flatten_triangulation(const Triangulation< dim, spacedim1 > &in_tria, Triangulation< dim, spacedim2 > &out_tria)
constexpr double PI
Definition numbers.h:240
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
void add_parameters(ParameterHandler &prm)