]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
New function GridGenerator::moebius() creating a loop of cells which are rotated...
authorleicht <leicht@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 9 Feb 2007 06:48:31 +0000 (06:48 +0000)
committerleicht <leicht@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 9 Feb 2007 06:48:31 +0000 (06:48 +0000)
git-svn-id: https://svn.dealii.org/trunk@14441 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/grid_generator.h
deal.II/deal.II/source/grid/grid_generator.cc

index 12c4f364e210fad7a98dc900647cc6e1daa82432..ad85d444afb9211424b83ea41a7d21708a14e2a4 100644 (file)
@@ -663,6 +663,24 @@ class GridGenerator
                                  const unsigned int  n_cells = 0);
 
 
+                                    /**
+                                     * Produce a ring of cells in 3D that is
+                                     * cut open, twisted and glued together
+                                     * again. This results in a kind of
+                                     * moebius-loop.
+                                     *
+                                     * @param tria        The triangulation to be worked on.
+                                     * @param n_cells     The number of cells in the loop. Must be greater than 4.
+                                     * @param n_rotations The number of rotations (Pi/2 each) to be performed before glueing the loop together.
+                                     * @param R           The radius of the circle, which forms the middle line of the torus containing the loop of cells. Must be greater than @p r.
+                                     * @param r           The radius of the cylinder bend together as loop.
+                                     */
+    static void moebius (Triangulation<3>&  tria,
+                        const unsigned int   n_cells,
+                        const unsigned int   n_rotations,
+                        const double         R,
+                        const double         r);
+    
                                     /**
                                      * This function transformes the
                                      * @p Triangulation @p tria
index 3d6c1e3463209503be193301261e8bae57879bd5..803b82ee8c27a4d5f7e6543179e99c034534b922 100644 (file)
@@ -150,6 +150,64 @@ void GridGenerator::hyper_cube (Triangulation<dim> &tria,
 }
 
 
+#if deal_II_dimension == 3
+
+void
+GridGenerator::moebius (
+  Triangulation<3>&  tria,
+  const unsigned int   n_cells,
+  const unsigned int   n_rotations,
+  const double         R,
+  const double         r)
+{
+  const unsigned int dim=3;
+  Assert (n_cells>4, ExcMessage("More than 4 cells are needed to create a moebius grid."));
+  Assert (r>0 && R>0, ExcMessage("Outer and inner radius must be positive."));
+  Assert (R>r, ExcMessage("Outer radius must be greater than inner radius."));
+    
+  
+  std::vector<Point<dim> > vertices (4*n_cells);
+  double beta_step=n_rotations*deal_II_numbers::PI/2.0/n_cells;
+  double alpha_step=2.0*deal_II_numbers::PI/n_cells;
+  
+  for (unsigned int i=0; i<n_cells; ++i)
+    for (unsigned int j=0; j<4; ++j)
+      {
+       vertices[4*i+j][0]=R*std::cos(i*alpha_step)+r*std::cos(i*beta_step+j*deal_II_numbers::PI/2.0)*std::cos(i*alpha_step);
+       vertices[4*i+j][1]=R*std::sin(i*alpha_step)+r*std::cos(i*beta_step+j*deal_II_numbers::PI/2.0)*std::sin(i*alpha_step);
+       vertices[4*i+j][2]=r*std::sin(i*beta_step+j*deal_II_numbers::PI/2.0);
+      }
+
+  unsigned int offset=0;
+  
+  std::vector<CellData<dim> > cells (n_cells);
+  for (unsigned int i=0; i<n_cells; ++i)
+    {
+      for (unsigned int j=0; j<2; ++j)
+       {
+         cells[i].vertices[0+4*j]=offset+0+4*j;
+         cells[i].vertices[1+4*j]=offset+3+4*j;
+         cells[i].vertices[2+4*j]=offset+2+4*j;
+         cells[i].vertices[3+4*j]=offset+1+4*j;
+       }
+      offset+=4;
+      cells[i].material_id=0;
+    }
+  
+                                  // now correct the last four vertices
+  cells[n_cells-1].vertices[4]=(0+n_rotations)%4;
+  cells[n_cells-1].vertices[5]=(3+n_rotations)%4;
+  cells[n_cells-1].vertices[6]=(2+n_rotations)%4;
+  cells[n_cells-1].vertices[7]=(1+n_rotations)%4;
+
+  GridReordering<dim>::invert_all_cells_of_negative_grid(vertices,cells);
+  tria.create_triangulation_compatibility (vertices, cells, SubCellData());
+}
+
+#endif
+
+
+
 // Implementation for 2D only
 template<int dim>
 void

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.