-//---------------------------- 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
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
// $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
}
+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,
}
-
+//TODO: Colorize edges as circumference, left and right radius
void
GridGenerator::hyper_L (Triangulation<2> &tria,
const double a,
+//TODO: Colorize edges as circumference and cut plane
void
GridGenerator::half_hyper_ball (Triangulation<2> &tria,
const Point<2> &p,
// 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>