]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Implement a cylinder shell in 3d.
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 19 Apr 2004 21:09:40 +0000 (21:09 +0000)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 19 Apr 2004 21:09:40 +0000 (21:09 +0000)
git-svn-id: https://svn.dealii.org/trunk@9066 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/grid_generator.h
deal.II/deal.II/source/grid/grid_generator.cc
tests/bits/Makefile
tests/bits/cylinder_shell_01.cc [new file with mode: 0644]

index 206d93254ad3f48834e6982b56554ba398283118..bd07bcf03e6744671e6f8f913b4580a659bff273 100644 (file)
@@ -536,6 +536,28 @@ class GridGenerator
                             const double        outer_radius,
                             const unsigned int  n_cells = 0);
 
+                                    /**
+                                     * Produce a domain that is the space
+                                     * between two cylinders in 3d, with
+                                     * given length, inner and outer radius
+                                     * and a given number of elements for
+                                     * this initial triangulation. If @p
+                                     * n_radial_cells is zero (as is the
+                                     * default), then it is computed
+                                     * adaptively such that the resulting
+                                     * elements have the least aspect
+                                     * ratio. The same holds for @p
+                                     * n_axial_cells.
+                                     *
+                                     * The triangulation needs to be void
+                                     * upon calling this function.
+                                     */
+    static void cylinder_shell (Triangulation<3>   &tria,
+                                const double        length,
+                                const double        inner_radius,
+                                const double        outer_radius,
+                                const unsigned int  n_radial_cells = 0,
+                                const unsigned int  n_axial_cells = 0);
 
                                     /**
                                      * This class produces a half hyper-ball,
index 55e392be2165c3bc92d8d1b9e2116d374ac057ef..8b1b7a10065b6526f6e686013d5dd595dfa71178 100644 (file)
@@ -1172,6 +1172,86 @@ void GridGenerator::hyper_shell (Triangulation<3>   &,
 }
 
 
+
+void GridGenerator::cylinder_shell (Triangulation<3>   &tria,
+                                    const double        length,
+                                    const double        inner_radius,
+                                    const double        outer_radius,
+                                    const unsigned int  n_radial_cells,
+                                    const unsigned int  n_axial_cells)
+{
+  Assert ((inner_radius > 0) && (inner_radius < outer_radius),
+         ExcInvalidRadii ());
+  
+  const double pi = deal_II_numbers::PI;
+  
+                                  // determine the number of cells
+                                  // for the grid. if not provided by
+                                  // the user determine it such that
+                                  // the length of each cell on the
+                                  // median (in the middle between
+                                  // the two circles) is equal to its
+                                  // radial extent (which is the
+                                  // difference between the two
+                                  // radii)
+  const unsigned int N_r = (n_radial_cells == 0 ?
+                            static_cast<unsigned int>
+                            (std::ceil((2*pi* (outer_radius + inner_radius)/2) /
+                                       (outer_radius - inner_radius))) :
+                            n_radial_cells);
+  const unsigned int N_z = (n_axial_cells == 0 ?
+                            static_cast<unsigned int>
+                            (std::ceil (length /
+                                        (2*pi*(outer_radius + inner_radius)/2/N_r))) :
+                            n_axial_cells);
+
+                                  // set up N vertices on the
+                                  // outer and N vertices on
+                                  // the inner circle. the
+                                  // first N ones are on the
+                                  // outer one, and all are
+                                  // numbered counter-clockwise
+  std::vector<Point<2> > vertices_2d(2*N_r);
+  for (unsigned int i=0; i<N_r; ++i)
+    {
+      vertices_2d[i] = Point<2>( std::cos(2*pi * i/N_r),
+                                 std::sin(2*pi * i/N_r)) * outer_radius;
+      vertices_2d[i+N_r] = vertices_2d[i] * (inner_radius/outer_radius);
+    };
+
+  std::vector<Point<3> > vertices_3d;
+  vertices_3d.reserve (2*N_r*(N_z+1));
+  for (unsigned int j=0; j<=N_z; ++j)
+    for (unsigned int i=0; i<2*N_r; ++i)
+      {
+        const Point<3> v (vertices_2d[i][0],
+                          vertices_2d[i][1],
+                          j*length/N_z);
+        vertices_3d.push_back (v);
+      }
+                            
+  std::vector<CellData<3> > cells (N_r*N_z, CellData<3>());
+  
+  for (unsigned int j=0; j<N_z; ++j)
+    for (unsigned int i=0; i<N_r; ++i) 
+      {
+        cells[i+j*N_r].vertices[0] = i + j*2*N_r;
+        cells[i+j*N_r].vertices[1] = (i+1)%N_r + j*2*N_r;
+        cells[i+j*N_r].vertices[2] = N_r+((i+1)%N_r) + j*2*N_r;
+        cells[i+j*N_r].vertices[3] = N_r+i + j*2*N_r;
+        
+        cells[i+j*N_r].vertices[4] = i + (j+1)*2*N_r;
+        cells[i+j*N_r].vertices[5] = (i+1)%N_r + (j+1)*2*N_r;
+        cells[i+j*N_r].vertices[6] = N_r+((i+1)%N_r) + (j+1)*2*N_r;
+        cells[i+j*N_r].vertices[7] = N_r+i + (j+1)*2*N_r;
+        
+        cells[i+j*N_r].material_id = 0;
+      }
+  
+  tria.create_triangulation (vertices_3d, cells, SubCellData());
+}
+
+
 #endif
 
 
index 37e954b547f5bccd9ead3471ff24b1433e802898..1726bb96c3f441b6b0c946c92b864a12404ec251 100644 (file)
@@ -57,7 +57,8 @@ tests_x = anna_? \
          subdomain_ids_* \
          dof_constraints_* \
          apply_boundary_values_* \
-         error_estimator_*
+         error_estimator_* \
+         cylinder_shell_*
 
 ifeq ($(USE_CONTRIB_PETSC),yes)
   tests_x += petsc_*
diff --git a/tests/bits/cylinder_shell_01.cc b/tests/bits/cylinder_shell_01.cc
new file mode 100644 (file)
index 0000000..9517325
--- /dev/null
@@ -0,0 +1,44 @@
+//----------------------------  cylinder_shell_01.cc  ---------------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2003, 2004 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------------------  cylinder_shell_01.cc  ---------------------------
+
+// test that the grid generated by GridGenerator::cylinder_shell<3> works as
+// expected
+
+
+#include "../tests.h"
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <grid/tria_iterator.h>
+#include <grid/tria_accessor.h>
+#include <grid/grid_generator.h>
+#include <fstream>
+
+    
+
+int main () 
+{
+  std::ofstream logfile("cylinder_shell_01.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  logfile.precision (2);
+
+                                   // generate a hyperball in 3d
+  Triangulation<3> tria;
+  GridGenerator::cylinder_shell (tria, 1, .8, 1);
+
+                                   // make sure that all cells have positive
+                                   // volume
+  for (Triangulation<3>::active_cell_iterator cell=tria.begin_active();
+       cell!=tria.end(); ++cell)
+    deallog << cell << ' ' << cell->measure () << std::endl;
+}

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.