22#ifndef APPLICATIONS_GRID_TOOLS_PERIODIC_BOX_H_
23#define APPLICATIONS_GRID_TOOLS_PERIODIC_BOX_H_
26#include <deal.II/distributed/tria.h>
27#include <deal.II/grid/grid_generator.h>
28#include <deal.II/grid/grid_tools.h>
31#include <exadg/grid/deformed_cube_manifold.h>
37create_periodic_box(dealii::Triangulation<dim> & triangulation,
38 unsigned int const n_refine_space,
39 std::vector<dealii::GridTools::PeriodicFacePair<
40 typename dealii::Triangulation<dim>::cell_iterator>> & periodic_faces,
41 unsigned int const n_subdivisions,
44 bool const curvilinear_mesh =
false,
45 double const deformation = 0.1)
47 dealii::GridGenerator::subdivided_hyper_cube(triangulation, n_subdivisions, left, right);
51 unsigned int const frequency = 2;
52 apply_deformed_cube_manifold(triangulation, left, right, deformation, frequency);
55 for(
auto const & cell : triangulation.cell_iterators())
57 for(
unsigned int const face_number : cell->face_indices())
60 if((std::fabs(cell->face(face_number)->center()(0) - left) < 1e-12))
61 cell->face(face_number)->set_all_boundary_ids(0);
62 else if((std::fabs(cell->face(face_number)->center()(0) - right) < 1e-12))
63 cell->face(face_number)->set_all_boundary_ids(1);
65 else if((std::fabs(cell->face(face_number)->center()(1) - left) < 1e-12))
66 cell->face(face_number)->set_all_boundary_ids(2);
67 else if((std::fabs(cell->face(face_number)->center()(1) - right) < 1e-12))
68 cell->face(face_number)->set_all_boundary_ids(3);
70 else if(dim == 3 and (std::fabs(cell->face(face_number)->center()(2) - left) < 1e-12))
71 cell->face(face_number)->set_all_boundary_ids(4);
72 else if(dim == 3 and (std::fabs(cell->face(face_number)->center()(2) - right) < 1e-12))
73 cell->face(face_number)->set_all_boundary_ids(5);
77 dealii::GridTools::collect_periodic_faces(triangulation, 0, 1, 0 , periodic_faces);
78 dealii::GridTools::collect_periodic_faces(triangulation, 2, 3, 1 , periodic_faces);
80 dealii::GridTools::collect_periodic_faces(
81 triangulation, 4, 5, 2 , periodic_faces);
83 triangulation.add_periodicity(periodic_faces);
86 triangulation.refine_global(n_refine_space);