]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
new mesh hyper_rectangle
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 16 Mar 2001 14:57:21 +0000 (14:57 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 16 Mar 2001 14:57:21 +0000 (14:57 +0000)
git-svn-id: https://svn.dealii.org/trunk@4237 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/grid_generator.h
deal.II/deal.II/source/grid/grid_generator.cc
deal.II/doc/news/2001/c-3-1.html

index dc548d9d1c40124be43b67791b593d86729afd6f..9efb8109733df7fe5ef904b6abedc03dbd17cb08 100644 (file)
@@ -33,6 +33,11 @@ template <int dim> class Triangulation;
  *       hypercube domain triangulated with exactly one element. You can
  *       get tensor product meshes by successive refinement of this cell.
  *
+ *    @item Rectangular coordinate parallel domains as a generalization
+ *      of hypercubes are generated by
+ *      @ref{GridGenerator}@p{::hyper_rectangle (tria, p1, p2)}, with two
+ *      opposite corner points @p{p1} and @p{p2}.
+ *
  *    @item Generalized L-shape domain:
  *      using the @ref{GridGenerator}@p{::hyper_L (tria, a,b)} function produces
  *      the hypercube with the interval $[a,b]$ without the hypercube
@@ -119,6 +124,21 @@ class GridGenerator
                            const double        left = 0.,
                            const double        right= 1.);
 
+                                    /**
+                                     * Create a coordinate-parallel
+                                     * parallelepiped from the two
+                                     * diagonally opposite corner
+                                     * points @p{p1} and @p{p2}.
+                                     *
+                                     * The triangulation needs to be
+                                     * void upon calling this
+                                     * function.
+                                     */
+    template <int dim>
+    static void hyper_rectangle (Triangulation<dim> &tria,
+                                const Point<dim>& p1,
+                                const Point<dim>& p2);
+
                                     /**
                                      * Hypercube with a layer of
                                      * hypercubes around it. The
index d3bb2715e3932160f0e9bfa4411b5023302d40d6..73b87fab774543540a1dcefb4da238d7c6d38c31 100644 (file)
 #include <cmath>
 
 
-#if deal_II_dimension == 1
-
-template <>
-void GridGenerator::hyper_cube<> (Triangulation<1> &tria,
-                                 const double left,
-                                 const double right) {
-  const Point<1> vertices[2] = { Point<1>(left), Point<1>(right) };
-  std::vector<CellData<1> > cells (1, CellData<1>());
-  cells[0].vertices[0] = 0;
-  cells[0].vertices[1] = 1;
-  cells[0].material_id = 0;
+template <int dim>
+void GridGenerator::hyper_rectangle (Triangulation<dim>& tria,
+                                    const Point<dim>& p_1,
+                                    const Point<dim>& p_2)
+{
+                                  // First, normalize input such that
+                                  // p1 is lower in all coordinate directions.
+  Point<dim> p1(p_1);
+  Point<dim> p2(p_2);
   
-  tria.create_triangulation (std::vector<Point<1> >(&vertices[0], &vertices[2]),
-                            cells,
-                            SubCellData());       // no boundary information
-};
+  for (unsigned int i=0;i<dim;++i)
+    if (p1(i) > p2(i))
+      {
+       double h = p1(i);
+       p1(i) = p2(i);
+       p2(i) = h;
+      }
+  
+  std::vector<Point<dim> > vertices (GeometryInfo<dim>::vertices_per_cell);
+  switch (dim)
+    {
+      case 1:
+       vertices[0](0) = p1(0);
+       vertices[1](0) = p2(0);
+       break;
+      case 2:
+       vertices[0](0) = p1(0);
+       vertices[1](0) = p2(0);
+       vertices[2](0) = p2(0);
+       vertices[3](0) = p1(0);
+       vertices[0](1) = p1(1);
+       vertices[1](1) = p1(1);
+       vertices[2](1) = p2(1);
+       vertices[3](1) = p2(1);
+       break;
+      case 3:
+       vertices[0](0) = p1(0);
+       vertices[1](0) = p2(0);
+       vertices[2](0) = p2(0);
+       vertices[3](0) = p1(0);
+       vertices[4](0) = p1(0);
+       vertices[5](0) = p2(0);
+       vertices[6](0) = p2(0);
+       vertices[7](0) = p1(0);
+       vertices[0](1) = p1(0);
+       vertices[1](1) = p1(0);
+       vertices[2](1) = p1(0);
+       vertices[3](1) = p1(0);
+       vertices[4](1) = p2(0);
+       vertices[5](1) = p2(0);
+       vertices[6](1) = p2(0);
+       vertices[7](1) = p2(0);
+       vertices[0](2) = p1(0);
+       vertices[1](2) = p1(0);
+       vertices[2](2) = p2(0);
+       vertices[3](2) = p2(0);
+       vertices[4](2) = p1(0);
+       vertices[5](2) = p1(0);
+       vertices[6](2) = p2(0);
+       vertices[7](2) = p2(0);
+       break;
+      default:
+       Assert (false, ExcNotImplemented ());
+    }
+
+                                  // Prepare cell data
+  std::vector<CellData<dim> > cells (1);
+  for (unsigned int i=0;i<GeometryInfo<dim>::vertices_per_cell;++i)
+    cells[0].vertices[i] = i;
+  cells[0].material_id = 0;
+
+  tria.create_triangulation (vertices, cells, SubCellData());
+}
+
+
+
+template <int dim>
+void GridGenerator::hyper_cube(Triangulation<dim> &tria,
+                              const double left,
+                              const double right)
+{
+  Point<dim> p1;
+  Point<dim> p2;
+  for (unsigned int i=0;i<dim;++i)
+    {
+      p1(i) = left;
+      p2(i) = right;
+    }
+  hyper_rectangle (tria, p1, p2);
+}
 
 
+#if deal_II_dimension == 1
 
 template <>
 void GridGenerator::hyper_cube_slit<> (Triangulation<1> &,
@@ -78,28 +153,6 @@ void GridGenerator::hyper_shell<> (Triangulation<1> &,
 
 #if deal_II_dimension == 2
 
-template <>
-void GridGenerator::hyper_cube<> (Triangulation<2> &tria,
-                                 const double left,
-                                 const double right)
-{
-  const Point<2> vertices[4] = { Point<2>(left,left),
-                                  Point<2>(right,left),
-                                  Point<2>(right,right),
-                                  Point<2>(left,right)  };
-  const int cell_vertices[1][4] = { { 0,1,2,3 } };
-  std::vector<CellData<2> > cells (1, CellData<2>());
-  for (unsigned int j=0; j<4; ++j)
-    cells[0].vertices[j] = cell_vertices[0][j];
-  cells[0].material_id = 0;
-  
-  tria.create_triangulation (std::vector<Point<2> >(&vertices[0], &vertices[4]),
-                            cells,
-                            SubCellData());       // no boundary information
-}
-
-
-
 template<>
 void GridGenerator::enclosed_hyper_cube (Triangulation<2> &tria,
                                         const double l,
@@ -383,32 +436,6 @@ void GridGenerator::half_hyper_shell<> (Triangulation<2>   &tria,
 
 #if deal_II_dimension == 3
 
-template <>
-void GridGenerator::hyper_cube<> (Triangulation<3> &tria,
-                                 const double left,
-                                 const double right) {
-  const Point<3> vertices[8] = { Point<3>(left,left,left),
-                                  Point<3>(right,left,left),
-                                  Point<3>(right,left,right),
-                                  Point<3>(left,left,right),
-
-                                  Point<3>(left,right,left),                            
-                                  Point<3>(right,right,left),
-                                  Point<3>(right,right,right),
-                                  Point<3>(left,right,right)};
-  const int cell_vertices[1][8] = { { 0,1,2,3,4,5,6,7 } };
-  std::vector<CellData<3> > cells (1, CellData<3>());
-  for (unsigned int j=0; j<8; ++j)
-    cells[0].vertices[j] = cell_vertices[0][j];
-  cells[0].material_id = 0;
-  
-  tria.create_triangulation (std::vector<Point<3> >(&vertices[0], &vertices[8]),
-                            cells,
-                            SubCellData());       // no boundary information
-};
-
-
-
 template <>
 void GridGenerator::hyper_cube_slit<> (Triangulation<3> &,
                                       const double,
@@ -630,9 +657,15 @@ void GridGenerator::hyper_shell<> (Triangulation<3>   &,
 
 
 // explicit instantiations
-// template void GridGenerator::hyper_cube (Triangulation<deal_II_dimension> &,
-//                                      const double,
-//                                      const double);
+template void
+GridGenerator::hyper_rectangle (Triangulation<deal_II_dimension> &,
+                               const Point<deal_II_dimension>&,
+                               const Point<deal_II_dimension>&);
+template void
+GridGenerator::hyper_cube (Triangulation<deal_II_dimension> &,
+                          const double,
+                          const double);
+
 // template void GridGenerator::hyper_L (Triangulation<deal_II_dimension> &,
 //                                   const double,
 //                                   const double);
index 882a396b1ea7fc315dc7b53e9c683ac640e58e3b..22b71f0effbd6505f225a509f787621bfb611ebc 100644 (file)
@@ -246,13 +246,23 @@ documentation, etc</a>.
 <h3>deal.II</h3>
 
 <ol>
+  <li> <p>
+       New: <code class="class">GridGenerator</code>::<code
+       class="member">hyper_rectangle</code> creates
+       coordinate-parallel rectangles in arbitrary dimension.
+       <br>
+       (GK 2001/03/16)
+       </p>
+
   <li> <p>
        Changed: The syntax of the <code
-       class="member">FiniteElemente<dim>::get_unit_support_points</code>
+       class="class">FiniteElement&lt;dim&gt;</code>::<code
+       class="member">get_unit_support_points</code>
        function is changed; it returns a reference to the vector of
        points. These points are now computed by the constructor of the
        <code class="class">FiniteElement</code> and not on each <code
-       class="member">FiniteElemente<dim>::get_unit_support_points</code>
+       class="class">FiniteElement&lt;dim&gt;</code>::<code
+       class="member">get_unit_support_points</code>
        function call as before.
        <br>
        (WB 2001/03/14)
@@ -260,7 +270,8 @@ documentation, etc</a>.
 
   <li> <p> 
        Removed: The function <code
-       class="member">FiniteElemente<dim>::get_support_points</code>
+       class="class">FiniteElement&lt;dim&gt;</code>::<code
+       class="member">get_support_points</code>
        is removed as the <code class="class">FiniteElement</code>
        cannot know the position of the support points by itself. This is
        because the support points depend on the unit support points

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.