From: Luca Heltai Date: Sun, 12 Apr 2015 16:14:18 +0000 (+0200) Subject: Made subdivided hyper cube/rectangle work with codimension one and two meshes. X-Git-Tag: v8.3.0-rc1~287^2~3 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=35adb35b502deb4a99c6f7d412359a1ec55ba80a;p=dealii.git Made subdivided hyper cube/rectangle work with codimension one and two meshes. --- diff --git a/doc/news/changes.h b/doc/news/changes.h index 2a5876482a..8f734286d3 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -428,6 +428,13 @@ inconvenience this causes.

Specific improvements

    +
  1. New: The GridGenerator::subdivided_hyper_cube() and + GridGenerator::subdivided_hyper_rectangle() now work also for codimension + one and two Triangulation(); +
    + (Luca Heltai, 2015/04/12) +
  2. +
  3. Changed: The cells of coarse meshes in parallel::distributed::Triangulation used to be ordered in a Cuthill-McKee numbering, which yields very high surface-to-volume ratios of the individual diff --git a/include/deal.II/grid/grid_generator.h b/include/deal.II/grid/grid_generator.h index c219ed3d0b..f5ec81f9ca 100644 --- a/include/deal.II/grid/grid_generator.h +++ b/include/deal.II/grid/grid_generator.h @@ -89,8 +89,8 @@ namespace GridGenerator * * @note The triangulation needs to be void upon calling this function. */ - template - void subdivided_hyper_cube (Triangulation &tria, + template + void subdivided_hyper_cube (Triangulation &tria, const unsigned int repetitions, const double left = 0., const double right= 1.); @@ -114,7 +114,7 @@ namespace GridGenerator void hyper_rectangle (Triangulation &tria, const Point &p1, const Point &p2, - const bool colorize = false); + const bool colorize = false); /** * Create a coordinate-parallel parallelepiped from the two diagonally @@ -147,13 +147,13 @@ namespace GridGenerator * @note For an example of the use of this function see the step-28 tutorial * program. */ - template + template void - subdivided_hyper_rectangle (Triangulation &tria, + subdivided_hyper_rectangle (Triangulation &tria, const std::vector &repetitions, - const Point &p1, - const Point &p2, - const bool colorize=false); + const Point &p1, + const Point &p2, + const bool colorize=false); /** * Like the previous function. However, here the second argument does not @@ -176,7 +176,7 @@ namespace GridGenerator const std::vector > &step_sizes, const Point &p_1, const Point &p_2, - const bool colorize); + const bool colorize); /** * Like the previous function, but with the following twist: the @p diff --git a/source/grid/grid_generator.cc b/source/grid/grid_generator.cc index 093feed06a..b4cb0578ad 100644 --- a/source/grid/grid_generator.cc +++ b/source/grid/grid_generator.cc @@ -99,11 +99,11 @@ namespace GridGenerator - template + template void - colorize_subdivided_hyper_rectangle (Triangulation &tria, - const Point &p1, - const Point &p2, + colorize_subdivided_hyper_rectangle (Triangulation &tria, + const Point &p1, + const Point &p2, const double epsilon) { @@ -114,13 +114,13 @@ namespace GridGenerator // should be smaller than the smallest cell // diameter. - typename Triangulation::face_iterator face = tria.begin_face(), - endface = tria.end_face(); + typename Triangulation::face_iterator face = tria.begin_face(), + endface = tria.end_face(); for (; face!=endface; ++face) { if (face->boundary_indicator() == 0) { - const Point center (face->center()); + const Point center (face->center()); if (std::abs(center(0)-p1[0]) < epsilon) face->set_boundary_indicator(0); else if (std::abs(center(0) - p2[0]) < epsilon) @@ -142,7 +142,7 @@ namespace GridGenerator } } - for (typename Triangulation::cell_iterator cell = tria.begin(); + for (typename Triangulation::cell_iterator cell = tria.begin(); cell != tria.end(); ++cell) { char id = 0; @@ -879,9 +879,9 @@ namespace GridGenerator } - template + template void - subdivided_hyper_cube (Triangulation &tria, + subdivided_hyper_cube (Triangulation &tria, const unsigned int repetitions, const double left, const double right) @@ -890,115 +890,26 @@ namespace GridGenerator Assert (left < right, ExcMessage ("Invalid left-to-right bounds of hypercube")); - // first generate the necessary - // points - const double delta = (right-left)/repetitions; - std::vector > points; - switch (dim) + Point p0,p1; + for (unsigned int i=0; i (left+x*delta)); - break; - - case 2: - for (unsigned int y=0; y<=repetitions; ++y) - for (unsigned int x=0; x<=repetitions; ++x) - points.push_back (Point (left+x*delta, - left+y*delta)); - break; - - case 3: - for (unsigned int z=0; z<=repetitions; ++z) - for (unsigned int y=0; y<=repetitions; ++y) - for (unsigned int x=0; x<=repetitions; ++x) - points.push_back (Point (left+x*delta, - left+y*delta, - left+z*delta)); - break; - - default: - Assert (false, ExcNotImplemented()); + p0[i] = left; + p1[i] = right; } - // next create the cells - // Prepare cell data - std::vector > cells; - // Define these as abbreviations - // for the step sizes below. The - // number of points in a single - // direction is repetitions+1 - const unsigned int dy = repetitions+1; - const unsigned int dz = dy*dy; - switch (dim) - { - case 1: - cells.resize (repetitions); - for (unsigned int x=0; x reps(dim, repetitions); + subdivided_hyper_rectangle(tria, reps, p0, p1); } - template + template void subdivided_hyper_rectangle ( - Triangulation &tria, + Triangulation &tria, const std::vector &repetitions, - const Point &p_1, - const Point &p_2, + const Point &p_1, + const Point &p_2, const bool colorize) { // contributed by Joerg R. Weimar @@ -1008,8 +919,8 @@ namespace GridGenerator // First, normalize input such that // p1 is lower in all coordinate // directions. - Point p1(p_1); - Point p2(p_2); + Point p1(p_1); + Point p2(p_2); for (unsigned int i=0; i p2(i)) @@ -1019,39 +930,38 @@ namespace GridGenerator // are >= 1, and calculate deltas // convert repetitions from double // to int by taking the ceiling. - Point delta; + std::vector > delta(dim); for (unsigned int i=0; i= 1, ExcInvalidRepetitions(repetitions[i])); - delta[i] = (p2[i]-p1[i])/repetitions[i]; + delta[i][i] = (p2[i]-p1[i])/repetitions[i]; } // then generate the necessary // points - std::vector > points; + std::vector > points; switch (dim) { case 1: for (unsigned int x=0; x<=repetitions[0]; ++x) - points.push_back (Point (p1[0]+x*delta[0])); + points.push_back (p1+(double)x*delta[0]); break; case 2: for (unsigned int y=0; y<=repetitions[1]; ++y) for (unsigned int x=0; x<=repetitions[0]; ++x) - points.push_back (Point (p1[0]+x*delta[0], - p1[1]+y*delta[1])); + points.push_back (p1+(double)x*delta[0] + +(double)y*delta[1]); break; case 3: for (unsigned int z=0; z<=repetitions[2]; ++z) for (unsigned int y=0; y<=repetitions[1]; ++y) for (unsigned int x=0; x<=repetitions[0]; ++x) - points.push_back (Point (p1[0]+x*delta[0], - p1[1]+y*delta[1], - p1[2]+z*delta[2])); + points.push_back (p1+(double)x*delta[0] + + (double)y*delta[1] + (double)z*delta[2]); break; default: @@ -1133,8 +1043,9 @@ namespace GridGenerator // use a large epsilon to // compare numbers to avoid // roundoff problems. - const double epsilon - = 0.01 * *std::min_element (&delta[0], &delta[0]+dim); + double epsilon = 10; + for (unsigned int i=0; i 0, ExcMessage ("The distance between corner points must be positive.")) diff --git a/source/grid/grid_generator.inst.in b/source/grid/grid_generator.inst.in index f17220b8d1..f7f6ef0ae8 100644 --- a/source/grid/grid_generator.inst.in +++ b/source/grid/grid_generator.inst.in @@ -29,6 +29,20 @@ namespace GridGenerator hyper_cube ( Triangulation &, const double, const double, const bool); + template void + subdivided_hyper_cube ( + Triangulation &, + const unsigned int, const double, const double); + + + template void + subdivided_hyper_rectangle + (Triangulation &, + const std::vector&, + const Point&, + const Point&, + bool); + template void merge_triangulations @@ -80,20 +94,6 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS for (deal_II_dimension : DIMENSIONS) { namespace GridGenerator \{ - template void - subdivided_hyper_cube ( - Triangulation &, - const unsigned int, const double, const double); - - - - template void - subdivided_hyper_rectangle - (Triangulation &, - const std::vector&, - const Point&, - const Point&, - bool); template void