]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Improve calculating boolean values 9936/head
authorReza Rastak <rastak@stanford.edu>
Wed, 22 Apr 2020 02:49:43 +0000 (19:49 -0700)
committerReza Rastak <rastak@stanford.edu>
Fri, 24 Apr 2020 22:04:50 +0000 (15:04 -0700)
12 files changed:
include/deal.II/fe/fe_tools.templates.h
include/deal.II/hp/fe_collection.h
include/deal.II/lac/affine_constraints.templates.h
include/deal.II/numerics/data_out_dof_data.templates.h
source/base/parameter_handler.cc
source/base/quadrature.cc
source/base/quadrature_lib.cc
source/distributed/tria.cc
source/dofs/dof_handler_policy.cc
source/dofs/dof_tools.cc
source/lac/trilinos_sparse_matrix.cc
source/multigrid/mg_transfer_internal.cc

index a20409b117aaf7f2292858b1c18e8c921aa38350..b631e1bd60148287eb99e60986c8fd231f766431 100644 (file)
@@ -1459,23 +1459,18 @@ namespace FETools
                                       fe2.dofs_per_cell,
                                       fe1.dofs_per_cell));
 
-    // first try the easy way: maybe
-    // the FE wants to implement things
-    // itself:
-    bool fe_implements_interpolation = true;
+    // first try the easy way: maybe the FE wants to implement things itself:
     try
       {
         internal::FEToolsGetInterpolationMatrixHelper::gim_forwarder(
           fe1, fe2, interpolation_matrix);
+        return;
       }
     catch (
       typename FiniteElement<dim, spacedim>::ExcInterpolationNotImplemented &)
       {
         // too bad....
-        fe_implements_interpolation = false;
       }
-    if (fe_implements_interpolation == true)
-      return;
 
     // uh, so this was not the
     // case. hm. then do it the hard
index 8b78569b85934ffc0162ff61ee06dff80127738f..bf06f99033af9688aa359fb4badac9c6eefbf641 100644 (file)
@@ -1006,13 +1006,12 @@ namespace hp
   FECollection<dim, spacedim>::hp_constraints_are_implemented() const
   {
     Assert(finite_elements.size() > 0, ExcNoFiniteElements());
-
-    bool hp_constraints = true;
-    for (unsigned int i = 0; i < finite_elements.size(); ++i)
-      hp_constraints =
-        hp_constraints && finite_elements[i]->hp_constraints_are_implemented();
-
-    return hp_constraints;
+    return std::all_of(
+      finite_elements.cbegin(),
+      finite_elements.cend(),
+      [](const std::shared_ptr<const FiniteElement<dim, spacedim>> &fe) {
+        return fe->hp_constraints_are_implemented();
+      });
   }
 
 
index 5712424ca784d98a95aa6e5598537bba1c40962e..f6ffb99b6665a118c2dec153ead122393cef540c 100644 (file)
@@ -4013,23 +4013,21 @@ AffineConstraints<number>::add_entries_local_to_global(
   Assert(sparsity_pattern.n_rows() == sparsity_pattern.n_cols(),
          ExcNotQuadratic());
 
-  const size_type n_local_dofs       = local_dof_indices.size();
-  bool            dof_mask_is_active = false;
-  if (dof_mask.n_rows() == n_local_dofs)
-    {
-      dof_mask_is_active = true;
-      AssertDimension(dof_mask.n_cols(), n_local_dofs);
-    }
-
+  const size_type n_local_dofs = local_dof_indices.size();
   typename internals::AffineConstraintsData<number>::ScratchDataAccessor
     scratch_data;
 
-  // if the dof mask is not active, all we have to do is to add some indices
-  // in a matrix format. To do this, we first create an array of all the
-  // indices that are to be added. these indices are the local dof indices
-  // plus some indices that come from constraints.
-  if (dof_mask_is_active == false)
+  const bool dof_mask_is_active = (dof_mask.n_rows() == n_local_dofs);
+  if (dof_mask_is_active == true)
+    {
+      AssertDimension(dof_mask.n_cols(), n_local_dofs);
+    }
+  else
     {
+      // if the dof mask is not active, all we have to do is to add some indices
+      // in a matrix format. To do this, we first create an array of all the
+      // indices that are to be added. these indices are the local dof indices
+      // plus some indices that come from constraints.
       std::vector<size_type> &actual_dof_indices = scratch_data->columns;
       actual_dof_indices.resize(n_local_dofs);
       make_sorted_row_list(local_dof_indices, actual_dof_indices);
@@ -4107,11 +4105,8 @@ AffineConstraints<number>::add_entries_local_to_global(
   const bool                    keep_constrained_entries,
   const Table<2, bool> &        dof_mask) const
 {
-  const size_type n_local_rows       = row_indices.size();
-  const size_type n_local_cols       = col_indices.size();
-  bool            dof_mask_is_active = false;
-  if (dof_mask.n_rows() == n_local_rows && dof_mask.n_cols() == n_local_cols)
-    dof_mask_is_active = true;
+  const size_type n_local_rows = row_indices.size();
+  const size_type n_local_cols = col_indices.size();
 
   // if constrained entries should be kept, need to add rows and columns of
   // those to the sparsity pattern
@@ -4131,6 +4126,8 @@ AffineConstraints<number>::add_entries_local_to_global(
   // in a matrix format. To do this, we first create an array of all the
   // indices that are to be added. these indices are the local dof indices
   // plus some indices that come from constraints.
+  const bool dof_mask_is_active =
+    dof_mask.n_rows() == n_local_rows && dof_mask.n_cols() == n_local_cols;
   if (dof_mask_is_active == false)
     {
       std::vector<size_type> actual_row_indices(n_local_rows);
@@ -4176,14 +4173,12 @@ AffineConstraints<number>::add_entries_local_to_global(
   typename internals::AffineConstraintsData<number>::ScratchDataAccessor
     scratch_data;
 
-  bool dof_mask_is_active = false;
-  if (dof_mask.n_rows() == n_local_dofs)
+  const bool dof_mask_is_active = (dof_mask.n_rows() == n_local_dofs);
+  if (dof_mask_is_active == true)
     {
-      dof_mask_is_active = true;
       AssertDimension(dof_mask.n_cols(), n_local_dofs);
     }
-
-  if (dof_mask_is_active == false)
+  else
     {
       std::vector<size_type> &actual_dof_indices = scratch_data->columns;
       actual_dof_indices.resize(n_local_dofs);
index ea183e48956ee48f31626070527b15368e2329e9..062fd6c7d91923fe50511ade97ae563242addd94 100644 (file)
@@ -227,15 +227,12 @@ namespace internal
     {
       for (unsigned int dataset = 0; dataset < dof_data.size(); ++dataset)
         {
-          bool duplicate = false;
-          for (unsigned int j = 0; j < dataset; ++j)
-            if (finite_elements[dataset].get() == finite_elements[j].get())
-              {
-                duplicate = true;
-                break;
-              }
-
-          if (duplicate == false)
+          const bool is_duplicate = std::any_of(
+            finite_elements.cbegin(),
+            finite_elements.cbegin() + dataset,
+            [&](const std::shared_ptr<dealii::hp::FECollection<dim, spacedim>>
+                  &fe) { return finite_elements[dataset].get() == fe.get(); });
+          if (is_duplicate == false)
             {
               if (cell->active())
                 {
index 0c98baea71255d593e9d3e8587959ad3025a1dc0..2f3f1576078278bf7063c40d2aa398d5f8ebaf4b 100644 (file)
@@ -1440,15 +1440,13 @@ ParameterHandler::recursively_print_parameters(
 
       // if there are any parameters in this section then print them
       // as an itemized list
-      bool parameters_exist_here = false;
-      for (const auto &p : current_section)
-        if ((is_parameter_node(p.second) == true) ||
-            (is_alias_node(p.second) == true))
-          {
-            parameters_exist_here = true;
-            break;
-          }
-
+      const bool parameters_exist_here =
+        std::any_of(current_section.begin(),
+                    current_section.end(),
+                    [](const boost::property_tree::ptree::value_type &p) {
+                      return is_parameter_node(p.second) ||
+                             is_alias_node(p.second);
+                    });
       if (parameters_exist_here)
         {
           out << "\\begin{itemize}" << '\n';
index 389f4a29c89e56c7c25e1fa36848a2172c98555c..bcfc43bea093da9774b9211d2549057cff4c3d18 100644 (file)
@@ -1536,15 +1536,14 @@ namespace internal
       bool
       uses_both_endpoints(const Quadrature<1> &base_quadrature)
       {
-        bool at_left = false, at_right = false;
-        for (unsigned int i = 0; i < base_quadrature.size(); ++i)
-          {
-            if (base_quadrature.point(i) == Point<1>(0.0))
-              at_left = true;
-            if (base_quadrature.point(i) == Point<1>(1.0))
-              at_right = true;
-          }
-
+        const bool at_left =
+          std::any_of(base_quadrature.get_points().cbegin(),
+                      base_quadrature.get_points().cend(),
+                      [](const Point<1> &p) { return p == Point<1>{0.}; });
+        const bool at_right =
+          std::any_of(base_quadrature.get_points().cbegin(),
+                      base_quadrature.get_points().cend(),
+                      [](const Point<1> &p) { return p == Point<1>{1.}; });
         return (at_left && at_right);
       }
     } // namespace
index 680486bb08f0fc4ba412bc2b80253f320a6418f2..cf7eb51b8a82984c73f1eb0c41b03b898c6f9032 100644 (file)
@@ -608,21 +608,22 @@ template <>
 unsigned int
 QGaussOneOverR<2>::quad_size(const Point<2> singularity, const unsigned int n)
 {
-  double eps       = 1e-8;
-  bool   on_edge   = false;
-  bool   on_vertex = false;
-  for (unsigned int i = 0; i < 2; ++i)
-    if ((std::abs(singularity[i]) < eps) ||
-        (std::abs(singularity[i] - 1) < eps))
-      on_edge = true;
-  if (on_edge &&
-      (std::abs((singularity - Point<2>(.5, .5)).norm_square() - .5) < eps))
-    on_vertex = true;
+  const double eps = 1e-8;
+  const bool   on_edge =
+    std::any_of(singularity.begin_raw(),
+                singularity.end_raw(),
+                [eps](double coord) {
+                  return std::abs(coord) < eps || std::abs(coord - 1.) < eps;
+                });
+  const bool on_vertex =
+    on_edge &&
+    std::abs((singularity - Point<2>(.5, .5)).norm_square() - .5) < eps;
   if (on_vertex)
-    return (2 * n * n);
-  if (on_edge)
-    return (4 * n * n);
-  return (8 * n * n);
+    return 2 * n * n;
+  else if (on_edge)
+    return 4 * n * n;
+  else
+    return 8 * n * n;
 }
 
 template <>
index 8aa5670873f3aa06a8fefb9ab80818502f2b395f..28ed1d7f3fc6d4dc92e235cabdd200408b8ec0e3 100644 (file)
@@ -2622,16 +2622,12 @@ namespace parallel
       //
       // The problem is that we cannot just ask for the first active cell, but
       // instead need to filter over locally owned cells.
-      bool have_coarser_cell = false;
-      for (typename Triangulation<dim, spacedim>::active_cell_iterator cell =
-             this->begin_active(this->n_global_levels() - 2);
-           cell != this->end(this->n_global_levels() - 2);
-           ++cell)
-        if (cell->is_locally_owned())
-          {
-            have_coarser_cell = true;
-            break;
-          }
+      const bool have_coarser_cell =
+        std::any_of(this->begin_active(this->n_global_levels() - 2),
+                    this->end_active(this->n_global_levels() - 2),
+                    [](const CellAccessor<dim, spacedim> &cell) {
+                      return cell.is_locally_owned();
+                    });
 
       // return true if at least one process has a coarser cell
       return 0 < Utilities::MPI::max(have_coarser_cell ? 1 : 0,
@@ -3632,13 +3628,13 @@ namespace parallel
           this->prepare_coarsening_and_refinement();
 
           // see if any flags are still set
-          mesh_changed = false;
-          for (const auto &cell : this->active_cell_iterators())
-            if (cell->refine_flag_set() || cell->coarsen_flag_set())
-              {
-                mesh_changed = true;
-                break;
-              }
+          mesh_changed =
+            std::any_of(this->begin_active(),
+                        active_cell_iterator{this->end()},
+                        [](const CellAccessor<dim, spacedim> &cell) {
+                          return cell.refine_flag_set() ||
+                                 cell.coarsen_flag_set();
+                        });
 
           // actually do the refinement to change the local mesh by
           // calling the base class refinement function directly
index a21950148c8056ee13069de0db342efa5a525fa6..ff7d43120c5e0c2a19c04d29cd4badc850d47d50 100644 (file)
@@ -4909,13 +4909,13 @@ namespace internal
         // ranks changed. In that case, we can apply the renumbering with some
         // local renumbering only (this is similar to the renumber_mg_dofs()
         // function below)
-        bool locally_owned_set_changes = false;
-        for (types::global_dof_index i : new_numbers)
-          if (dof_handler->locally_owned_dofs().is_element(i) == false)
-            {
-              locally_owned_set_changes = true;
-              break;
-            }
+        const bool locally_owned_set_changes =
+          std::any_of(new_numbers.cbegin(),
+                      new_numbers.cend(),
+                      [this](const types::global_dof_index i) {
+                        return dof_handler->locally_owned_dofs().is_element(
+                                 i) == false;
+                      });
 
         if (Utilities::MPI::sum(static_cast<unsigned int>(
                                   locally_owned_set_changes),
index 05b451da9395c09e22a0078908d5831bac77a13c..63e5dc05188721f3eb2a39c7634bd114c464a2f3 100644 (file)
@@ -1729,26 +1729,14 @@ namespace DoFTools
                                 dof_handler.get_fe(0).n_components()));
     std::fill(n_dofs_on_subdomain.begin(), n_dofs_on_subdomain.end(), 0);
 
-    // in debug mode, make sure that there are some cells at least with
-    // this subdomain id
-#ifdef DEBUG
-    {
-      bool found = false;
-      for (typename Triangulation<
-             DoFHandlerType::dimension,
-             DoFHandlerType::space_dimension>::active_cell_iterator cell =
-             dof_handler.get_triangulation().begin_active();
-           cell != dof_handler.get_triangulation().end();
-           ++cell)
-        if (cell->subdomain_id() == subdomain)
-          {
-            found = true;
-            break;
-          }
-      Assert(found == true,
-             ExcMessage("There are no cells for the given subdomain!"));
-    }
-#endif
+    // Make sure there are at least some cells with this subdomain id
+    Assert(std::any_of(
+             dof_handler.begin_active(),
+             typename DoFHandlerType::active_cell_iterator{dof_handler.end()},
+             [subdomain](const typename DoFHandlerType::cell_accessor &cell) {
+               return cell.subdomain_id() == subdomain;
+             }),
+           ExcMessage("There are no cells for the given subdomain!"));
 
     std::vector<types::subdomain_id> subdomain_association(
       dof_handler.n_dofs());
index 2d8031a1b05077cb9f7010eb0bc6f09a7c6ced69..014e3c6b84c32c5a32b0c1f770d966d42568acc3 100644 (file)
@@ -648,8 +648,7 @@ namespace TrilinosWrappers
       // check whether the relevant rows correspond to exactly the same map as
       // the owned rows. In that case, do not create the nonlocal graph and
       // fill the columns by demand
-      bool have_ghost_rows = false;
-      {
+      const bool have_ghost_rows = [&]() {
         std::vector<dealii::types::global_dof_index> indices;
         relevant_rows.fill_index_vector(indices);
         Epetra_Map relevant_map(
@@ -661,11 +660,8 @@ namespace TrilinosWrappers
                indices.data())),
           0,
           row_space_map.Comm());
-        if (relevant_map.SameAs(row_space_map))
-          have_ghost_rows = false;
-        else
-          have_ghost_rows = true;
-      }
+        return !relevant_map.SameAs(row_space_map);
+      }();
 
       const unsigned int n_rows = relevant_rows.n_elements();
       std::vector<TrilinosWrappers::types::int_type> ghost_rows;
index 20179c8eb34130c64be1a4f7300453e79253e87c..4ba4c140d0b2df7be3ef6d2eeaa6336ec1762d7f 100644 (file)
@@ -764,16 +764,15 @@ namespace internal
               if (!cell->has_children())
                 continue;
 
-              bool consider_cell = false;
-              if (tria.locally_owned_subdomain() ==
-                    numbers::invalid_subdomain_id ||
-                  cell->level_subdomain_id() == tria.locally_owned_subdomain())
-                consider_cell = true;
+              bool consider_cell =
+                (tria.locally_owned_subdomain() ==
+                   numbers::invalid_subdomain_id ||
+                 cell->level_subdomain_id() == tria.locally_owned_subdomain());
 
               // due to the particular way we store DoF indices (via children),
               // we also need to add the DoF indices for coarse cells where we
               // own at least one child
-              bool cell_is_remote = !consider_cell;
+              const bool cell_is_remote = !consider_cell;
               for (unsigned int c = 0;
                    c < GeometryInfo<dim>::max_children_per_cell;
                    ++c)

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.