From 60791f33919c449675921f90a7fda992a3ddad5d Mon Sep 17 00:00:00 2001
From: brian <brian@0785d39b-7218-0410-832d-ea1e28bc413d>
Date: Mon, 16 Dec 2002 18:01:14 +0000
Subject: [PATCH] added half_hyper_ball class

git-svn-id: https://svn.dealii.org/trunk@6826 0785d39b-7218-0410-832d-ea1e28bc413d
---
 deal.II/deal.II/include/grid/grid_generator.h | 33 +++++++++++++++-
 deal.II/deal.II/source/grid/grid_generator.cc | 39 +++++++++++++++++++
 2 files changed, 71 insertions(+), 1 deletion(-)

diff --git a/deal.II/deal.II/include/grid/grid_generator.h b/deal.II/deal.II/include/grid/grid_generator.h
index d27009c435..b196539aad 100644
--- a/deal.II/deal.II/include/grid/grid_generator.h
+++ b/deal.II/deal.II/include/grid/grid_generator.h
@@ -64,7 +64,21 @@ template <typename number> class SparseMatrix;
  *      to the triangulation object you passed to this function if you later want
  *      the triangulation to be refined at the outer boundaries.
  *
- *    @item Hyper shell: A hyper shell is the region between two hyper
+ *    @item Half Hyper ball:
+ *      You get half of the circle generated by Hyper ball.
+ *      with center @p{p} and with radius @p{r} by calling
+ *      @ref{GridGenerator}@p{::half_hyper_ball (tria, p, r)}. The half-circle is 
+ *      triangulated by four cells. The diameter of the center cell is
+ *      chosen to be the same as for the Hyper ball class. 
+ *      To create a half-hyperball in one dimension results in
+ *      an error.
+ *
+ *      Do not forget to also attach a suitable boundary approximation object
+ *      to the triangulation object you passed to this function if you later want
+ *      the triangulation to be refined at the outer boundaries. The class
+ *      @ref{HalfHyperBallBoundary} will provide a boundary object.
+ *
+*    @item Hyper shell: A hyper shell is the region between two hyper
  *      sphere with the same origin. Therefore, it is a ring in two
  *      spatial dimensions. To triangulation it, call the function
  *      @ref{GridGenerator}@p{::hyper_shell (tria, origin, inner_radius, outer_radius, N)},
@@ -448,6 +462,23 @@ class GridGenerator
 			     const unsigned int  n_cells = 0);
 
 
+				     /**
+				      * This class produces a half hyper-ball,
+				      * which contains four elements.
+				      *
+				      * The triangulation needs to be
+				      * void upon calling this
+				      * function.
+				      *
+				      * Currently, only a two-dimensional 
+				      * version is implemented. The appropriate
+				      * boundary class is 
+				      * @ref{HalfHyperBallBoundary}.
+				      */
+    static void half_hyper_ball (Triangulation<2> &tria,
+				 const Point<2>   &center = Point<2>(),
+				 const double      radius = 1.);
+
 				     /**
 				      * Produce a half hyper-shell,
 				      * i.e. the space between two
diff --git a/deal.II/deal.II/source/grid/grid_generator.cc b/deal.II/deal.II/source/grid/grid_generator.cc
index b3a67b781c..204b400ef7 100644
--- a/deal.II/deal.II/source/grid/grid_generator.cc
+++ b/deal.II/deal.II/source/grid/grid_generator.cc
@@ -464,6 +464,45 @@ GridGenerator::cylinder (Triangulation<2> &tria,
 
 
 
+void
+GridGenerator::half_hyper_ball (Triangulation<2> &tria,
+				const Point<2>   &p,
+				const double      radius)
+{
+				   // equilibrate cell sizes at
+				   // transition from the inner part
+				   // to the radial cells
+  const double a = 1./(1+std::sqrt(2.0));
+  const Point<2> vertices[8] = { p+Point<2>(0,-1)*radius,
+				   p+Point<2>(+1,-1)*(radius/std::sqrt(2.0)),
+				   p+Point<2>(0,-1)*(radius/std::sqrt(2.0)*a),
+				   p+Point<2>(+1,-1)*(radius/std::sqrt(2.0)*a),
+				   p+Point<2>(0,+1)*(radius/std::sqrt(2.0)*a),
+				   p+Point<2>(+1,+1)*(radius/std::sqrt(2.0)*a),
+				   p+Point<2>(0,+1)*radius,
+				   p+Point<2>(+1,+1)*(radius/std::sqrt(2.0)) };
+  
+  const int cell_vertices[5][4] = {{0, 1, 3, 2},
+				   {2, 3, 5, 4},
+				   {1, 7, 5, 3},
+				   {6, 4, 5, 7}};
+
+  std::vector<CellData<2> > cells (4, CellData<2>());
+  
+  for (unsigned int i=0; i<4; ++i) 
+    {
+      for (unsigned int j=0; j<4; ++j)
+	cells[i].vertices[j] = cell_vertices[i][j];
+      cells[i].material_id = 0;
+    };
+  
+  tria.create_triangulation (std::vector<Point<2> >(&vertices[0], &vertices[8]),
+			     cells,
+			     SubCellData());       // no boundary information
+};
+
+
+
 void
 GridGenerator::half_hyper_shell (Triangulation<2>   &tria,
 				 const Point<2>     &center,
-- 
2.39.5