]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fix remaining issues for #9195 9432/head
authorPeter Munch <peterrmuench@gmail.com>
Sun, 26 Jan 2020 19:44:03 +0000 (20:44 +0100)
committerPeter Munch <peterrmuench@gmail.com>
Mon, 27 Jan 2020 19:15:08 +0000 (20:15 +0100)
include/deal.II/grid/tria.h
include/deal.II/grid/tria_description.h
source/grid/tria_description.cc

index 97b4572d525a8cd1c5768ee94f23131f3cb83e46..621dab4caa5e81def84f726aa62996499be8ecb3 100644 (file)
@@ -1836,7 +1836,7 @@ public:
    * Create a triangulation from the provided
    * TriangulationDescription::Description.
    *
-   * @note The namespace construction::Utilities contains functions
+   * @note The namespace TriangulationDescription::Utilities contains functions
    *       to create TriangulationDescription::Description.
    *
    * @param construction_data The data needed for this process.
index d927f680b5d298ccc603ca3e29fca59bf3fe3bc3..9d414fad70c714f2f9a789f2b43bb36160e137aa 100644 (file)
@@ -271,16 +271,16 @@ namespace TriangulationDescription
    * the cell id, the subdomain_id and the level_subdomain_id as well as
    * information related to manifold_id and boundary_id.
    *
+   * @note Similarly to dealii::CellData, this structure stores information
+   * about a cell. However, in contrast to dealii::CellData, it also stores
+   * a unique id, partitioning information, and information related to cell
+   * faces and edges.
+   *
    * @author Peter Munch, 2019
    */
   template <int dim>
   struct CellData
   {
-    /**
-     * Constructor.
-     */
-    CellData() = default;
-
     /**
      * Boost serialization function
      */
@@ -312,7 +312,7 @@ namespace TriangulationDescription
     /**
      * Manifold id of the cell.
      */
-    types::material_id manifold_id;
+    types::manifold_id manifold_id;
 
     /**
      * Manifold id of all lines of the cell.
@@ -419,22 +419,22 @@ namespace TriangulationDescription
      *
      * @param tria Partitioned input triangulation.
      * @param comm MPI_Communicator to be used. In the case
-     *             of dealii::distributed::Triangulation, the communicators
-     *             have to match.
+     *   of dealii::distributed::Triangulation, the communicators have to match.
      * @param construct_multilevel_hierarchy Signal if the multigrid levels
      *        should be constructed.
-     * @param my_rank_in Construct Description for this rank (only
-     *                   working for serial triangulations).
-     * @return Description to be used to setup a Triangulation.
+     * @param my_rank_in Construct Description for the specified rank (only
+     *   working for serial triangulations that have been partitioned by
+     *   functions like GridToold::partition_triangulation()).
+     * @return Description to be used to set up a Triangulation.
      *
      * @note Multilevel hierarchies are supported if it is enabled in
-     *       parallel::fullydistributed::Triangulation.
+     *   parallel::fullydistributed::Triangulation.
      *
      * @note Hanging nodes in the input triangulation are supported. However,
-     *       to be able to use this
-     *       feature in the case of parallel::fullydistributed::Triangulation,
-     *       the user has to enable multilevel hierarchy support in
-     *       parallel::fullydistributed::Triangulation.
+     *   to be able to use this feature in the case of
+     *   parallel::fullydistributed::Triangulation, the user has to enable
+     *   multilevel hierarchy support in
+     *   parallel::fullydistributed::Triangulation.
      *
      * @author Peter Munch, 2019
      */
@@ -448,27 +448,27 @@ namespace TriangulationDescription
 
 
     /**
-     * Construct a construction::Description. In contrast
+     * Construct a TriangulationDescription::Description. In contrast
      * to the function above, this function is also responsible for creating
      * a serial triangulation and for its partitioning (by calling the
-     * provided std::functions). Internally only selected processes (every
-     * n-th/each root of a group of size group_size) create a serial
+     * provided `std::functions' objects). Internally only selected processes (
+     * every n-th/each root of a group of size group_size) create a serial
      * triangulation and the ConstructionData for all processes in its group,
      * which is communicated.
      *
      * @note A reasonable group size is the size of a NUMA domain or the
      * size of a compute node.
      *
-     * @param serial_grid_generator A function, which creates a serial triangulation.
-     * @param serial_grid_partitioner A function, which can partition a serial
-     *        triangulation, i.e., sets the sudomain_ids of the active cells.
-     *        The function takes as the first argument a serial triangulation,
-     *        as the second argument the MPI communicator, and as the third
-     *        argument the group size.
-     * @param comm MPI communicator
+     * @param serial_grid_generator A function which creates a serial triangulation.
+     * @param serial_grid_partitioner A function which can partition a serial
+     *   triangulation, i.e., sets the sudomain_ids of the active cells.
+     *   The function takes as the first argument a serial triangulation,
+     *   as the second argument the MPI communicator, and as the third
+     *   argument the group size.
+     * @param comm MPI communicator.
      * @param group_size The size of each group.
      * @param construct_multilevel_hierarchy Construct multigrid levels.
-     * @return Description to be used to setup a Triangulation.
+     * @return Description to be used to set up a Triangulation.
      *
      * @author Peter Munch, 2019
      */
index da532e6ae43e81fae8d8b6c9ed06812d9da59084..0905de30172d97cb1b825ce305f046010a19bcbc 100644 (file)
@@ -38,11 +38,12 @@ namespace TriangulationDescription
        */
       template <int dim, int spacedim>
       void
-      set_user_flag_reverse(TriaIterator<CellAccessor<dim, spacedim>> cell)
+      set_user_flag_and_of_its_parents(
+        const TriaIterator<CellAccessor<dim, spacedim>> &cell)
       {
         cell->set_user_flag();
         if (cell->level() != 0)
-          set_user_flag_reverse(cell->parent());
+          set_user_flag_and_of_its_parents(cell->parent());
       }
 
 
@@ -108,21 +109,20 @@ namespace TriangulationDescription
     {
       if (auto tria_pdt = dynamic_cast<
             const parallel::distributed::Triangulation<dim, spacedim> *>(&tria))
-        AssertThrow(comm == tria_pdt->get_communicator(),
-                    ExcMessage("MPI communicators do not match."));
+        Assert(comm == tria_pdt->get_communicator(),
+               ExcMessage("MPI communicators do not match."));
 
       // First, figure out for what rank we are supposed to build the
       // ConstructionData object
       unsigned int my_rank = my_rank_in;
-      AssertThrow(my_rank == numbers::invalid_unsigned_int ||
-                    my_rank < dealii::Utilities::MPI::n_mpi_processes(comm),
-                  ExcMessage(
-                    "Rank has to be smaller than available processes."));
+      Assert(my_rank == numbers::invalid_unsigned_int ||
+               my_rank < dealii::Utilities::MPI::n_mpi_processes(comm),
+             ExcMessage("Rank has to be smaller than available processes."));
 
       if (auto tria_pdt = dynamic_cast<
             const parallel::distributed::Triangulation<dim, spacedim> *>(&tria))
         {
-          AssertThrow(
+          Assert(
             my_rank == numbers::invalid_unsigned_int ||
               my_rank == dealii::Utilities::MPI::this_mpi_process(comm),
             ExcMessage(
@@ -139,8 +139,8 @@ namespace TriangulationDescription
         }
       else
         {
-          AssertThrow(
-            false, ExcMessage("This type of triangulation is not supported!"));
+          Assert(false,
+                 ExcMessage("This type of triangulation is not supported!"));
         }
 
       Description<dim, spacedim> construction_data;
@@ -182,7 +182,7 @@ namespace TriangulationDescription
       // check if multilevel hierarchy should be constructed
       if (construct_multilevel_hierarchy == false)
         {
-          AssertThrow(
+          Assert(
             tria.has_hanging_nodes() == false,
             ExcMessage(
               "Hanging nodes are only supported if multilevel hierarchy is constructed!"));
@@ -324,10 +324,9 @@ namespace TriangulationDescription
 
           // 1b) loop over levels (from fine to coarse) and mark on each level
           //     the locally relevant cells
-          for (unsigned int level =
-                 tria.get_triangulation().n_global_levels() - 1;
-               level != numbers::invalid_unsigned_int;
-               level--)
+          for (int level = tria.get_triangulation().n_global_levels() - 1;
+               level >= 0;
+               --level)
             {
               // collect vertices connected to a (on any level) locally owned
               // cell
@@ -361,7 +360,7 @@ namespace TriangulationDescription
               // mark all locally relevant cells
               for (auto cell : tria.cell_iterators_on_level(level))
                 if (is_locally_relevant_on_level(cell))
-                  set_user_flag_reverse(cell);
+                  set_user_flag_and_of_its_parents(cell);
             }
 
           // 2) set_up coarse-grid triangulation
@@ -589,9 +588,12 @@ namespace TriangulationDescription
             std::min(group_root + group_size,
                      dealii::Utilities::MPI::n_mpi_processes(comm));
 
-          // 3) create ConstructionData for the other processes in group
+          // 3) create Description for the other processes in group; since
+          // we expect that this function is called for huge meshes, one
+          // Description is created at a time and send away; only once the
+          // Description has been sent away, the next rank is processed.
           for (unsigned int other_rank = group_root + 1; other_rank < end_group;
-               other_rank++)
+               ++other_rank)
             {
               // 3a) create construction data for other ranks
               const auto construction_data =
@@ -631,7 +633,7 @@ namespace TriangulationDescription
                           len,
                           MPI_CHAR,
                           status.MPI_SOURCE,
-                          status.MPI_TAG,
+                          mpi_tag,
                           comm,
                           &status);
           AssertThrowMPI(ierr);

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.