]> https://gitweb.dealii.org/ - dealii.git/commitdiff
GridGenerator::flatten_triangulation. 95/head
authorLuca Heltai <luca.heltai@sissa.it>
Tue, 19 Aug 2014 12:37:08 +0000 (14:37 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Wed, 20 Aug 2014 14:03:02 +0000 (16:03 +0200)
Added GridTools::flatten_triangulation + one test.

Fixed WB comments.

Moved flatten triangulation to GridGenerator.

Fixed TH and WB comments.

Fixed indentation issue.

doc/news/changes.h
include/deal.II/grid/grid_generator.h
include/deal.II/grid/grid_tools.h
source/grid/grid_generator.cc
source/grid/grid_generator.inst.in
source/grid/grid_out.inst.in
source/grid/grid_tools.cc
source/grid/grid_tools.inst.in
tests/grid/flatten_grid_01.cc [new file with mode: 0644]
tests/grid/flatten_grid_01.output [new file with mode: 0644]

index 2e6d6f72b8c54380229fdaa01d3411bbab32864f..859da6936a4bbf3c69131b38f937e24c02b7ea2f 100644 (file)
@@ -97,6 +97,24 @@ inconvenience this causes.
 
 
 <ol>
+  <li> New: There is now a GridGenerator::flatten_triangulation()
+  taking a Triangulation<dim, spacedim_1> as input and returning
+  a Triangulation<dim, spacedim_2> as output. The output
+  triangulation will contain a single level with all active
+  cells of the input triangulation, and will be topologically
+  equivalent to the input triangulation. If the two spacedimensions
+  are equal, then this function will copy the triangulation
+  removing all levels, e.g., flattening it. If the two spacedimensions
+  are different, then this function will copy the vertices only
+  up to the smallest spacedimension parameter. <br>
+  Using this function, you can create a Triangulation<2,3> from
+  a Triangulation<2,2> or project to the plane z=0 your
+  Triangulation<2,3>. No checks are performed on the validity of
+  the resulting Triangulation.
+  <br>
+  (Luca Heltai, 2014/08/19)
+  </li>
+
   <li> Fixed: Newer versions of GCC (e.g. 4.9.x) are no longer compatible
   with BOOST 1.46. Consequently, the CMake scripts now require at least
   BOOST 1.48 in order to use a BOOST version found on the system. If no
index 44a45c0edf664bd14d2bc0fe03c5be52d2a5d404..14ea1a5924e3bcfdf4ea1e930d6095e1d519f311 100644 (file)
@@ -690,6 +690,39 @@ namespace GridGenerator
                          const double               height,
                          Triangulation<3,3>        &result);
 
+  /**
+   * Given an input triangulation @p in_tria, this function makes a new
+   * flat triangulation @p out_tria which contains a single level with
+   * all active cells of the input triangulation. If spacedim1 and
+   * spacedim2 are different, only the smallest spacedim components of
+   * the vertices are copied over. This is useful to create a
+   * Triangulation<2,3> out of a Triangulation<2,2>, or to project a
+   * Triangulation<2,3> into a Triangulation<2,2>, by neglecting the z
+   * components of the vertices.
+   *
+   * No internal checks are performed on the vertices, which are
+   * assumed to make sense topologically in the target #spacedim2
+   * dimensional space. If this is not the case, you will encounter
+   * problems when using the triangulation later on.
+   *
+   * All informations about cell manifold_ids and material ids are
+   * copied from one triangulation to the other, and only the boundary
+   * manifold_ids and boundary_ids are copied over from the faces of
+   * @p in_tria to the faces of @p out_tria. If you need to specify
+   * manifold ids on interior faces, they have to be specified
+   * manually after the triangulation is created.
+   *
+   * This function will fail if the input Triangulation is of type
+   * parallel::distributed::Triangulation, as well as when the input
+   * Triangulation contains hanging nodes.
+   *
+   * @author Luca Heltai, 2014
+   */
+  template <int dim, int spacedim1, int spacedim2>
+  void flatten_triangulation(const Triangulation<dim,spacedim1> &in_tria,
+                            Triangulation<dim,spacedim2> &out_tria);
+
+
   /**
    * This function transforms the @p Triangulation @p tria smoothly to a
    * domain that is described by the boundary points in the map @p
index 76dbe802e342728d4d68c300331e5741af5d6814..4355d8984a6d19b290fb1188041650945388e3b3 100644 (file)
@@ -1001,6 +1001,7 @@ namespace GridTools
   std::vector<typename Container::active_cell_iterator>
   get_patch_around_cell(const typename Container::active_cell_iterator &cell);
 
+
   /*@}*/
   /**
    *  @name Lower-dimensional meshes for parts of higher-dimensional meshes 
index 4f4cf102bfe71742f340cebfc2210cba41a62428..3eca78a3f7b3db5f9550cbbcbd94a6b2d0f61d3d 100644 (file)
@@ -38,6 +38,8 @@
 #include <deal.II/fe/fe_q.h>
 #include <deal.II/numerics/matrix_tools.h>
 
+#include <deal.II/distributed/tria.h>
+
 #include <iostream>
 #include <cmath>
 #include <limits>
@@ -3866,6 +3868,87 @@ namespace GridGenerator
             }
       }
   }
+  
+  template <int dim, int spacedim1, int spacedim2>
+  void flatten_triangulation(const Triangulation<dim, spacedim1> &in_tria,
+                            Triangulation<dim,spacedim2> &out_tria) 
+  {
+    const parallel::distributed::Triangulation<dim, spacedim1> * pt =
+      dynamic_cast<const parallel::distributed::Triangulation<dim, spacedim1> *>(&in_tria);
+    
+    Assert (pt == NULL, 
+           ExcMessage("Cannot use this function on parallel::distributed::Triangulation."));
+
+    std::vector<Point<spacedim2> > v;
+    std::vector<CellData<dim> > cells;
+    SubCellData subcelldata;
+  
+    const unsigned int spacedim = std::min(spacedim1,spacedim2);
+    const std::vector<Point<spacedim1> > &in_vertices = in_tria.get_vertices();
+  
+    v.resize(in_vertices.size());
+    for(unsigned int i=0; i<in_vertices.size(); ++i)
+      for(unsigned int d=0; d<spacedim; ++d)
+       v[i][d] = in_vertices[i][d];
+
+    cells.resize(in_tria.n_active_cells());
+    typename Triangulation<dim,spacedim1>::active_cell_iterator 
+      cell = in_tria.begin_active(),
+      endc = in_tria.end();
+  
+    for(unsigned int id=0; cell != endc; ++cell, ++id) 
+      {
+       for(unsigned int i=0; i<GeometryInfo<dim>::vertices_per_cell; ++i) 
+         cells[id].vertices[i] = cell->vertex_index(i);
+       cells[id].material_id = cell->material_id();
+       cells[id].manifold_id = cell->manifold_id();
+      }
+
+    if(dim>1) {
+      typename Triangulation<dim,spacedim1>::active_face_iterator
+       face = in_tria.begin_active_face(),
+       endf = in_tria.end_face();
+      
+      // Face counter for both dim == 2 and dim == 3
+      unsigned int f=0;
+      switch(dim) {
+      case 2:
+       {
+         subcelldata.boundary_lines.resize(in_tria.n_active_faces());
+         for(; face != endf; ++face)
+           if(face->at_boundary()) 
+             {
+               for(unsigned int i=0; i<GeometryInfo<dim>::vertices_per_face; ++i) 
+                 subcelldata.boundary_lines[f].vertices[i] = face->vertex_index(i);
+               subcelldata.boundary_lines[f].boundary_id = face->boundary_indicator();
+               subcelldata.boundary_lines[f].manifold_id = face->manifold_id();
+               ++f;
+             }
+         subcelldata.boundary_lines.resize(f);
+       }
+       break;
+      case 3:
+       {
+         subcelldata.boundary_quads.resize(in_tria.n_active_faces());
+         for(; face != endf; ++face) 
+           if(face->at_boundary()) 
+             {
+               for(unsigned int i=0; i<GeometryInfo<dim>::vertices_per_face; ++i) 
+                 subcelldata.boundary_quads[f].vertices[i] = face->vertex_index(i);
+               subcelldata.boundary_quads[f].boundary_id = face->boundary_indicator();
+               subcelldata.boundary_quads[f].manifold_id = face->manifold_id();
+               ++f;
+             }
+         subcelldata.boundary_quads.resize(f);
+       }
+       break;
+      default:
+       Assert(false, ExcInternalError());
+      }
+    }
+    out_tria.create_triangulation(v, cells, subcelldata);
+  }
+
 }
 
 // explicit instantiations
index 4101f47e36e4137f59d749ad4fcda33fb3d91bea..74600bcddad495f9e4caf1fc8f3ad0e5420da39c 100644 (file)
@@ -99,3 +99,16 @@ namespace GridGenerator \{
 \}  
  }
 
+
+for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS; deal_II_space_dimension_2 : SPACE_DIMENSIONS)
+{
+namespace GridGenerator \{
+#if (deal_II_dimension <= deal_II_space_dimension) && (deal_II_dimension <= deal_II_space_dimension_2)
+     template 
+     void 
+     flatten_triangulation<>(const Triangulation<deal_II_dimension,deal_II_space_dimension> &,
+                             Triangulation<deal_II_dimension,deal_II_space_dimension_2>&);
+#endif
+\}    
+}
+
index c9cbe7f4bbe153f639cfe5f32742c76edc76a38a..bb3390e29efb6d47af093c46a8477981b824e05d 100644 (file)
@@ -96,6 +96,9 @@ for (deal_II_dimension : DIMENSIONS)
    template void GridOut::write_ucd
      (const Triangulation<1, 3> &,
       std::ostream &) const;
+   template void GridOut::write_msh
+     (const Triangulation<1, 3> &,
+      std::ostream &) const;
 #endif
   }
 
index c2ba9e2d4a4e3b54681219d54374e8546cb38043..01d3e4e9fcdcd9962b2233e8a8c44bd1ea40cad8 100644 (file)
@@ -2595,8 +2595,6 @@ next_cell:
 #endif
   }
 
-
-
 } /* namespace GridTools */
 
 
index 6191916d47eec55b7e9048be9d5ce6cdc46d200f..6198f016cb884bd86f02354ce6a2efb261caaa99 100644 (file)
@@ -373,3 +373,4 @@ for (deal_II_dimension : DIMENSIONS ; deal_II_space_dimension : SPACE_DIMENSIONS
 #endif
 }
 
+
diff --git a/tests/grid/flatten_grid_01.cc b/tests/grid/flatten_grid_01.cc
new file mode 100644 (file)
index 0000000..c795eec
--- /dev/null
@@ -0,0 +1,63 @@
+// ---------------------------------------------------------------------
+// $Id$
+//
+// Copyright (C) 2004 - 2013 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.
+//
+// ---------------------------------------------------------------------
+
+// Generate a grid, refine it once, flatten it and output the result.
+
+#include "../tests.h"
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_generator.h>
+
+template <int dim, int spacedim1, int spacedim2> 
+void test() {
+  deallog << "Testing <" << dim << "," << spacedim1
+         << "> VS <" << dim << "," << spacedim2 
+         << ">" << std::endl;
+  
+  Triangulation<dim, spacedim1> tria1;
+  GridGenerator::hyper_cube(tria1);
+  tria1.refine_global(1);
+  
+  Triangulation<dim, spacedim2> tria2;
+  GridGenerator::flatten_triangulation(tria1, tria2);
+  GridOut go;
+  go.write_msh(tria2, deallog.get_file_stream());
+}
+
+int main()
+{
+  initlog();
+  test<1,1,1>();
+  test<1,1,2>();
+  test<1,1,3>();
+  // 
+  test<1,2,1>();
+  test<1,2,2>();
+  test<1,2,3>();
+  // 
+  test<1,3,1>();
+  test<1,3,2>();
+  test<1,3,3>();
+  //
+  test<2,2,2>();
+  test<2,2,3>();
+  //
+  test<2,3,2>();
+  test<2,3,3>();
+  //
+  test<3,3,3>();
+}
diff --git a/tests/grid/flatten_grid_01.output b/tests/grid/flatten_grid_01.output
new file mode 100644 (file)
index 0000000..7304f8a
--- /dev/null
@@ -0,0 +1,231 @@
+
+DEAL::Testing <1,1> VS <1,1>
+$NOD
+3
+1  0.00000 0 0
+2  1.00000 0 0
+3  0.500000 0 0
+$ENDNOD
+$ELM
+2
+1 1 0 0 2 1 3 
+2 1 0 0 2 3 2 
+$ENDELM
+DEAL::Testing <1,1> VS <1,2>
+$NOD
+3
+1  0.00000 0.00000 0
+2  1.00000 0.00000 0
+3  0.500000 0.00000 0
+$ENDNOD
+$ELM
+2
+1 1 0 0 2 1 3 
+2 1 0 0 2 3 2 
+$ENDELM
+DEAL::Testing <1,1> VS <1,3>
+$NOD
+3
+1  0.00000 0.00000 0.00000
+2  1.00000 0.00000 0.00000
+3  0.500000 0.00000 0.00000
+$ENDNOD
+$ELM
+2
+1 1 0 0 2 1 3 
+2 1 0 0 2 3 2 
+$ENDELM
+DEAL::Testing <1,2> VS <1,1>
+$NOD
+3
+1  0.00000 0 0
+2  1.00000 0 0
+3  0.500000 0 0
+$ENDNOD
+$ELM
+2
+1 1 0 0 2 1 3 
+2 1 0 0 2 3 2 
+$ENDELM
+DEAL::Testing <1,2> VS <1,2>
+$NOD
+3
+1  0.00000 0.00000 0
+2  1.00000 0.00000 0
+3  0.500000 0.00000 0
+$ENDNOD
+$ELM
+2
+1 1 0 0 2 1 3 
+2 1 0 0 2 3 2 
+$ENDELM
+DEAL::Testing <1,2> VS <1,3>
+$NOD
+3
+1  0.00000 0.00000 0.00000
+2  1.00000 0.00000 0.00000
+3  0.500000 0.00000 0.00000
+$ENDNOD
+$ELM
+2
+1 1 0 0 2 1 3 
+2 1 0 0 2 3 2 
+$ENDELM
+DEAL::Testing <1,3> VS <1,1>
+$NOD
+3
+1  0.00000 0 0
+2  1.00000 0 0
+3  0.500000 0 0
+$ENDNOD
+$ELM
+2
+1 1 0 0 2 1 3 
+2 1 0 0 2 3 2 
+$ENDELM
+DEAL::Testing <1,3> VS <1,2>
+$NOD
+3
+1  0.00000 0.00000 0
+2  1.00000 0.00000 0
+3  0.500000 0.00000 0
+$ENDNOD
+$ELM
+2
+1 1 0 0 2 1 3 
+2 1 0 0 2 3 2 
+$ENDELM
+DEAL::Testing <1,3> VS <1,3>
+$NOD
+3
+1  0.00000 0.00000 0.00000
+2  1.00000 0.00000 0.00000
+3  0.500000 0.00000 0.00000
+$ENDNOD
+$ELM
+2
+1 1 0 0 2 1 3 
+2 1 0 0 2 3 2 
+$ENDELM
+DEAL::Testing <2,2> VS <2,2>
+$NOD
+9
+1  0.00000 0.00000 0
+2  1.00000 0.00000 0
+3  0.00000 1.00000 0
+4  1.00000 1.00000 0
+5  0.500000 0.00000 0
+6  0.00000 0.500000 0
+7  1.00000 0.500000 0
+8  0.500000 1.00000 0
+9  0.500000 0.500000 0
+$ENDNOD
+$ELM
+4
+1 3 0 0 4 1 5 9 6 
+2 3 0 0 4 5 2 7 9 
+3 3 0 0 4 6 9 8 3 
+4 3 0 0 4 9 7 4 8 
+$ENDELM
+DEAL::Testing <2,2> VS <2,3>
+$NOD
+9
+1  0.00000 0.00000 0.00000
+2  1.00000 0.00000 0.00000
+3  0.00000 1.00000 0.00000
+4  1.00000 1.00000 0.00000
+5  0.500000 0.00000 0.00000
+6  0.00000 0.500000 0.00000
+7  1.00000 0.500000 0.00000
+8  0.500000 1.00000 0.00000
+9  0.500000 0.500000 0.00000
+$ENDNOD
+$ELM
+4
+1 3 0 0 4 1 5 9 6 
+2 3 0 0 4 5 2 7 9 
+3 3 0 0 4 6 9 8 3 
+4 3 0 0 4 9 7 4 8 
+$ENDELM
+DEAL::Testing <2,3> VS <2,2>
+$NOD
+9
+1  0.00000 0.00000 0
+2  1.00000 0.00000 0
+3  0.00000 1.00000 0
+4  1.00000 1.00000 0
+5  0.500000 0.00000 0
+6  0.00000 0.500000 0
+7  1.00000 0.500000 0
+8  0.500000 1.00000 0
+9  0.500000 0.500000 0
+$ENDNOD
+$ELM
+4
+1 3 0 0 4 1 5 9 6 
+2 3 0 0 4 5 2 7 9 
+3 3 0 0 4 6 9 8 3 
+4 3 0 0 4 9 7 4 8 
+$ENDELM
+DEAL::Testing <2,3> VS <2,3>
+$NOD
+9
+1  0.00000 0.00000 0.00000
+2  1.00000 0.00000 0.00000
+3  0.00000 1.00000 0.00000
+4  1.00000 1.00000 0.00000
+5  0.500000 0.00000 0.00000
+6  0.00000 0.500000 0.00000
+7  1.00000 0.500000 0.00000
+8  0.500000 1.00000 0.00000
+9  0.500000 0.500000 0.00000
+$ENDNOD
+$ELM
+4
+1 3 0 0 4 1 5 9 6 
+2 3 0 0 4 5 2 7 9 
+3 3 0 0 4 6 9 8 3 
+4 3 0 0 4 9 7 4 8 
+$ENDELM
+DEAL::Testing <3,3> VS <3,3>
+$NOD
+27
+1  0.00000 0.00000 0.00000
+2  1.00000 0.00000 0.00000
+3  0.00000 1.00000 0.00000
+4  1.00000 1.00000 0.00000
+5  0.00000 0.00000 1.00000
+6  1.00000 0.00000 1.00000
+7  0.00000 1.00000 1.00000
+8  1.00000 1.00000 1.00000
+9  0.500000 0.00000 0.00000
+10  0.00000 0.500000 0.00000
+11  0.00000 0.00000 0.500000
+12  1.00000 0.500000 0.00000
+13  1.00000 0.00000 0.500000
+14  0.500000 1.00000 0.00000
+15  0.00000 1.00000 0.500000
+16  1.00000 1.00000 0.500000
+17  0.500000 0.00000 1.00000
+18  0.00000 0.500000 1.00000
+19  1.00000 0.500000 1.00000
+20  0.500000 1.00000 1.00000
+21  0.500000 0.00000 0.500000
+22  0.500000 0.500000 0.00000
+23  0.00000 0.500000 0.500000
+24  1.00000 0.500000 0.500000
+25  0.500000 1.00000 0.500000
+26  0.500000 0.500000 1.00000
+27  0.500000 0.500000 0.500000
+$ENDNOD
+$ELM
+8
+1 5 0 0 8 1 9 21 11 10 22 27 23 
+2 5 0 0 8 9 2 13 21 22 12 24 27 
+3 5 0 0 8 10 22 27 23 3 14 25 15 
+4 5 0 0 8 22 12 24 27 14 4 16 25 
+5 5 0 0 8 11 21 17 5 23 27 26 18 
+6 5 0 0 8 21 13 6 17 27 24 19 26 
+7 5 0 0 8 23 27 26 18 15 25 20 7 
+8 5 0 0 8 27 24 19 26 25 16 8 20 
+$ENDELM

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.