Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.4.1
\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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
grid_generator.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 1999 - 2022 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
25#include <deal.II/base/point.h>
26#include <deal.II/base/table.h>
27
29#include <deal.II/grid/tria.h>
30
32
33#include <array>
34#include <map>
35
37
61{
66
95 template <int dim, int spacedim>
96 void
98 const double left = 0.,
99 const double right = 1.,
100 const bool colorize = false);
101
129 template <int dim>
130 void
132 const std::vector<Point<dim>> &vertices);
133
139 template <int dim, int spacedim>
140 void
143
144
171 template <int dim, int spacedim>
172 void
174 const unsigned int repetitions,
175 const double left = 0.,
176 const double right = 1.,
177 const bool colorize = false);
178
214 template <int dim, int spacedim>
215 void
217 const Point<dim> & p1,
218 const Point<dim> & p2,
219 const bool colorize = false);
220
270 template <int dim, int spacedim>
271 void
273 const std::vector<unsigned int> &repetitions,
274 const Point<dim> & p1,
275 const Point<dim> & p2,
276 const bool colorize = false);
277
293 template <int dim>
294 void
296 const std::vector<std::vector<double>> &step_sizes,
297 const Point<dim> & p_1,
298 const Point<dim> & p_2,
299 const bool colorize = false);
300
311 template <int dim>
312 void
314 const std::vector<std::vector<double>> &spacing,
315 const Point<dim> & p,
316 const Table<dim, types::material_id> &material_id,
317 const bool colorize = false);
318
343 template <int dim, int spacedim>
344 void
346 const std::vector<unsigned int> &holes);
347
391 template <int dim>
392 void
394 const double inner_radius = 0.4,
395 const double outer_radius = 1.,
396 const double pad_bottom = 2.,
397 const double pad_top = 2.,
398 const double pad_left = 1.,
399 const double pad_right = 1.,
400 const Point<dim> & center = Point<dim>(),
401 const types::manifold_id polar_manifold_id = 0,
402 const types::manifold_id tfi_manifold_id = 1,
403 const double L = 1.,
404 const unsigned int n_slices = 2,
405 const bool colorize = false);
406
493 template <int dim>
494 void
496 const double shell_region_width = 0.03,
497 const unsigned int n_shells = 2,
498 const double skewness = 2.0,
499 const bool colorize = false);
500
521 template <int dim, int spacedim>
522 void
524 const std::vector<Point<spacedim>> &vertices,
525 const bool colorize = false);
526
543 template <int dim>
544 void
546 const Point<dim> (&corners)[dim],
547 const bool colorize = false);
548
564 template <int dim>
565 void
567 const Point<dim> (&corners)[dim],
568 const bool colorize = false);
569
581 template <int dim>
582 void
584 const unsigned int n_subdivisions,
585 const Point<dim> (&corners)[dim],
586 const bool colorize = false);
587
596 template <int dim>
597 void
599#ifndef _MSC_VER
600 const unsigned int (&n_subdivisions)[dim],
601#else
602 const unsigned int *n_subdivisions,
603#endif
604 const Point<dim> (&corners)[dim],
605 const bool colorize = false);
606
630 template <int dim, int spacedim>
631 void
633 const Point<spacedim> & origin,
634 const std::array<Tensor<1, spacedim>, dim> &edges,
635 const std::vector<unsigned int> &subdivisions = {},
636 const bool colorize = false);
637
654 template <int dim>
655 void
657 const double left = 0.,
658 const double right = 1.,
659 const double thickness = 1.,
660 const bool colorize = false);
661
738 template <int dim>
739 void
741 const Point<dim> & center = Point<dim>(),
742 const double radius = 1.,
743 const bool attach_spherical_manifold_on_boundary_cells = false);
744
777 template <int dim>
778 void
780 const Point<dim> & center = Point<dim>(),
781 const double radius = 1.);
782
801 void
803 const unsigned int n_rotate_middle_square);
804
826 void
828 const bool face_orientation,
829 const bool face_flip,
830 const bool face_rotation,
831 const bool manipulate_left_cube);
832
833
859 template <int spacedim>
860 void
863 const double radius = 1.);
864
898 template <int dim>
899 void
901 const Point<dim> & center = Point<dim>(),
902 const double radius = 1.);
903
916 template <int dim>
917 void
919 const Point<dim> & center = Point<dim>(),
920 const double radius = 1.);
921
945 template <int dim>
946 void
948 const double radius = 1.,
949 const double half_length = 1.);
950
951
985 template <int dim>
986 void
988 const unsigned int x_subdivisions,
989 const double radius = 1.,
990 const double half_length = 1.);
991
992
1023 template <int dim>
1024 void
1026 const double radius_0 = 1.0,
1027 const double radius_1 = 0.5,
1028 const double half_length = 1.0);
1029
1153 template <int dim, int spacedim>
1154 void
1156 const std::vector<std::pair<Point<spacedim>, double>> &openings,
1157 const std::pair<Point<spacedim>, double> &bifurcation,
1158 const double aspect_ratio = 0.5);
1159
1186 template <int dim, int spacedim>
1187 void
1189 const std::vector<unsigned int> &sizes,
1190 const bool colorize_cells = false);
1191
1227 template <int dim>
1228 void
1230 const double left = -1.,
1231 const double right = 1.,
1232 const bool colorize = false);
1233
1265 template <int dim, int spacedim>
1266 void
1268 const std::vector<unsigned int> &repetitions,
1269 const Point<dim> & bottom_left,
1270 const Point<dim> & top_right,
1271 const std::vector<int> & n_cells_to_remove);
1272
1292 template <int dim>
1293 void
1295 const double left = 0.,
1296 const double right = 1.,
1297 const bool colorize = false);
1298
1356 template <int dim>
1357 void
1359 const Point<dim> & center,
1360 const double inner_radius,
1361 const double outer_radius,
1362 const unsigned int n_cells = 0,
1363 bool colorize = false);
1364
1398 template <int dim>
1399 void
1401 const Point<dim> & inner_center,
1402 const Point<dim> & outer_center,
1403 const double inner_radius,
1404 const double outer_radius,
1405 const unsigned int n_cells);
1406
1436 template <int dim>
1437 void
1439 const Point<dim> & center,
1440 const double inner_radius,
1441 const double outer_radius,
1442 const unsigned int n_cells = 0,
1443 const bool colorize = false);
1444
1445
1471 template <int dim>
1472 void
1474 const Point<dim> & center,
1475 const double inner_radius,
1476 const double outer_radius,
1477 const unsigned int n_cells = 0,
1478 const bool colorize = false);
1479
1506 template <int dim>
1507 void
1509 const double length,
1510 const double inner_radius,
1511 const double outer_radius,
1512 const unsigned int n_radial_cells = 0,
1513 const unsigned int n_axial_cells = 0);
1514
1564 template <int dim, int spacedim>
1565 void
1567 const double R,
1568 const double r,
1569 const unsigned int n_cells_toroidal = 6,
1570 const double phi = 2.0 * numbers::PI);
1571
1601 template <int dim>
1602 void
1604 const double inner_radius = .25,
1605 const double outer_radius = .5,
1606 const double L = .5,
1607 const unsigned int repetitions = 1,
1608 const bool colorize = false);
1609
1680 template <int dim>
1681 void
1683 const Point<dim> & center,
1684 const double inner_radius = 0.125,
1685 const double outer_radius = 0.25,
1686 const unsigned int n_shells = 1,
1687 const double skewness = 0.1,
1688 const unsigned int n_cells_per_shell = 0,
1689 const bool colorize = false);
1690
1704 void
1706 const unsigned int n_cells,
1707 const unsigned int n_rotations,
1708 const double R,
1709 const double r);
1710
1746 template <int dim, int spacedim>
1747 void
1750 const std::string & grid_generator_function_name,
1751 const std::string & grid_generator_function_arguments);
1752
1802 template <int dim>
1803 void
1808 const Point<3> &interior_point = Point<3>(),
1809 const double & outer_ball_radius = 1.0);
1810
1825 void
1827 Triangulation<3> & vol_tria,
1831
1836
1900 template <int dim, int spacedim>
1901 void
1903 const Triangulation<dim, spacedim> &triangulation_2,
1905 const double duplicated_vertex_tolerance = 1.0e-12,
1906 const bool copy_manifold_ids = false);
1907
1922 template <int dim, int spacedim>
1923 void
1925 const std::vector<const Triangulation<dim, spacedim> *> &triangulations,
1927 const double duplicated_vertex_tolerance = 1.0e-12,
1928 const bool copy_manifold_ids = false);
1929
1977 template <int dim, int spacedim = dim>
1978 void
1980 const std::vector<unsigned int> & extents,
1982
2020 template <int dim, int spacedim>
2021 void
2023 const Triangulation<dim, spacedim> &triangulation_1,
2024 const Triangulation<dim, spacedim> &triangulation_2,
2026
2063 template <int dim, int spacedim>
2064 void
2066 const Triangulation<dim, spacedim> &input_triangulation,
2068 & cells_to_remove,
2070
2127 void
2129 const Triangulation<2, 2> & input,
2130 const unsigned int n_slices,
2131 const double height,
2132 Triangulation<3, 3> & result,
2133 const bool copy_manifold_ids = false,
2134 const std::vector<types::manifold_id> &manifold_priorities = {});
2135
2142 void
2144 const Triangulation<2, 2> & input,
2145 const unsigned int n_slices,
2146 const double height,
2147 Triangulation<2, 2> & result,
2148 const bool copy_manifold_ids = false,
2149 const std::vector<types::manifold_id> &manifold_priorities = {});
2150
2151
2170 void
2172 const Triangulation<2, 2> & input,
2173 const std::vector<double> & slice_coordinates,
2174 Triangulation<3, 3> & result,
2175 const bool copy_manifold_ids = false,
2176 const std::vector<types::manifold_id> &manifold_priorities = {});
2177
2184 void
2186 const Triangulation<2, 2> & input,
2187 const std::vector<double> & slice_coordinates,
2188 Triangulation<2, 2> & result,
2189 const bool copy_manifold_ids = false,
2190 const std::vector<types::manifold_id> &manifold_priorities = {});
2191
2192
2193
2230 template <int dim, int spacedim1, int spacedim2>
2231 void
2234
2268 template <int dim, int spacedim>
2269 void
2272
2276 template <int spacedim>
2277 void
2279 Triangulation<1, spacedim> & out_tria);
2280
2281
2288 namespace Airfoil
2289 {
2295 {
2300 std::string airfoil_type;
2301
2309 std::string naca_id;
2310
2317
2323
2328 double height;
2329
2334
2351
2359
2363 unsigned int refinements;
2364
2368 unsigned int n_subdivision_x_0;
2369
2373 unsigned int n_subdivision_x_1;
2374
2378 unsigned int n_subdivision_x_2;
2379
2383 unsigned int n_subdivision_y;
2384
2397
2402
2408 void
2410 };
2411
2441 template <int dim>
2442 void
2445 const AdditionalData & additional_data = AdditionalData());
2446
2447
2448
2461 template <int dim>
2462 void
2465 std::vector<GridTools::PeriodicFacePair<
2466 typename Triangulation<dim, dim>::cell_iterator>> &periodic_faces,
2467 const AdditionalData &additional_data = AdditionalData());
2468
2469 } // namespace Airfoil
2470
2484 template <int dim, int spacedim>
2485 void
2486 subdivided_hyper_rectangle_with_simplices(
2488 const std::vector<unsigned int> &repetitions,
2489 const Point<dim> & p1,
2490 const Point<dim> & p2,
2491 const bool colorize = false);
2492
2507 template <int dim, int spacedim>
2508 void
2509 subdivided_hyper_cube_with_simplices(Triangulation<dim, spacedim> &tria,
2510 const unsigned int repetitions,
2511 const double p1 = 0.0,
2512 const double p2 = 1.0,
2513 const bool colorize = false);
2514
2516
2523
2524#ifdef _MSC_VER
2525 // Microsoft's VC++ has a bug where it doesn't want to recognize that
2526 // an implementation (definition) of the extract_boundary_mesh function
2527 // matches a declaration. This can apparently only be avoided by
2528 // doing some contortion with the return type using the following
2529 // intermediate type. This is only used when using MS VC++ and uses
2530 // the direct way of doing it otherwise
2531 template <template <int, int> class MeshType, int dim, int spacedim>
2532 struct ExtractBoundaryMesh
2533 {
2534 using return_type =
2535 std::map<typename MeshType<dim - 1, spacedim>::cell_iterator,
2536 typename MeshType<dim, spacedim>::face_iterator>;
2537 };
2538#endif
2539
2634 template <template <int, int> class MeshType, int dim, int spacedim>
2635#ifdef DOXYGEN
2636 return_type
2637#else
2638# ifndef _MSC_VER
2639 std::map<typename MeshType<dim - 1, spacedim>::cell_iterator,
2640 typename MeshType<dim, spacedim>::face_iterator>
2641# else
2642 typename ExtractBoundaryMesh<MeshType, dim, spacedim>::return_type
2643# endif
2644#endif
2645 extract_boundary_mesh(const MeshType<dim, spacedim> & volume_mesh,
2646 MeshType<dim - 1, spacedim> & surface_mesh,
2647 const std::set<types::boundary_id> &boundary_ids =
2648 std::set<types::boundary_id>());
2649
2651
2652
2657
2658
2667 int,
2668 << "The number of repetitions " << arg1 << " must be >=1.");
2673 int,
2674 << "The vector of repetitions must have " << arg1
2675 << " elements.");
2676
2681 "The input to this function is oriented in a way that will"
2682 " cause all cells to have negative measure.");
2684
2685#ifndef DOXYGEN
2686 // These functions are only implemented with specializations; declare them
2687 // here
2688 template <>
2689 void
2691 const double,
2692 const double,
2693 const double,
2694 const unsigned int,
2695 const bool);
2696
2697 template <>
2698 void
2700 const double,
2701 const double,
2702 const double,
2703 const unsigned int,
2704 const bool);
2705
2706 template <>
2707 void
2709 const double,
2710 const double,
2711 const double,
2712 const unsigned int,
2713 const bool);
2714
2715 template <>
2716 void
2718 const double,
2719 const unsigned int,
2720 const double,
2721 const bool);
2722
2723 template <>
2724 void
2726 const double,
2727 const unsigned int,
2728 const double,
2729 const bool);
2730
2731 template <>
2732 void
2734 const double,
2735 const unsigned int,
2736 const double,
2737 const bool);
2738
2739
2740
2741#endif
2742} // namespace GridGenerator
2743
2744
2745
2747
2748#endif
Definition: point.h:111
Definition: tensor.h:503
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:442
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:443
Point< 3 > center
Point< 3 > vertices[4]
bool colorize
Definition: grid_out.cc:4605
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcInvalidInputOrientation()
static ::ExceptionBase & ExcInvalidRepetitions(int arg1)
#define DeclExceptionMsg(Exception, defaulttext)
Definition: exceptions.h:487
static ::ExceptionBase & ExcInvalidRadii()
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:509
static ::ExceptionBase & ExcInvalidRepetitionsDimension(int arg1)
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_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 convert_hypercube_to_simplex_mesh(const Triangulation< 1, spacedim > &in_tria, Triangulation< 1, spacedim > &out_tria)
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_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_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 hyper_ball(Triangulation< dim > &tria, const Point< dim > &center=Point< dim >(), const double radius=1., const bool attach_spherical_manifold_on_boundary_cells=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 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 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 cheese(Triangulation< dim, spacedim > &tria, const std::vector< unsigned int > &holes)
Rectangular domain with rectangular pattern of holes.
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)
void torus(Triangulation< dim, spacedim > &tria, const double R, const double r, const unsigned int n_cells_toroidal=6, const double phi=2.0 *numbers::PI)
void create_union_triangulation(const Triangulation< dim, spacedim > &triangulation_1, const Triangulation< dim, spacedim > &triangulation_2, Triangulation< dim, spacedim > &result)
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 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 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 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 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)
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 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)
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 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 flatten_triangulation(const Triangulation< dim, spacedim1 > &in_tria, Triangulation< dim, spacedim2 > &out_tria)
static constexpr double PI
Definition: numbers.h:233
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
void add_parameters(ParameterHandler &prm)
const ::Triangulation< dim, spacedim > & tria