]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Improve some checks for distributed Triangulations. 12495/head
authorDavid Wells <drwells@email.unc.edu>
Thu, 24 Jun 2021 20:26:02 +0000 (16:26 -0400)
committerDavid Wells <drwells@email.unc.edu>
Fri, 25 Jun 2021 18:32:11 +0000 (14:32 -0400)
include/deal.II/numerics/error_estimator.templates.h
source/dofs/dof_tools_sparsity.cc
source/numerics/error_estimator_1d.cc

index f0bdf443231d53128aec530f9c055efdbac267f5..10b0f655509fb68f744f216744c5e618a1190daf 100644 (file)
@@ -24,7 +24,7 @@
 #include <deal.II/base/quadrature_lib.h>
 #include <deal.II/base/work_stream.h>
 
-#include <deal.II/distributed/tria.h>
+#include <deal.II/distributed/tria_base.h>
 
 #include <deal.II/dofs/dof_accessor.h>
 #include <deal.II/dofs/dof_handler.h>
@@ -1160,30 +1160,23 @@ KellyErrorEstimator<dim, spacedim>::estimate(
   const types::material_id  material_id,
   const Strategy            strategy)
 {
-#ifdef DEAL_II_WITH_P4EST
-  if (dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim> *>(
-        &dof_handler.get_triangulation()) != nullptr)
-    Assert((subdomain_id_ == numbers::invalid_subdomain_id) ||
-             (subdomain_id_ ==
-              dynamic_cast<
-                const parallel::distributed::Triangulation<dim, spacedim> &>(
-                dof_handler.get_triangulation())
-                .locally_owned_subdomain()),
-           ExcMessage(
-             "For parallel distributed triangulations, the only "
-             "valid subdomain_id that can be passed here is the "
-             "one that corresponds to the locally owned subdomain id."));
-
-  const types::subdomain_id subdomain_id =
-    ((dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim> *>(
-        &dof_handler.get_triangulation()) != nullptr) ?
-       dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim>
-                      &>(dof_handler.get_triangulation())
-         .locally_owned_subdomain() :
-       subdomain_id_);
-#else
-  const types::subdomain_id subdomain_id = subdomain_id_;
-#endif
+  types::subdomain_id subdomain_id = numbers::invalid_subdomain_id;
+  if (const auto *triangulation = dynamic_cast<
+        const parallel::DistributedTriangulationBase<dim, spacedim> *>(
+        &dof_handler.get_triangulation()))
+    {
+      Assert((subdomain_id_ == numbers::invalid_subdomain_id) ||
+               (subdomain_id_ == triangulation->locally_owned_subdomain()),
+             ExcMessage(
+               "For distributed Triangulation objects and associated "
+               "DoFHandler objects, asking for any subdomain other than the "
+               "locally owned one does not make sense."));
+      subdomain_id = triangulation->locally_owned_subdomain();
+    }
+  else
+    {
+      subdomain_id = subdomain_id_;
+    }
 
   const unsigned int n_components = dof_handler.get_fe(0).n_components();
   (void)n_components;
index 887f7d9adfd09506578e0843f699df76355dabec..7f049f6e3be7274e3f9b958831a82f63b10d58fe 100644 (file)
@@ -75,18 +75,20 @@ namespace DoFTools
     Assert(sparsity.n_cols() == n_dofs,
            ExcDimensionMismatch(sparsity.n_cols(), n_dofs));
 
-    // If we have a distributed::Triangulation only allow locally_owned
+    // If we have a distributed Triangulation only allow locally_owned
     // subdomain. Not setting a subdomain is also okay, because we skip
     // ghost cells in the loop below.
-    Assert((dof.get_triangulation().locally_owned_subdomain() ==
-            numbers::invalid_subdomain_id) ||
-             (subdomain_id == numbers::invalid_subdomain_id) ||
-             (subdomain_id ==
-              dof.get_triangulation().locally_owned_subdomain()),
-           ExcMessage(
-             "For parallel::distributed::Triangulation objects and "
-             "associated DoF handler objects, asking for any subdomain other "
-             "than the locally owned one does not make sense."));
+    if (const auto *triangulation = dynamic_cast<
+          const parallel::DistributedTriangulationBase<dim, spacedim> *>(
+          &dof.get_triangulation()))
+      {
+        Assert((subdomain_id == numbers::invalid_subdomain_id) ||
+                 (subdomain_id == triangulation->locally_owned_subdomain()),
+               ExcMessage(
+                 "For distributed Triangulation objects and associated "
+                 "DoFHandler objects, asking for any subdomain other than the "
+                 "locally owned one does not make sense."));
+      }
 
     std::vector<types::global_dof_index> dofs_on_this_cell;
     dofs_on_this_cell.reserve(dof.get_fe_collection().max_dofs_per_cell());
@@ -143,18 +145,20 @@ namespace DoFTools
            ExcDimensionMismatch(couplings.n_cols(),
                                 dof.get_fe(0).n_components()));
 
-    // If we have a distributed::Triangulation only allow locally_owned
+    // If we have a distributed Triangulation only allow locally_owned
     // subdomain. Not setting a subdomain is also okay, because we skip
     // ghost cells in the loop below.
-    Assert((dof.get_triangulation().locally_owned_subdomain() ==
-            numbers::invalid_subdomain_id) ||
-             (subdomain_id == numbers::invalid_subdomain_id) ||
-             (subdomain_id ==
-              dof.get_triangulation().locally_owned_subdomain()),
-           ExcMessage(
-             "For parallel::distributed::Triangulation objects and "
-             "associated DoF handler objects, asking for any subdomain other "
-             "than the locally owned one does not make sense."));
+    if (const auto *triangulation = dynamic_cast<
+          const parallel::DistributedTriangulationBase<dim, spacedim> *>(
+          &dof.get_triangulation()))
+      {
+        Assert((subdomain_id == numbers::invalid_subdomain_id) ||
+                 (subdomain_id == triangulation->locally_owned_subdomain()),
+               ExcMessage(
+                 "For distributed Triangulation objects and associated "
+                 "DoFHandler objects, asking for any subdomain other than the "
+                 "locally owned one does not make sense."));
+      }
 
     const hp::FECollection<dim, spacedim> &fe_collection =
       dof.get_fe_collection();
@@ -541,18 +545,20 @@ namespace DoFTools
     AssertDimension(sparsity.n_rows(), n_dofs);
     AssertDimension(sparsity.n_cols(), n_dofs);
 
-    // If we have a distributed::Triangulation only allow locally_owned
+    // If we have a distributed Triangulation only allow locally_owned
     // subdomain. Not setting a subdomain is also okay, because we skip
     // ghost cells in the loop below.
-    Assert((dof.get_triangulation().locally_owned_subdomain() ==
-            numbers::invalid_subdomain_id) ||
-             (subdomain_id == numbers::invalid_subdomain_id) ||
-             (subdomain_id ==
-              dof.get_triangulation().locally_owned_subdomain()),
-           ExcMessage(
-             "For parallel::distributed::Triangulation objects and "
-             "associated DoF handler objects, asking for any subdomain other "
-             "than the locally owned one does not make sense."));
+    if (const auto *triangulation = dynamic_cast<
+          const parallel::DistributedTriangulationBase<dim, spacedim> *>(
+          &dof.get_triangulation()))
+      {
+        Assert((subdomain_id == numbers::invalid_subdomain_id) ||
+                 (subdomain_id == triangulation->locally_owned_subdomain()),
+               ExcMessage(
+                 "For distributed Triangulation objects and associated "
+                 "DoFHandler objects, asking for any subdomain other than the "
+                 "locally owned one does not make sense."));
+      }
 
     std::vector<types::global_dof_index> dofs_on_this_cell;
     std::vector<types::global_dof_index> dofs_on_other_cell;
@@ -1428,18 +1434,20 @@ namespace DoFTools
     Assert(flux_mask.n_cols() == n_comp,
            ExcDimensionMismatch(flux_mask.n_cols(), n_comp));
 
-    // If we have a distributed::Triangulation only allow locally_owned
+    // If we have a distributed Triangulation only allow locally_owned
     // subdomain. Not setting a subdomain is also okay, because we skip
     // ghost cells in the loop below.
-    Assert((dof.get_triangulation().locally_owned_subdomain() ==
-            numbers::invalid_subdomain_id) ||
-             (subdomain_id == numbers::invalid_subdomain_id) ||
-             (subdomain_id ==
-              dof.get_triangulation().locally_owned_subdomain()),
-           ExcMessage(
-             "For parallel::distributed::Triangulation objects and "
-             "associated DoF handler objects, asking for any subdomain other "
-             "than the locally owned one does not make sense."));
+    if (const auto *triangulation = dynamic_cast<
+          const parallel::DistributedTriangulationBase<dim, spacedim> *>(
+          &dof.get_triangulation()))
+      {
+        Assert((subdomain_id == numbers::invalid_subdomain_id) ||
+                 (subdomain_id == triangulation->locally_owned_subdomain()),
+               ExcMessage(
+                 "For distributed Triangulation objects and associated "
+                 "DoFHandler objects, asking for any subdomain other than the "
+                 "locally owned one does not make sense."));
+      }
 
     Assert(
       face_has_flux_coupling,
index 5a098cc9f664bc5cf6a717905fbb42acea9b0aec..8e11d0f6d8c8d15f733cf7447c1709eb7aabc909 100644 (file)
@@ -18,7 +18,7 @@
 #include <deal.II/base/quadrature_lib.h>
 #include <deal.II/base/work_stream.h>
 
-#include <deal.II/distributed/tria.h>
+#include <deal.II/distributed/tria_base.h>
 
 #include <deal.II/dofs/dof_accessor.h>
 #include <deal.II/dofs/dof_handler.h>
@@ -311,31 +311,24 @@ KellyErrorEstimator<1, spacedim>::estimate(
   const Strategy            strategy)
 {
   AssertThrow(strategy == cell_diameter_over_24, ExcNotImplemented());
-  using number = typename InputVector::value_type;
-#ifdef DEAL_II_WITH_P4EST
-  if (dynamic_cast<const parallel::distributed::Triangulation<1, spacedim> *>(
-        &dof_handler.get_triangulation()) != nullptr)
-    Assert((subdomain_id_ == numbers::invalid_subdomain_id) ||
-             (subdomain_id_ ==
-              dynamic_cast<
-                const parallel::distributed::Triangulation<1, spacedim> &>(
-                dof_handler.get_triangulation())
-                .locally_owned_subdomain()),
-           ExcMessage(
-             "For parallel distributed triangulations, the only "
-             "valid subdomain_id that can be passed here is the "
-             "one that corresponds to the locally owned subdomain id."));
-
-  const types::subdomain_id subdomain_id =
-    ((dynamic_cast<const parallel::distributed::Triangulation<1, spacedim> *>(
-        &dof_handler.get_triangulation()) != nullptr) ?
-       dynamic_cast<const parallel::distributed::Triangulation<1, spacedim> &>(
-         dof_handler.get_triangulation())
-         .locally_owned_subdomain() :
-       subdomain_id_);
-#else
-  const types::subdomain_id subdomain_id = subdomain_id_;
-#endif
+  using number                     = typename InputVector::value_type;
+  types::subdomain_id subdomain_id = numbers::invalid_subdomain_id;
+  if (const auto *triangulation = dynamic_cast<
+        const parallel::DistributedTriangulationBase<1, spacedim> *>(
+        &dof_handler.get_triangulation()))
+    {
+      Assert((subdomain_id_ == numbers::invalid_subdomain_id) ||
+               (subdomain_id_ == triangulation->locally_owned_subdomain()),
+             ExcMessage(
+               "For distributed Triangulation objects and associated "
+               "DoFHandler objects, asking for any subdomain other than the "
+               "locally owned one does not make sense."));
+      subdomain_id = triangulation->locally_owned_subdomain();
+    }
+  else
+    {
+      subdomain_id = subdomain_id_;
+    }
 
   const unsigned int n_components       = dof_handler.get_fe(0).n_components();
   const unsigned int n_solution_vectors = solutions.size();

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.