]> https://gitweb.dealii.org/ - dealii.git/commitdiff
new parallelogram grid in 2d
authorGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 14 Mar 2005 22:49:58 +0000 (22:49 +0000)
committerGuido Kanschat <dr.guido.kanschat@gmail.com>
Mon, 14 Mar 2005 22:49:58 +0000 (22:49 +0000)
git-svn-id: https://svn.dealii.org/trunk@10135 0785d39b-7218-0410-832d-ea1e28bc413d

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

index c45842f9afd81c401112a6e2e948f37066d0297b..9d6440d1854072f739549e4929346ad686a9d1ea 100644 (file)
@@ -1,15 +1,15 @@
-//----------------------------  grid_generator.h  ---------------------------
+//---------------------------------------------------------------------------
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+//    Copyright (C) 1998 - 2005 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.
 //
-//----------------------------  grid_generator.h  ---------------------------
+//---------------------------------------------------------------------------
 #ifndef __deal2__grid_generator_h
 #define __deal2__grid_generator_h
 
@@ -256,6 +256,28 @@ class GridGenerator
                                const Point<dim>                &p2,
                                const bool                       colorize=false);
 
+                                    /**
+                                     * A parallelogram. The first
+                                     * corner point is the
+                                     * origin. The <tt>dim</tt>
+                                     * adjacent points are the
+                                     * one-dimensional subtensors of
+                                     * the tensor provided and
+                                     * additional points will be sums
+                                     * of these two vectors.
+                                     * Colorizing is done according
+                                     * to hyper_rectangle().
+                                     *
+                                     * @note This function is
+                                     * implemented in 2d only.
+                                     */
+    template <int dim>
+    static void
+    parallelogram(Triangulation<dim>&  tria,
+                 const Tensor<2,dim>& corners,
+                 const bool           colorize=false);
+                 
+                             
                                     /**
                                      * Hypercube with a layer of
                                      * hypercubes around it. The
index d87aa11418a390c962b2874f92b730475116140a..f1787843fbdcef2216213d19a5cb5c770cc7f56d 100644 (file)
@@ -2,7 +2,7 @@
 //    $Id$
 //    Version: $Name$
 //
-//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+//    Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -159,6 +159,35 @@ void GridGenerator::hyper_cube (Triangulation<dim> &tria,
 }
 
 
+template<int dim>
+void
+GridGenerator::parallelogram (
+  Triangulation<dim>&  tria,
+  const Tensor<2,dim>& corners,
+  const bool      colorize)
+{
+  Assert (dim==2, ExcNotImplemented());
+  
+  std::vector<Point<dim> > vertices (GeometryInfo<dim>::vertices_per_cell);
+
+  vertices[1] = corners[0];
+  vertices[3] = corners[1];
+  vertices[2] = vertices[1];
+  vertices[2] += vertices[3];
+                                  // 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());
+
+                                  // Assign boundary indicators
+  if (colorize)
+    colorize_hyper_rectangle (tria);
+}
+
+
 template <int dim>
 void
 GridGenerator::subdivided_hyper_cube (Triangulation<dim> &tria,
@@ -598,7 +627,7 @@ GridGenerator::hyper_cube_slit (Triangulation<2> &tria,
 }
 
 
-
+//TODO: Colorize edges as circumference, left and right radius
 void
 GridGenerator::hyper_L (Triangulation<2> &tria,
                        const double a,
@@ -765,6 +794,7 @@ GridGenerator::cylinder (Triangulation<2> &tria,
 
 
 
+//TODO: Colorize edges as circumference and cut plane
 void
 GridGenerator::half_hyper_ball (Triangulation<2> &tria,
                                const Point<2>   &p,
@@ -1377,20 +1407,23 @@ void GridGenerator::laplace_solve (const SparseMatrix<double> &S,
 
 // explicit instantiations
 template void
-GridGenerator::hyper_rectangle<deal_II_dimension> (Triangulation<deal_II_dimension> &,
-                                                  const Point<deal_II_dimension>&,
-                                                  const Point<deal_II_dimension>&,
-                                                  const bool);
+GridGenerator::hyper_rectangle<deal_II_dimension> (
+  Triangulation<deal_II_dimension> &,
+  const Point<deal_II_dimension>&, const Point<deal_II_dimension>&,
+  const bool);
+template void
+GridGenerator::hyper_cube<deal_II_dimension> (
+  Triangulation<deal_II_dimension> &, const double, const double);
 template void
-GridGenerator::hyper_cube<deal_II_dimension> (Triangulation<deal_II_dimension> &,
-                                             const double,
-                                             const double);
+GridGenerator::parallelogram<deal_II_dimension> (
+  Triangulation<deal_II_dimension> &,
+  const Tensor<2,deal_II_dimension>&,
+  const bool);
 
 template void
-GridGenerator::subdivided_hyper_cube<deal_II_dimension> (Triangulation<deal_II_dimension> &,
-                                                         const unsigned int,
-                                                         const double,
-                                                         const double);
+GridGenerator::subdivided_hyper_cube<deal_II_dimension> (
+  Triangulation<deal_II_dimension> &,
+  const unsigned int, const double, const double);
 
 template void
 GridGenerator::subdivided_hyper_rectangle<deal_II_dimension> 

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.