From 4de22c5b84b3df77a36e8b2d0de7423fb33f7c1d Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Wed, 22 Nov 2023 16:03:30 +0100 Subject: [PATCH] New test case --- tests/distributed_grids/create_torus_mesh.cc | 59 +++++++++++++++++++ .../create_torus_mesh.output | 2 + 2 files changed, 61 insertions(+) create mode 100644 tests/distributed_grids/create_torus_mesh.cc create mode 100644 tests/distributed_grids/create_torus_mesh.output diff --git a/tests/distributed_grids/create_torus_mesh.cc b/tests/distributed_grids/create_torus_mesh.cc new file mode 100644 index 0000000000..c84f0a7c18 --- /dev/null +++ b/tests/distributed_grids/create_torus_mesh.cc @@ -0,0 +1,59 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 2023 - 2023 by the deal.II authors +// +// This file is part of the deal.II library. +// +// The deal.II library is free software; you can use it, redistribute +// it, and/or modify it under the terms of the GNU Lesser General +// Public License as published by the Free Software Foundation; either +// version 2.1 of the License, or (at your option) any later version. +// The full text of the license can be found in the file LICENSE.md at +// the top level directory of deal.II. +// +// --------------------------------------------------------------------- + + + +// Test GridGenerator::torus, which used to run into an assertion with p4est +// due to bad orientations before 2023, see +// https://github.com/dealii/dealii/issues/16272 + +#include + +#include + +#include "../tests.h" + + +int +main(int argc, char *argv[]) +{ + Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1); + + initlog(); + + // create mesh + const int dim = 3; + + parallel::distributed::Triangulation triangulation( + MPI_COMM_WORLD, + typename Triangulation::MeshSmoothing( + Triangulation::smoothing_on_refinement | + Triangulation::smoothing_on_coarsening)); + + const double R = 2.; + const double r = 1.; + const unsigned int n_cells_toroidal = 1; + const double phi = 0.5 * dealii::numbers::PI; + dealii::GridGenerator::torus( + triangulation, R, r, n_cells_toroidal, phi); + triangulation.reset_all_manifolds(); + triangulation.refine_global(1); + + deallog << "Tria info: " << triangulation.n_vertices() << " " + << triangulation.n_cells() << " " << triangulation.n_faces() << " " + << triangulation.n_raw_lines() << std::endl; + + return 0; +} diff --git a/tests/distributed_grids/create_torus_mesh.output b/tests/distributed_grids/create_torus_mesh.output new file mode 100644 index 0000000000..3cd4c0231a --- /dev/null +++ b/tests/distributed_grids/create_torus_mesh.output @@ -0,0 +1,2 @@ + +DEAL::Tria info: 75 45 170 214 -- 2.39.5