(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.
*/
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.
*/
+#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
||
(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;
+ }
}
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// 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;
+}
--- /dev/null
+
+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