]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added move constructor to Triangulation 2769/head
authordanshapero <shapero.daniel@gmail.com>
Tue, 9 Aug 2016 18:00:14 +0000 (11:00 -0700)
committerdanshapero <shapero.daniel@gmail.com>
Tue, 9 Aug 2016 18:00:14 +0000 (11:00 -0700)
doc/news/changes.h
include/deal.II/grid/tria.h
source/grid/tria.cc
tests/grid/move_constructor.cc [new file with mode: 0644]
tests/grid/move_constructor.with_cxx11=on.output [new file with mode: 0644]

index e5e4c546f6b37cd00e8dd3e013d32d066c89735e..e892bfc8d2ed7879e415d5eebacbfdc6f46db29e 100644 (file)
@@ -558,6 +558,11 @@ SparsityPattern.
  (Wolfgang Bangerth, 2016/07/08)
  </li>
 
+ <li> New: A move constructor has been added to Triangulation.
+ <br>
+ (Daniel Shapero, 2016/07/07)
+ </li>
+
  <li> Fixed: The function DoFTools::dof_couplings_from_component_couplings
  for hp::FECollection arguments was compiled but not exported from the
  object file. This is now fixed.
index c0991e1eb77ce9d5a61401cf8673b4bab5687a90..64f910d8c216698050d1b882c3cdb3f57a9b1376 100644 (file)
@@ -1516,6 +1516,16 @@ public:
    */
   Triangulation (const Triangulation<dim, spacedim> &t);
 
+#ifdef DEAL_II_WITH_CXX11
+  /**
+   * Move constructor.
+   *
+   * Create a new triangulation by stealing the internal data of another
+   * triangulation.
+   */
+  Triangulation (Triangulation<dim, spacedim> &&tria);
+#endif
+
   /**
    * Delete the object and all levels of the hierarchy.
    */
index 1f3629c7e0afd5a4af3b15d376fca707519a3ab0..9f7a9c3b33f50dce877d3a3366f6870d810d18dd 100644 (file)
@@ -8996,13 +8996,56 @@ Triangulation (const Triangulation<dim, spacedim> &other)
 
 
 
+#ifdef DEAL_II_WITH_CXX11
+template <int dim, int spacedim>
+Triangulation<dim, spacedim>::
+Triangulation (Triangulation<dim, spacedim> &&tria)
+  :
+  Subscriptor(tria),
+  smooth_grid(tria.smooth_grid),
+  periodic_face_pairs_level_0(std::move(tria.periodic_face_pairs_level_0)),
+  periodic_face_map(std::move(tria.periodic_face_map)),
+  levels(std::move(tria.levels)),
+  faces(std::move(tria.faces)),
+  vertices(std::move(tria.vertices)),
+  vertices_used(std::move(tria.vertices_used)),
+  manifold(std::move(tria.manifold)),
+  anisotropic_refinement(tria.anisotropic_refinement),
+  check_for_distorted_cells(tria.check_for_distorted_cells),
+  number_cache(tria.number_cache),
+  vertex_to_boundary_id_map_1d(tria.vertex_to_boundary_id_map_1d),
+  vertex_to_manifold_id_map_1d(tria.vertex_to_manifold_id_map_1d)
+{
+  for (unsigned int i=0; i<tria.levels.size(); ++i)
+    tria.levels[i] = nullptr;
+
+  tria.faces = nullptr;
+
+  tria.number_cache = internal::Triangulation::NumberCache<dim>();
+
+  tria.vertex_to_boundary_id_map_1d = nullptr;
+  tria.vertex_to_manifold_id_map_1d = nullptr;
+}
+#endif
+
+
+
 template <int dim, int spacedim>
 Triangulation<dim, spacedim>::~Triangulation ()
 {
   for (unsigned int i=0; i<levels.size(); ++i)
-    delete levels[i];
+    if (levels[i])
+      {
+        delete levels[i];
+        levels[i] = 0;
+      }
   levels.clear ();
-  delete faces;
+
+  if (faces)
+    {
+      delete faces;
+      faces = 0;
+    }
 
   // the vertex_to_boundary_id_map_1d field
   // should be unused except in 1d
@@ -9010,14 +9053,23 @@ Triangulation<dim, spacedim>::~Triangulation ()
           ||
           (vertex_to_boundary_id_map_1d == 0),
           ExcInternalError());
-  delete vertex_to_boundary_id_map_1d;
+  if (vertex_to_boundary_id_map_1d)
+    {
+      delete vertex_to_boundary_id_map_1d;
+      vertex_to_boundary_id_map_1d = 0;
+    }
+
   // the vertex_to_manifold_id_map_1d field
   // should be unused except in 1d
   Assert ((dim == 1)
           ||
           (vertex_to_manifold_id_map_1d == 0),
           ExcInternalError());
-  delete vertex_to_manifold_id_map_1d;
+  if (vertex_to_manifold_id_map_1d)
+    {
+      delete vertex_to_manifold_id_map_1d;
+      vertex_to_manifold_id_map_1d = 0;
+    }
 }
 
 
diff --git a/tests/grid/move_constructor.cc b/tests/grid/move_constructor.cc
new file mode 100644 (file)
index 0000000..6ee0c49
--- /dev/null
@@ -0,0 +1,122 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2003 - 2015 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/manifold_lib.h>
+
+#include <fstream>
+
+
+template <int dim>
+void print_tria_info(const Triangulation<dim> &tria)
+{
+  deallog << (tria.n_active_cells() != 0) << ", "
+          << (tria.n_active_hexs() != 0) << ", "
+          << (tria.n_active_quads() != 0) << ", "
+          << (tria.n_active_lines() != 0) << ", "
+          << (tria.n_levels() != 0) << ", "
+          << (tria.n_vertices() != 0) << ", "
+          << (tria.get_periodic_face_map().size() != 0) << ", "
+          << (&tria.get_manifold(0)
+              == (const Manifold<dim> *)&Triangulation<dim>::straight_boundary)
+          << std::endl;
+}
+
+
+template <int dim>
+void test_hyper_cube()
+{
+  deallog << "Dimension: " << dim << std::endl;
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria, -1.0, 1.0);
+  tria.refine_global(2);
+
+  print_tria_info(tria);
+
+  Triangulation<dim> new_tria = std::move(tria);
+  print_tria_info(new_tria);
+  print_tria_info(tria);
+}
+
+
+template <int dim>
+void test_hyper_shell()
+{
+  deallog << "Dimension: " << dim << std::endl;
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_ball(tria);
+
+  const SphericalManifold<dim> boundary;
+  tria.set_all_manifold_ids_on_boundary(0);
+  tria.set_manifold(0, boundary);
+
+  tria.refine_global(3);
+  print_tria_info(tria);
+
+  Triangulation<dim> new_tria = std::move(tria);
+  print_tria_info(new_tria);
+  print_tria_info(tria);
+}
+
+
+template <int dim>
+void test_periodic_cube()
+{
+  deallog << "Dimension: " << dim << std::endl;
+
+  Triangulation<dim> tria;
+  GridGenerator::hyper_cube(tria, -1.0, 1.0, true);
+
+  using periodic_face_pair =
+    GridTools::PeriodicFacePair<typename Triangulation<dim>::cell_iterator>;
+  std::vector<periodic_face_pair> periodicity_vector;
+  GridTools::collect_periodic_faces(tria, 0, 1, 0, periodicity_vector);
+  tria.add_periodicity(periodicity_vector);
+  tria.refine_global(3);
+
+  print_tria_info(tria);
+
+  Triangulation<dim> new_tria = std::move(tria);
+  print_tria_info(new_tria);
+  print_tria_info(tria);
+}
+
+
+int main()
+{
+  initlog();
+
+  deallog << "Hyper-cube:" << std::endl;
+  test_hyper_cube<2>();
+  test_hyper_cube<3>();
+
+  deallog << std::endl << "Hyper-shell:" << std::endl;
+  test_hyper_shell<2>();
+  test_hyper_shell<3>();
+
+  deallog << std::endl << "Periodic hyper-cube:" << std::endl;
+  test_periodic_cube<2>();
+  test_periodic_cube<3>();
+
+  return 0;
+}
diff --git a/tests/grid/move_constructor.with_cxx11=on.output b/tests/grid/move_constructor.with_cxx11=on.output
new file mode 100644 (file)
index 0000000..76860e1
--- /dev/null
@@ -0,0 +1,30 @@
+
+DEAL::Hyper-cube:
+DEAL::Dimension: 2
+DEAL::1, 0, 1, 1, 1, 1, 0, 1
+DEAL::1, 0, 1, 1, 1, 1, 0, 1
+DEAL::0, 0, 0, 0, 0, 0, 0, 1
+DEAL::Dimension: 3
+DEAL::1, 1, 1, 1, 1, 1, 0, 1
+DEAL::1, 1, 1, 1, 1, 1, 0, 1
+DEAL::0, 0, 0, 0, 0, 0, 0, 1
+DEAL::
+DEAL::Hyper-shell:
+DEAL::Dimension: 2
+DEAL::1, 0, 1, 1, 1, 1, 0, 0
+DEAL::1, 0, 1, 1, 1, 1, 0, 0
+DEAL::0, 0, 0, 0, 0, 0, 0, 1
+DEAL::Dimension: 3
+DEAL::1, 1, 1, 1, 1, 1, 0, 0
+DEAL::1, 1, 1, 1, 1, 1, 0, 0
+DEAL::0, 0, 0, 0, 0, 0, 0, 1
+DEAL::
+DEAL::Periodic hyper-cube:
+DEAL::Dimension: 2
+DEAL::1, 0, 1, 1, 1, 1, 1, 1
+DEAL::1, 0, 1, 1, 1, 1, 1, 1
+DEAL::0, 0, 0, 0, 0, 0, 0, 1
+DEAL::Dimension: 3
+DEAL::1, 1, 1, 1, 1, 1, 1, 1
+DEAL::1, 1, 1, 1, 1, 1, 1, 1
+DEAL::0, 0, 0, 0, 0, 0, 0, 1

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.