]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Test for p::d::Tria::load(forest). 11721/head
authorMarc Fehling <mafehling.git@gmail.com>
Wed, 10 Feb 2021 02:02:22 +0000 (19:02 -0700)
committerMarc Fehling <marc.fehling@gmx.net>
Thu, 11 Feb 2021 04:31:16 +0000 (21:31 -0700)
doc/news/changes/minor/20210209Munch
include/deal.II/distributed/tria.h
source/distributed/tria.cc
tests/mpi/p4est_load_01.cc [new file with mode: 0644]
tests/mpi/p4est_load_01.with_p4est=true.mpirun=2.output [new file with mode: 0644]

index e8e7e83aaa9a866d2b98e9c3adfafee206cd620b..74879eeeeec7f2194d1baa077bd675892c807129 100644 (file)
@@ -1,6 +1,6 @@
-New: The method parallel::distributed::Triagnulation::load() 
+New: The method parallel::distributed::Triangulation::load()
 can now also accept a p4est/p8est forest, which can be queried
-from an existing triangulation via parallel::distributed::Triagnulation::get_p4est().
+from an existing triangulation via parallel::distributed::Triangulation::get_p4est().
 <br>
 (Marc Fehling, Peter Munch, 2021/02/09)
 
index 238db08b55c627375bd8552e59d94ba86a7ebf38..dda7535d94fb1f6ebf15e101175bf0bd427c838f 100644 (file)
@@ -607,7 +607,7 @@ namespace parallel
       /**
        * Load the refinement information from a given parallel forest. This
        * forest might be obtained from the function call to
-       * parallel::distributed::Triangulation::get_p4est();
+       * parallel::distributed::Triangulation::get_p4est().
        */
       void
       load(const typename dealii::internal::p4est::types<dim>::forest *forest);
index cf2b51eae344df8074dafb893665599bf7ba9eb6..ceb9d84b5fcaccfdb16229bff81c4cde92581a81 100644 (file)
@@ -1608,11 +1608,22 @@ namespace parallel
     }
 
 
+
     template <int dim, int spacedim>
     void
     Triangulation<dim, spacedim>::load(
       const typename dealii::internal::p4est::types<dim>::forest *forest)
     {
+      Assert(this->n_cells() > 0,
+             ExcMessage(
+               "load() only works if the Triangulation already contains "
+               "a coarse mesh!"));
+      Assert(this->n_cells() == forest->trees->elem_count,
+             ExcMessage(
+               "Coarse mesh of the Triangulation does not match the one "
+               "of the provided forest!"));
+
+      // clear the old forest
       if (parallel_ghost != nullptr)
         {
           dealii::internal::p4est::functions<dim>::ghost_destroy(
@@ -1625,7 +1636,7 @@ namespace parallel
       // note: we can keep the connectivity, since the coarse grid does not
       // change
 
-      // create deep copy
+      // create deep copy of the new forest
       typename dealii::internal::p4est::types<dim>::forest *temp =
         const_cast<typename dealii::internal::p4est::types<dim>::forest *>(
           forest);
diff --git a/tests/mpi/p4est_load_01.cc b/tests/mpi/p4est_load_01.cc
new file mode 100644 (file)
index 0000000..7e4119d
--- /dev/null
@@ -0,0 +1,62 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 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.
+//
+// ---------------------------------------------------------------------
+
+
+
+// load a triangulation from an existing forest
+
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include "../tests.h"
+
+
+template <int dim>
+void
+test()
+{
+  // prepare two independent triangulations with the same coarse mesh
+  parallel::distributed::Triangulation<dim> tria1(MPI_COMM_WORLD),
+    tria2(MPI_COMM_WORLD);
+  GridGenerator::subdivided_hyper_cube(tria1, 2);
+  GridGenerator::subdivided_hyper_cube(tria2, 2);
+
+  // refine the first, and transfer the forest to the other
+  tria1.refine_global(2);
+  tria2.load(tria1.get_p4est());
+
+  // verify that global information matches
+  AssertDimension(tria1.n_global_levels(), tria2.n_global_levels());
+  AssertDimension(tria1.n_global_active_cells(), tria2.n_global_active_cells());
+
+  deallog << "OK" << std::endl;
+}
+
+
+int
+main(int argc, char *argv[])
+{
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv, 1);
+  MPILogInitAll                    log;
+
+  deallog.push("2d");
+  test<2>();
+  deallog.pop();
+  deallog.push("3d");
+  test<3>();
+  deallog.pop();
+}
diff --git a/tests/mpi/p4est_load_01.with_p4est=true.mpirun=2.output b/tests/mpi/p4est_load_01.with_p4est=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..4aa9062
--- /dev/null
@@ -0,0 +1,7 @@
+
+DEAL:0:2d::OK
+DEAL:0:3d::OK
+
+DEAL:1:2d::OK
+DEAL:1:3d::OK
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.