]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Implement transformations of grids. Extend test case.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 26 Jul 2002 07:52:53 +0000 (07:52 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 26 Jul 2002 07:52:53 +0000 (07:52 +0000)
git-svn-id: https://svn.dealii.org/trunk@6290 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/grid_tools.h
deal.II/deal.II/source/grid/grid_tools.cc
deal.II/doc/news/2002/c-3-4.html
tests/deal.II/grid_tools.cc

index ec70cb5a5342b66983b9ef3486c8092a0a1fddfe..7992e0f0b041f69bce69f676a34fe14a570d10dc 100644 (file)
 
 
 #include <base/config.h>
+#include <grid/tria.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
 
-template <int dim> class Triangulation;
 
 
 
@@ -51,9 +53,120 @@ class GridTools
                                      */
     static
     double diameter (const Triangulation<1> &tria);
+
+                                    /**
+                                     * Transform the vertices of the
+                                     * given triangulation by
+                                     * applying the predicate to all
+                                     * its vertices. Since the
+                                     * internal consistency of a
+                                     * triangulation can only be
+                                     * guaranteed if the
+                                     * transformation is applied to
+                                     * the vertices of only one level
+                                     * of a hierarchically refined
+                                     * cells, this function may only
+                                     * be used on coarse grids,
+                                     * i.e. before any refinement of
+                                     * it has taken place.
+                                     *
+                                     * The predicate given as
+                                     * argument is used to transform
+                                     * each vertex. Its respective
+                                     * type has to offer a
+                                     * function-like syntax, i.e. the
+                                     * predicate is either an object
+                                     * of a type that has an
+                                     * @p{operator()}, or it is a
+                                     * pointer to the function. In
+                                     * either case, argument and
+                                     * return value have to be of
+                                     * type @p{Point<dim>}.
+                                     */
+    template <int dim, typename Predicate>
+    static
+    void transform (const Predicate    &predicate,
+                   Triangulation<dim> &triangulation);
+
+                                    /**
+                                     * Shift each vertex of the
+                                     * triangulation by the given
+                                     * shift vector. This function
+                                     * uses the @ref{transform}
+                                     * function above, so the
+                                     * requirements on the
+                                     * triangulation stated there
+                                     * hold for this function as
+                                     * well.
+                                     */
+    template <int dim>
+    static
+    void shift (const Point<dim>   &shift_vector,
+               Triangulation<dim> &triangulation);
+
+
+                                    /**
+                                     * Rotate all vertices of the
+                                     * given two-dimensional
+                                     * triangulation in
+                                     * counter-clockwise sense around
+                                     * the origin of the coordinate
+                                     * system by the given angle
+                                     * (given in radians, rather than
+                                     * degrees). This function uses
+                                     * the @ref{transform} function
+                                     * above, so the requirements on
+                                     * the triangulation stated there
+                                     * hold for this function as
+                                     * well.
+                                     */
+    static
+    void rotate (const double      angle,
+                Triangulation<2> &triangulation);
+    
+                                    /**
+                                     * Exception
+                                     */
+    DeclException0 (ExcTriangulationHasBeenRefined);
+};
+
+
+
+/* ----------------- Template function --------------- */
+
+template <int dim, typename Predicate>
+void GridTools::transform (const Predicate    &predicate,
+                          Triangulation<dim> &triangulation)
+{
+  Assert (triangulation.n_levels() == 1,
+         ExcTriangulationHasBeenRefined());
+  
+  std::vector<bool> treated_vertices (triangulation.n_vertices(),
+                                     false);
+
+                                  // loop over all active cells, and
+                                  // transform those vertices that
+                                  // have not yet been touched. note
+                                  // that we get to all vertices in
+                                  // the triangulation by only
+                                  // visiting the active cells.
+  typename Triangulation<dim>::active_cell_iterator
+    cell = triangulation.begin_active (),
+    endc = triangulation.end ();
+  for (; cell!=endc; ++cell)
+    for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_cell; ++v)
+      if (treated_vertices[cell->vertex_index(v)] == false)
+       {
+                                          // transform this vertex
+         cell->vertex(v) = predicate(cell->vertex(v));
+                                          // and mark it as treated
+         treated_vertices[cell->vertex_index(v)] = true;
+       };
 };
 
 
+
+
 /*----------------------------   grid_tools.h     ---------------------------*/
 /* end of #ifndef __deal2__grid_tools_H */
 #endif
index d8ee977a568aa95cd080c9ed849e6f47658203cc..605158cd738df85385008ce626b7006b78fd5474 100644 (file)
@@ -17,7 +17,9 @@
 #include <grid/tria.h>
 #include <grid/tria_accessor.h>
 #include <grid/tria_iterator.h>
+
 #include <cmath>
+#include <functional>
 
 #if deal_II_dimension != 1
 
@@ -93,8 +95,61 @@ GridTools::diameter (const Triangulation<1> &tria)
 
 
 
+// define some transformations in an anonymous namespace
+namespace 
+{
+  template <int dim>
+  inline
+  Point<dim> shift_point (const Point<dim> p,
+                         const Point<dim> shift)
+  {
+    return p+shift;
+  };
+
+
+
+  inline
+  Point<2> rotate2d (const Point<2> p,
+                    const double     angle)
+  {
+    return Point<2> (std::cos(angle)*p(0) - std::sin(angle) * p(1),
+                    std::sin(angle)*p(0) + std::cos(angle) * p(1));
+  };
+};
+
+
+template <int dim>
+void
+GridTools::shift (const Point<dim>   &shift_vector,
+                 Triangulation<dim> &triangulation)
+{
+  transform (std::bind2nd(ptr_fun(&shift_point<dim>),
+                         shift_vector),
+            triangulation);
+};
+
+
+#if deal_II_dimension == 2
+
+void
+GridTools::rotate (const double      angle,
+                  Triangulation<2> &triangulation)
+{
+  transform (std::bind2nd(ptr_fun(&rotate2d),
+                         angle),
+            triangulation);
+};
+
+#endif
+
+
+
 #if deal_II_dimension != 1
 template
 double
 GridTools::diameter<deal_II_dimension> (const Triangulation<deal_II_dimension> &);
 #endif
+
+template
+void GridTools::shift<deal_II_dimension> (const Point<deal_II_dimension> &,
+                                         Triangulation<deal_II_dimension> &);
index 3cfcd15f6ab12cf150f13e0ac91c874315cd865e..bc6e371e8e7d18d34afebc2cd751a7cc2d02cda1 100644 (file)
@@ -153,6 +153,14 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p> New: The <code class="class">GridTools</code> class now
+       offers functions to apply general transformations to all
+       vertices of a triangulation, and also more specialized
+       functions that shift and rotate triangulations.
+       <br>
+       (RH 2002/07/25)
+       </p>
+
   <li> <p> New: The existing <code
        class="member">FETools::extrapolate</code> functions does not
        take into account hanging node constraints. Therefore, it works
index 4b5b565265e8517fe7ed25f7e2f55cc57e1cf689..ed002b3700883239d00aebcb07198f99df092362 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 2001 by the deal.II authors
+//    Copyright (C) 2001, 2002 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -16,6 +16,7 @@
 #include <grid/tria.h>
 #include <grid/grid_generator.h>
 #include <grid/grid_tools.h>
+#include <grid/grid_out.h>
 
 #include <fstream>
 
@@ -24,8 +25,9 @@ std::ofstream logfile("grid_tools.output");
 
 
 
+// check GridTools::diameter
 template <int dim>
-void test ()
+void test1 ()
 {
                                   // test 1: hypercube
   if (true)
@@ -65,15 +67,38 @@ void test ()
 };
 
 
+// GridTools::transform
+void test2 ()
+{
+  Triangulation<2> tria;
+  GridGenerator::hyper_cube(tria);
+
+  logfile << "Unchanged grid:" << std::endl;
+  GridOut().write_gnuplot (tria, logfile);
+  
+  logfile << "Shifted grid:" << std::endl;
+  const Point<2> shift(1,2);
+  GridTools::shift (shift, tria);
+  GridOut().write_gnuplot (tria, logfile);
+
+  logfile << "Rotated grid:" << std::endl;
+  GridTools::rotate (3.14159265258/4, tria);
+  GridOut().write_gnuplot (tria, logfile);
+};
+
+
 int main ()
 {
   logfile.precision(4);
   deallog.attach(logfile);
   deallog.depth_console(0);
 
-  test<1> ();
-  test<2> ();
-  test<3> ();
+  test1<1> ();
+  test1<2> ();
+  test1<3> ();
+
+  test2 ();
   
   return 0;
 };
+

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.