Add GridGenerator::eccentric_hyper_shell

Add an eccentric hyper shell geometry that is formed by two spheres
with different center points.

diff --git a/source/grid/grid_generator.cc b/source/grid/grid_generator.cc
index 83ef28df0882626ded266d78e4f641eb36d3d2ac..c42eb2acbfc7ca4cf3e06e397a3331b8a5525f5e 100644 (file)
--- a/source/grid/grid_generator.cc
+++ b/source/grid/grid_generator.cc
@@ -3003,6 +3003,58 @@ namespace GridGenerator
     tria.set_manifold(0, SphericalManifold<2>(center));
   }
 
+  template <int dim>
+  void
+  eccentric_hyper_shell(Triangulation<dim> &tria,
+                        const Point<dim> &  inner_center,
+                        const Point<dim> &  outer_center,
+                        const double        inner_radius,
+                        const double        outer_radius,
+                        const unsigned int  n_cells)
+  {
+    GridGenerator::hyper_shell(
+      tria, outer_center, inner_radius, outer_radius, n_cells, true);
+
+    // check the consistency of the dimensions provided
+    Assert(
+      outer_radius - inner_radius > outer_center.distance(inner_center),
+      ExcInternalError(
+        "Inner radius superior or equal to outer radius plus eccentricity."));
+
+    // shift nodes along the inner boundary according to the position of
+    // inner_circle
+    std::set<Point<dim> *> vertices_to_move;
+
+    for (const auto &cell : tria.active_cell_iterators())
+      for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
+        if (cell->face(f)->boundary_id() == 0)
+          for (unsigned int v = 0; v < GeometryInfo<dim>::vertices_per_face;
+               ++v)
+            vertices_to_move.insert(&cell->face(f)->vertex(v));
+
+    const auto shift = inner_center - outer_center;
+    for (const auto &p : vertices_to_move)
+      (*p) += shift;
+
+    // the original hyper_shell function assigns the same manifold id
+    // to all cells and faces. Set all manifolds ids to a different
+    // value (2), then use boundary ids to assign different manifolds to
+    // the inner (0) and outer manifolds (1). Use a transfinite manifold
+    // for all faces and cells aside from the boundaries.
+    tria.set_all_manifold_ids(2);
+    GridTools::copy_boundary_to_manifold_id(tria);
+
+    SphericalManifold<dim> inner_manifold(inner_center);
+    SphericalManifold<dim> outer_manifold(outer_center);
+
+    TransfiniteInterpolationManifold<dim> transfinite;
+    transfinite.initialize(tria);
+
+    tria.set_manifold(0, inner_manifold);
+    tria.set_manifold(1, outer_manifold);
+    tria.set_manifold(2, transfinite);
+  }
+
 
 
   template <int dim>