From 4a4c30fc9b99443cde84db64264a90ceb79cf8b4 Mon Sep 17 00:00:00 2001 From: wolf Date: Fri, 26 Jul 2002 07:52:53 +0000 Subject: [PATCH] Implement transformations of grids. Extend test case. git-svn-id: https://svn.dealii.org/trunk@6290 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/grid/grid_tools.h | 115 +++++++++++++++++++++- deal.II/deal.II/source/grid/grid_tools.cc | 55 +++++++++++ deal.II/doc/news/2002/c-3-4.html | 8 ++ tests/deal.II/grid_tools.cc | 35 ++++++- 4 files changed, 207 insertions(+), 6 deletions(-) diff --git a/deal.II/deal.II/include/grid/grid_tools.h b/deal.II/deal.II/include/grid/grid_tools.h index ec70cb5a53..7992e0f0b0 100644 --- a/deal.II/deal.II/include/grid/grid_tools.h +++ b/deal.II/deal.II/include/grid/grid_tools.h @@ -15,8 +15,10 @@ #include +#include +#include +#include -template 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}. + */ + template + static + void transform (const Predicate &predicate, + Triangulation &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 + static + void shift (const Point &shift_vector, + Triangulation &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 +void GridTools::transform (const Predicate &predicate, + Triangulation &triangulation) +{ + Assert (triangulation.n_levels() == 1, + ExcTriangulationHasBeenRefined()); + + std::vector 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::active_cell_iterator + cell = triangulation.begin_active (), + endc = triangulation.end (); + for (; cell!=endc; ++cell) + for (unsigned int v=0; v::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 diff --git a/deal.II/deal.II/source/grid/grid_tools.cc b/deal.II/deal.II/source/grid/grid_tools.cc index d8ee977a56..605158cd73 100644 --- a/deal.II/deal.II/source/grid/grid_tools.cc +++ b/deal.II/deal.II/source/grid/grid_tools.cc @@ -17,7 +17,9 @@ #include #include #include + #include +#include #if deal_II_dimension != 1 @@ -93,8 +95,61 @@ GridTools::diameter (const Triangulation<1> &tria) +// define some transformations in an anonymous namespace +namespace +{ + template + inline + Point shift_point (const Point p, + const Point 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 +void +GridTools::shift (const Point &shift_vector, + Triangulation &triangulation) +{ + transform (std::bind2nd(ptr_fun(&shift_point), + 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 (const Triangulation &); #endif + +template +void GridTools::shift (const Point &, + Triangulation &); diff --git a/deal.II/doc/news/2002/c-3-4.html b/deal.II/doc/news/2002/c-3-4.html index 3cfcd15f6a..bc6e371e8e 100644 --- a/deal.II/doc/news/2002/c-3-4.html +++ b/deal.II/doc/news/2002/c-3-4.html @@ -153,6 +153,14 @@ contributor's names are abbreviated by WB (Wolfgang Bangerth), GK

deal.II

    +
  1. New: The GridTools class now + offers functions to apply general transformations to all + vertices of a triangulation, and also more specialized + functions that shift and rotate triangulations. +
    + (RH 2002/07/25) +

    +
  2. New: The existing FETools::extrapolate functions does not take into account hanging node constraints. Therefore, it works diff --git a/tests/deal.II/grid_tools.cc b/tests/deal.II/grid_tools.cc index 4b5b565265..ed002b3700 100644 --- a/tests/deal.II/grid_tools.cc +++ b/tests/deal.II/grid_tools.cc @@ -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 #include #include +#include #include @@ -24,8 +25,9 @@ std::ofstream logfile("grid_tools.output"); +// check GridTools::diameter template -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; }; + -- 2.39.5