]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Remove DataOut::first_cell() and DataOut::next_cell() 12382/head
authorDaniel Arndt <arndtd@ornl.gov>
Wed, 2 Jun 2021 15:56:17 +0000 (11:56 -0400)
committerDaniel Arndt <arndtd@ornl.gov>
Wed, 2 Jun 2021 17:10:54 +0000 (13:10 -0400)
doc/news/changes/incompatibilities/20210602DanielArndt-1 [new file with mode: 0644]
include/deal.II/numerics/data_out.h
source/numerics/data_out.cc
tests/data_out/data_out_08.cc [deleted file]
tests/data_out/data_out_08.output [deleted file]
tests/mpi/p4est_2d_constraintmatrix_04.cc
tests/mpi/p4est_3d_constraintmatrix_03.cc

diff --git a/doc/news/changes/incompatibilities/20210602DanielArndt-1 b/doc/news/changes/incompatibilities/20210602DanielArndt-1
new file mode 100644 (file)
index 0000000..b6999e2
--- /dev/null
@@ -0,0 +1,4 @@
+Removed: The deprecated member functions DataOut::first_cell() and DataOut::next_cell()
+have been removed.
+<br>
+(Daniel Arndt, 2021/06/02)
index cdd7a1471767db995285dc60b71b0a54389d94b6..39c20cd0ff4d09de6fd84eedd17c83315253f0d8 100644 (file)
@@ -118,16 +118,16 @@ namespace internal
  * small to be seen individually) or because you only want to see a certain
  * region of the domain (for example only in the fluid part of the domain in
  * step-46), or for some other reason.
- *
- * For this, internally build_patches() does not generate the sequence of
- * cells to be converted into patches itself, but relies on the two function
- * that we'll call first_cell() and next_cell(). By default, they return the
- * first active cell, and the next active cell, respectively. But this can
- * be changed using the set_cell_selection() function that allows you to
- * replace this behavior. What set_cell_selection() wants to know is how
- * you want to pick out the first cell on which output should be generated,
- * and how given one cell on which output is generated you want to pick the
- * next cell.
+
+ * For this, internally build_patches() does not generate the sequence of cells
+ * to be converted into patches itself, but relies on the two private
+ * std::function objects first_cell_function() and next_cell_function(). By
+ * default, they return the first active cell, and the next active cell,
+ * respectively. But this can be changed using the set_cell_selection() function
+ * that allows you to replace this behavior. What set_cell_selection() wants to
+ * know is how you want to pick out the first cell on which output should be
+ * generated, and how given one cell on which output is generated you want to
+ * pick the next cell.
  *
  * This may,
  * for example, include only cells that are in parts of a domain (e.g., if you
@@ -416,35 +416,6 @@ public:
   const std::pair<FirstCellFunctionType, NextCellFunctionType>
   get_cell_selection() const;
 
-  /**
-   * Return the first cell which we want output for. The default
-   * implementation returns the first active cell, but you might want to
-   * return other cells in a derived class.
-   *
-   * @deprecated Use the set_cell_selection() function instead.
-   */
-  DEAL_II_DEPRECATED
-  virtual cell_iterator
-  first_cell();
-
-  /**
-   * Return the next cell after @p cell which we want output for.  If there
-   * are no more cells, any implementation of this function should return
-   * <tt>dof_handler->end()</tt>.
-   *
-   * The default implementation returns the next active cell, but you might
-   * want to return other cells in a derived class. Note that the default
-   * implementation assumes that the given @p cell is active, which is
-   * guaranteed as long as first_cell() is also used from the default
-   * implementation. Overloading only one of the two functions might not be a
-   * good idea.
-   *
-   * @deprecated Use the set_cell_selection() function instead.
-   */
-  DEAL_II_DEPRECATED
-  virtual cell_iterator
-  next_cell(const cell_iterator &cell);
-
 private:
   /**
    * A function object that is used to select what the first cell is going to
@@ -463,28 +434,6 @@ private:
                               const cell_iterator &)>
     next_cell_function;
 
-  /**
-   * Return the first cell produced by the first_cell()/next_cell() function
-   * pair that is locally owned. If this object operates on a non-distributed
-   * triangulation, the result equals what first_cell() returns.
-   *
-   * @deprecated Use the set_cell_selection() function instead.
-   */
-  DEAL_II_DEPRECATED
-  virtual cell_iterator
-  first_locally_owned_cell();
-
-  /**
-   * Return the next cell produced by the next_cell() function that is locally
-   * owned. If this object operates on a non-distributed triangulation, the
-   * result equals what first_cell() returns.
-   *
-   * @deprecated Use the set_cell_selection() function instead.
-   */
-  DEAL_II_DEPRECATED
-  virtual cell_iterator
-  next_locally_owned_cell(const cell_iterator &cell);
-
   /**
    * Build one patch. This function is called in a WorkStream context.
    *
index 2507c7d819b7980d60fcc398d4d9420524f5d839..11bb33f4086144287a0b109276aa5c89334f8904 100644 (file)
@@ -68,16 +68,25 @@ namespace internal
 template <int dim, int spacedim>
 DataOut<dim, spacedim>::DataOut()
 {
-  // For the moment, just call the existing virtual functions. This
-  // preserves backward compatibility. When these deprecated functions are
-  // removed, we'll just inline their functionality into the lambda
-  // functions created here.
   set_cell_selection(
     [this](const Triangulation<dim, spacedim> &) {
-      return this->first_locally_owned_cell();
+      typename Triangulation<dim, spacedim>::active_cell_iterator cell =
+        this->triangulation->begin_active();
+
+      // skip cells if the current one has no children (is active) and is a
+      // ghost or artificial cell
+      while ((cell != this->triangulation->end()) && !cell->is_locally_owned())
+        ++cell;
+      return cell;
     },
-    [this](const Triangulation<dim, spacedim> &, const cell_iterator &cell) {
-      return this->next_locally_owned_cell(cell);
+    [this](const Triangulation<dim, spacedim> &,
+           const cell_iterator &old_cell) {
+      typename Triangulation<dim, spacedim>::active_cell_iterator cell =
+        old_cell;
+      ++cell;
+      while ((cell != this->triangulation->end()) && !cell->is_locally_owned())
+        ++cell;
+      return cell;
     });
 }
 
@@ -1333,60 +1342,6 @@ DataOut<dim, spacedim>::get_cell_selection() const
 
 
 
-template <int dim, int spacedim>
-typename DataOut<dim, spacedim>::cell_iterator
-DataOut<dim, spacedim>::first_cell()
-{
-  return this->triangulation->begin_active();
-}
-
-
-
-template <int dim, int spacedim>
-typename DataOut<dim, spacedim>::cell_iterator
-DataOut<dim, spacedim>::next_cell(
-  const typename DataOut<dim, spacedim>::cell_iterator &cell)
-{
-  // convert the iterator to an active_iterator and advance this to the next
-  // active cell
-  typename Triangulation<dim, spacedim>::active_cell_iterator active_cell =
-    cell;
-  ++active_cell;
-  return active_cell;
-}
-
-
-
-template <int dim, int spacedim>
-typename DataOut<dim, spacedim>::cell_iterator
-DataOut<dim, spacedim>::first_locally_owned_cell()
-{
-  typename DataOut<dim, spacedim>::cell_iterator cell = first_cell();
-
-  // skip cells if the current one has no children (is active) and is a ghost
-  // or artificial cell
-  while ((cell != this->triangulation->end()) && cell->is_active() &&
-         !cell->is_locally_owned())
-    cell = next_cell(cell);
-
-  return cell;
-}
-
-
-
-template <int dim, int spacedim>
-typename DataOut<dim, spacedim>::cell_iterator
-DataOut<dim, spacedim>::next_locally_owned_cell(
-  const typename DataOut<dim, spacedim>::cell_iterator &old_cell)
-{
-  typename DataOut<dim, spacedim>::cell_iterator cell = next_cell(old_cell);
-  while ((cell != this->triangulation->end()) && cell->is_active() &&
-         !cell->is_locally_owned())
-    cell = next_cell(cell);
-  return cell;
-}
-
-
 // explicit instantiations
 #include "data_out.inst"
 
diff --git a/tests/data_out/data_out_08.cc b/tests/data_out/data_out_08.cc
deleted file mode 100644 (file)
index eec6a3e..0000000
+++ /dev/null
@@ -1,145 +0,0 @@
-// ---------------------------------------------------------------------
-//
-// Copyright (C) 2003 - 2020 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.
-//
-// ---------------------------------------------------------------------
-
-
-// This test documents two unrelated bugs in DataOut when used with a Filter (by
-// deriving from DataOut):
-// 1. The patch index computation in data_out.cc is wrong and causes an SIGV (or
-// an Assert after adding that):
-/*
-466: --------------------------------------------------------
-466: An error occurred in line <306> of file
-</ssd/branch_port_the_testsuite/deal.II/source/numerics/data_out.cc> in function
-466:     void dealii::DataOut<dim, DoFHandlerType>::build_one_patch(const
-std::pair<typename dealii::DataOut_DoFData<DoFHandlerType, DoFHandlerType::
-dimension, DoFHandlerType:: space_dimension>::cell_iterator, unsigned int>*,
-dealii::internal::DataOut::ParallelData<DoFHandlerType:: dimension,
-DoFHandlerType:: space_dimension>&, dealii::DataOutBase::Patch<DoFHandlerType::
-dimension, DoFHandlerType:: space_dimension>&, dealii::DataOut<dim,
-DoFHandlerType>::CurvedCellRegion,
-std::vector<dealii::DataOutBase::Patch<DoFHandlerType:: dimension,
-DoFHandlerType:: space_dimension> >&) [with int dim = 2, DoFHandlerType =
-dealii::DoFHandler<2>, typename dealii::DataOut_DoFData<DoFHandlerType,
-DoFHandlerType:: dimension, DoFHandlerType:: space_dimension>::cell_iterator =
-dealii::TriaIterator<dealii::CellAccessor<2, 2> >] 466: The violated condition
-was: 466:     cell_and_index->second < patches.size() 466: The name and call
-sequence of the exception was: 466:     ExcInternalError()
-*/
-// 2. DataOut used begin_active() instead of first_cell() in two places which
-// caused a wrong patch to be generated when the first active cell is not picked
-// by the filter.
-
-#include <deal.II/dofs/dof_accessor.h>
-#include <deal.II/dofs/dof_handler.h>
-#include <deal.II/dofs/dof_tools.h>
-
-#include <deal.II/fe/fe_dgq.h>
-
-#include <deal.II/grid/filtered_iterator.h>
-#include <deal.II/grid/grid_generator.h>
-#include <deal.II/grid/tria.h>
-#include <deal.II/grid/tria_iterator.h>
-
-#include <deal.II/lac/vector.h>
-
-#include <deal.II/numerics/data_out.h>
-
-#include "../tests.h"
-
-
-template <int dim>
-class FilteredDataOut : public DataOut<dim>
-{
-public:
-  FilteredDataOut(const unsigned int subdomain_id)
-    : subdomain_id(subdomain_id)
-  {}
-
-  virtual typename DataOut<dim>::cell_iterator
-  first_cell()
-  {
-    auto cell = this->dofs->begin_active();
-    while ((cell != this->dofs->end()) &&
-           (cell->subdomain_id() != subdomain_id))
-      ++cell;
-
-    return cell;
-  }
-
-  virtual typename DataOut<dim>::cell_iterator
-  next_cell(const typename DataOut<dim>::cell_iterator &old_cell)
-  {
-    if (old_cell != this->dofs->end())
-      {
-        const IteratorFilters::SubdomainEqualTo predicate(subdomain_id);
-
-        return ++(
-          FilteredIterator<typename Triangulation<dim>::active_cell_iterator>(
-            predicate, old_cell));
-      }
-    else
-      return old_cell;
-  }
-
-private:
-  const unsigned int subdomain_id;
-};
-
-
-template <int dim>
-void
-check()
-{
-  Triangulation<dim> tria;
-  GridGenerator::hyper_cube(tria, 0., 1.);
-  tria.refine_global(1);
-
-  Vector<double> cell_data(4);
-  for (unsigned int i = 0; i < 4; ++i)
-    cell_data(i) = i * 1.0;
-
-  // this should skip the first cell
-  typename Triangulation<dim>::active_cell_iterator it = tria.begin_active();
-  it->set_subdomain_id(1);
-
-  FE_DGQ<dim>     fe(0);
-  DoFHandler<dim> dof_handler(tria);
-  dof_handler.distribute_dofs(fe);
-
-  // we pick only subdomain==0 which will
-  // skip the first of the four cells
-  FilteredDataOut<dim> data_out(0);
-  data_out.attach_dof_handler(dof_handler);
-
-  data_out.add_data_vector(cell_data,
-                           "cell_data",
-                           DataOut<dim>::type_cell_data);
-  data_out.build_patches();
-
-  data_out.write_deal_II_intermediate(deallog.get_file_stream());
-
-  deallog << "OK" << std::endl;
-}
-
-int
-main(int argc, char **argv)
-{
-  initlog();
-  deallog << std::setprecision(2);
-
-  check<2>();
-
-  return 0;
-}
diff --git a/tests/data_out/data_out_08.output b/tests/data_out/data_out_08.output
deleted file mode 100644 (file)
index bd4f3f7..0000000
+++ /dev/null
@@ -1,38 +0,0 @@
-
-2 2
-[deal.II intermediate format graphics data]
-[written by deal.II x.y.z]
-[Version: 3]
-1
-cell_data
-3
-[deal.II intermediate Patch<2,2>]
-0.50 0.0 1.0 0.0 1.0 0.50 0.50 0.50 
-4294967295 4294967295 4294967295 2 
-0 1
-0
-1 4
-1.0 1.0 1.0 1.0 
-
-
-[deal.II intermediate Patch<2,2>]
-0.0 0.50 0.50 0.50 0.50 1.0 0.0 1.0 
-4294967295 2 4294967295 4294967295 
-1 1
-0
-1 4
-2.0 2.0 2.0 2.0 
-
-
-[deal.II intermediate Patch<2,2>]
-0.50 0.50 1.0 0.50 1.0 1.0 0.50 1.0 
-1 4294967295 0 4294967295 
-2 1
-0
-1 4
-3.0 3.0 3.0 3.0 
-
-
-0
-
-DEAL::OK
index 63db8eac75afb086847d3220d594b448f81e63c0..e9c99c3b664eecb4a02936b89b5e1daad7a81ce5 100644 (file)
@@ -56,44 +56,6 @@ const double T1 = 2.0;
 
 
 
-template <int dim>
-class FilteredDataOut : public DataOut<dim>
-{
-public:
-  FilteredDataOut(const unsigned int subdomain_id)
-    : subdomain_id(subdomain_id)
-  {}
-
-  virtual typename DataOut<dim>::cell_iterator
-  first_cell()
-  {
-    auto cell = this->triangulation->begin_active();
-    while ((cell != this->triangulation->end()) &&
-           (cell->subdomain_id() != subdomain_id))
-      ++cell;
-
-    return cell;
-  }
-
-  virtual typename DataOut<dim>::cell_iterator
-  next_cell(const typename DataOut<dim>::cell_iterator &old_cell)
-  {
-    if (old_cell != this->triangulation->end())
-      {
-        const IteratorFilters::SubdomainEqualTo predicate(subdomain_id);
-
-        return ++(
-          FilteredIterator<typename Triangulation<dim>::active_cell_iterator>(
-            predicate, old_cell));
-      }
-    else
-      return old_cell;
-  }
-
-private:
-  const unsigned int subdomain_id;
-};
-
 template <int dim>
 class TemperatureInitialValues : public Function<dim>
 {
@@ -287,7 +249,32 @@ test()
 
   std::vector<std::string> solution_names(1, "T");
 
-  FilteredDataOut<dim> data_out(tr.locally_owned_subdomain());
+  DataOut<dim> data_out;
+  data_out.set_cell_selection(
+    [](const Triangulation<dim> &t) {
+      auto cell = t.begin_active();
+      while ((cell != t.end()) &&
+             (cell->subdomain_id() != t.locally_owned_subdomain()))
+        ++cell;
+
+      return cell;
+    },
+
+    [](const Triangulation<dim> &                        t,
+       const typename Triangulation<dim>::cell_iterator &old_cell) ->
+    typename Triangulation<dim>::cell_iterator {
+      if (old_cell != t.end())
+        {
+          const IteratorFilters::SubdomainEqualTo predicate(
+            t.locally_owned_subdomain());
+
+          return ++(
+            FilteredIterator<typename Triangulation<dim>::active_cell_iterator>(
+              predicate, old_cell));
+        }
+      else
+        return old_cell;
+    });
   data_out.attach_dof_handler(dofh);
 
   data_out.add_data_vector(x_rel, solution_names);
index 2ad23eddef6b25f18406726575d9546295b6bd11..e715ff158a02e4017dbbf537a43df7466fb53ee2 100644 (file)
@@ -56,44 +56,6 @@ const double T1 = 2.0;
 
 
 
-template <int dim>
-class FilteredDataOut : public DataOut<dim>
-{
-public:
-  FilteredDataOut(const unsigned int subdomain_id)
-    : subdomain_id(subdomain_id)
-  {}
-
-  virtual typename DataOut<dim>::cell_iterator
-  first_cell()
-  {
-    auto cell = this->triangulation->begin_active();
-    while ((cell != this->triangulation->end()) &&
-           (cell->subdomain_id() != subdomain_id))
-      ++cell;
-
-    return cell;
-  }
-
-  virtual typename DataOut<dim>::cell_iterator
-  next_cell(const typename DataOut<dim>::cell_iterator &old_cell)
-  {
-    if (old_cell != this->triangulation->end())
-      {
-        const IteratorFilters::SubdomainEqualTo predicate(subdomain_id);
-
-        return ++(
-          FilteredIterator<typename Triangulation<dim>::active_cell_iterator>(
-            predicate, old_cell));
-      }
-    else
-      return old_cell;
-  }
-
-private:
-  const unsigned int subdomain_id;
-};
-
 template <int dim>
 class TemperatureInitialValues : public Function<dim>
 {
@@ -287,7 +249,32 @@ test()
 
   std::vector<std::string> solution_names(1, "T");
 
-  FilteredDataOut<dim> data_out(tr.locally_owned_subdomain());
+  DataOut<dim> data_out;
+  data_out.set_cell_selection(
+    [](const Triangulation<dim> &t) {
+      auto cell = t.begin_active();
+      while ((cell != t.end()) &&
+             (cell->subdomain_id() != t.locally_owned_subdomain()))
+        ++cell;
+
+      return cell;
+    },
+
+    [](const Triangulation<dim> &                        t,
+       const typename Triangulation<dim>::cell_iterator &old_cell) ->
+    typename Triangulation<dim>::cell_iterator {
+      if (old_cell != t.end())
+        {
+          const IteratorFilters::SubdomainEqualTo predicate(
+            t.locally_owned_subdomain());
+
+          return ++(
+            FilteredIterator<typename Triangulation<dim>::active_cell_iterator>(
+              predicate, old_cell));
+        }
+      else
+        return old_cell;
+    });
   data_out.attach_dof_handler(dofh);
 
   data_out.add_data_vector(x_rel, solution_names);

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.