]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added copy boundary to manifold id function.
authorLuca Heltai <luca.heltai@sissa.it>
Fri, 9 Jan 2015 14:41:07 +0000 (15:41 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Sat, 10 Jan 2015 16:15:05 +0000 (17:15 +0100)
doc/news/changes.h
include/deal.II/grid/grid_tools.h
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/manifold/copy_boundary_to_manifold_id_01.cc [new file with mode: 0644]
tests/manifold/copy_boundary_to_manifold_id_01.output [new file with mode: 0644]

index a75bd50a396005c09013eccf0eea555929f451b5..52be64b6d50c1825d989d63b15eaaee9df1762f1 100644 (file)
@@ -134,6 +134,13 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+  <li> New: GridTools::copy_boundary_to_manifold_id() copies 
+  boundary_ids to manifold_ids for faces on the boundary.
+  <br>
+  (Luca Heltai, 2015/01/09)
+  </li>
+
+
   <li> Fixed: Utilities::int_to_string() produced wrong results if
   the number of digits specified was ten or greater.
   <br>
index 7f83820c19e2d9f5295c117603ee085e36ccc913..6208c6480ba70a4719fc83f611e5f5884fe7730e 100644 (file)
@@ -1145,6 +1145,40 @@ namespace GridTools
    const FullMatrix<double>                                          &matrix = FullMatrix<double>(),
    const std::vector<unsigned int>                                   &first_vector_components = std::vector<unsigned int>());
 
+  /*@}*/
+  /**
+   * @name Dealing with boundary and manifold ids
+   */
+  /*@{*/
+
+  /**
+   * Copy boundary ids to manifold ids. The default manifold_id for
+   * new Triangulation objects is numbers::invalid_manifold_id. When
+   * refinements occurs, the Triangulation asks where to locate new
+   * points to the underlying manifold, and if the manifold_id of a
+   * boundary face is set to numbers::invalid_manifold_id, then
+   * Triangulation reverts back boundary_id (this was the behavior of
+   * deal.II < 8.2). This function copies the boundary_ids of the
+   * boundary faces to the manifold_ids of the same faces, allowing
+   * the user to change the boundary_ids and use them for boundary
+   * conditions regardless of the geometry, which will use
+   * manifold_ids to create new points. Only active cells will be
+   * iterated over. This is a function you'd typically call when there
+   * is only one active level on your Triangulation.
+   *
+   * The optional parameter @p reset_boundary_ids, indicates wether
+   * this function should reset the boundary_ids of the Triangulation
+   * to its default value 0 after copying its value to the
+   * manifold_id. By default, boundary_ids are left untouched.
+   *
+   * @ingroup manifold
+   * @ingroup boundary
+   *
+   * @author Luca Heltai, 2015
+   */
+  template <int dim, int spacedim>
+  void copy_boundary_to_manifold_id(Triangulation<dim, spacedim> &tria,
+                                    bool reset_boundary_ids=false);
 
   /*@}*/
   /**
index 0e7b73969395b8152275ca6588351335454fc578..38b91b6504868dcffd7b853c3689f99575137a8b 100644 (file)
@@ -2894,6 +2894,25 @@ next_cell:
 #endif
   }
 
+  template <int dim, int spacedim>
+  void copy_boundary_to_manifold_id(Triangulation<dim, spacedim> &tria,
+                                    bool reset_boundary_ids)
+  {
+
+    typename Triangulation<dim,spacedim>::active_cell_iterator
+    cell=tria.begin_active(), endc=tria.end();
+
+    for (; cell != endc; ++cell)
+      for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+        if (cell->face(f)->at_boundary())
+          {
+            cell->face(f)->set_manifold_id
+            (static_cast<types::manifold_id>(cell->face(f)->boundary_indicator()));
+            if (reset_boundary_ids == true)
+              cell->face(f)->set_boundary_indicator(0);
+          }
+  }
+
 } /* namespace GridTools */
 
 
index d31b93d5eca25ed58964232ad4366aa400791616..cef43b58c7fdaf4c92f0ade4d2e537d6daed659e 100644 (file)
@@ -388,4 +388,13 @@ for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS
 #endif
 }
 
+for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS)
+{
+#if deal_II_space_dimension >= deal_II_dimension
+ namespace GridTools \{
+      template void copy_boundary_to_manifold_id<deal_II_dimension, deal_II_space_dimension>
+      (Triangulation<deal_II_dimension, deal_II_space_dimension> &, bool);
+\}
+#endif
+}
 
diff --git a/tests/manifold/copy_boundary_to_manifold_id_01.cc b/tests/manifold/copy_boundary_to_manifold_id_01.cc
new file mode 100644 (file)
index 0000000..e45e723
--- /dev/null
@@ -0,0 +1,74 @@
+//------------------------------------------------------------
+//    Copyright (C) 2015, the deal.II authors
+//
+//    This file is subject to LGPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//------------------------------------------------------------
+
+
+// Copy from boundary ids to manifold ids
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+
+
+// all include files you need here
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_out.h>
+
+template<int dim, int spacedim>
+void print_info(Triangulation<dim,spacedim> &tria) {
+  typename Triangulation<dim,spacedim>::active_cell_iterator cell;
+
+  for (cell = tria.begin_active(); cell != tria.end(); ++cell)
+    {
+      for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+       if(cell->face(f)->at_boundary()) 
+         deallog << "face: " << cell->face(f)
+                 << ", boundary_id: " 
+                 << (int)cell->face(f)->boundary_indicator()
+                 << ", manifold_id: " 
+                 << (int)cell->face(f)->manifold_id() << std::endl;
+    }
+}
+
+
+// Helper function
+template <int dim, int spacedim>
+void test()
+{
+  deallog << "Testing dim=" << dim
+          << ", spacedim="<< spacedim << std::endl;
+
+  Triangulation<dim,spacedim> tria;
+  GridGenerator::hyper_cube (tria, 0., 1., true);
+
+  print_info(tria);
+  GridTools::copy_boundary_to_manifold_id(tria);
+  print_info(tria);
+  GridTools::copy_boundary_to_manifold_id(tria, true);
+  print_info(tria);
+}
+
+int main ()
+{
+  initlog(true);
+  
+  test<1,1>();
+  test<1,2>();
+  test<2,2>();
+  test<2,3>();
+  test<3,3>();
+
+  return 0;
+}
+
diff --git a/tests/manifold/copy_boundary_to_manifold_id_01.output b/tests/manifold/copy_boundary_to_manifold_id_01.output
new file mode 100644 (file)
index 0000000..c3b3c64
--- /dev/null
@@ -0,0 +1,60 @@
+
+DEAL::Testing dim=1, spacedim=1
+DEAL::face: 0, boundary_id: 0, manifold_id: -1
+DEAL::face: 1, boundary_id: 1, manifold_id: -1
+DEAL::face: 0, boundary_id: 0, manifold_id: 0
+DEAL::face: 1, boundary_id: 1, manifold_id: 1
+DEAL::face: 0, boundary_id: 0, manifold_id: 0
+DEAL::face: 1, boundary_id: 0, manifold_id: 1
+DEAL::Testing dim=1, spacedim=2
+DEAL::face: 0, boundary_id: 0, manifold_id: -1
+DEAL::face: 1, boundary_id: 1, manifold_id: -1
+DEAL::face: 0, boundary_id: 0, manifold_id: 0
+DEAL::face: 1, boundary_id: 1, manifold_id: 1
+DEAL::face: 0, boundary_id: 0, manifold_id: 0
+DEAL::face: 1, boundary_id: 0, manifold_id: 1
+DEAL::Testing dim=2, spacedim=2
+DEAL::face: 1, boundary_id: 0, manifold_id: -1
+DEAL::face: 2, boundary_id: 1, manifold_id: -1
+DEAL::face: 0, boundary_id: 2, manifold_id: -1
+DEAL::face: 3, boundary_id: 3, manifold_id: -1
+DEAL::face: 1, boundary_id: 0, manifold_id: 0
+DEAL::face: 2, boundary_id: 1, manifold_id: 1
+DEAL::face: 0, boundary_id: 2, manifold_id: 2
+DEAL::face: 3, boundary_id: 3, manifold_id: 3
+DEAL::face: 1, boundary_id: 0, manifold_id: 0
+DEAL::face: 2, boundary_id: 0, manifold_id: 1
+DEAL::face: 0, boundary_id: 0, manifold_id: 2
+DEAL::face: 3, boundary_id: 0, manifold_id: 3
+DEAL::Testing dim=2, spacedim=3
+DEAL::face: 1, boundary_id: 0, manifold_id: -1
+DEAL::face: 2, boundary_id: 1, manifold_id: -1
+DEAL::face: 0, boundary_id: 2, manifold_id: -1
+DEAL::face: 3, boundary_id: 3, manifold_id: -1
+DEAL::face: 1, boundary_id: 0, manifold_id: 0
+DEAL::face: 2, boundary_id: 1, manifold_id: 1
+DEAL::face: 0, boundary_id: 2, manifold_id: 2
+DEAL::face: 3, boundary_id: 3, manifold_id: 3
+DEAL::face: 1, boundary_id: 0, manifold_id: 0
+DEAL::face: 2, boundary_id: 0, manifold_id: 1
+DEAL::face: 0, boundary_id: 0, manifold_id: 2
+DEAL::face: 3, boundary_id: 0, manifold_id: 3
+DEAL::Testing dim=3, spacedim=3
+DEAL::face: 2, boundary_id: 0, manifold_id: -1
+DEAL::face: 3, boundary_id: 1, manifold_id: -1
+DEAL::face: 0, boundary_id: 2, manifold_id: -1
+DEAL::face: 4, boundary_id: 3, manifold_id: -1
+DEAL::face: 1, boundary_id: 4, manifold_id: -1
+DEAL::face: 5, boundary_id: 5, manifold_id: -1
+DEAL::face: 2, boundary_id: 0, manifold_id: 0
+DEAL::face: 3, boundary_id: 1, manifold_id: 1
+DEAL::face: 0, boundary_id: 2, manifold_id: 2
+DEAL::face: 4, boundary_id: 3, manifold_id: 3
+DEAL::face: 1, boundary_id: 4, manifold_id: 4
+DEAL::face: 5, boundary_id: 5, manifold_id: 5
+DEAL::face: 2, boundary_id: 0, manifold_id: 0
+DEAL::face: 3, boundary_id: 0, manifold_id: 1
+DEAL::face: 0, boundary_id: 0, manifold_id: 2
+DEAL::face: 4, boundary_id: 0, manifold_id: 3
+DEAL::face: 1, boundary_id: 0, manifold_id: 4
+DEAL::face: 5, boundary_id: 0, manifold_id: 5

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.