]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add Boundary::normal_vector. Implement it. Test it.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 25 Aug 2011 15:51:59 +0000 (15:51 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 25 Aug 2011 15:51:59 +0000 (15:51 +0000)
git-svn-id: https://svn.dealii.org/trunk@24197 0785d39b-7218-0410-832d-ea1e28bc413d

16 files changed:
deal.II/include/deal.II/grid/tria_boundary.h
deal.II/include/deal.II/grid/tria_boundary_lib.h
deal.II/source/grid/tria_boundary.cc
deal.II/source/grid/tria_boundary_lib.cc
tests/deal.II/normal_vector_01.cc [new file with mode: 0644]
tests/deal.II/normal_vector_01/cmp/generic [new file with mode: 0644]
tests/deal.II/normal_vector_01_2d.cc [new file with mode: 0644]
tests/deal.II/normal_vector_01_2d/cmp/generic [new file with mode: 0644]
tests/deal.II/normal_vector_02.cc [new file with mode: 0644]
tests/deal.II/normal_vector_02/cmp/generic [new file with mode: 0644]
tests/deal.II/normal_vector_02_2d.cc [new file with mode: 0644]
tests/deal.II/normal_vector_02_2d/cmp/generic [new file with mode: 0644]
tests/deal.II/normal_vector_03.cc [new file with mode: 0644]
tests/deal.II/normal_vector_03/cmp/generic [new file with mode: 0644]
tests/deal.II/normal_vector_03_2d.cc [new file with mode: 0644]
tests/deal.II/normal_vector_03_2d/cmp/generic [new file with mode: 0644]

index 046fdf3066a67d8306022eb08197243f8022eed1..1b2871ca11bd1b0d520e4ae8169f2859c474e1e6 100644 (file)
@@ -122,7 +122,8 @@ class Boundary : public Subscriptor
                                      * lines therefore is also on the
                                      * boundary).
                                      */
-    virtual Point<spacedim>
+    virtual 
+    Point<spacedim>
     get_new_point_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line) const = 0;
 
                                     /**
@@ -151,7 +152,8 @@ class Boundary : public Subscriptor
                                      * default implementation throws
                                      * an error in any case, however.
                                      */
-    virtual Point<spacedim>
+    virtual 
+    Point<spacedim>
     get_new_point_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad) const;
 
                                     /**
@@ -193,7 +195,8 @@ class Boundary : public Subscriptor
                                      * default implementation throws
                                      * an error in any case, however.
                                      */
-    virtual void
+    virtual
+    void
     get_intermediate_points_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line,
                                     std::vector<Point<spacedim> > &points) const;
     
@@ -227,7 +230,8 @@ class Boundary : public Subscriptor
                                      * default implementation throws
                                      * an error in any case, however.
                                      */
-    virtual void
+    virtual
+    void
     get_intermediate_points_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad,
                                     std::vector<Point<spacedim> > &points) const;
 
@@ -248,6 +252,46 @@ class Boundary : public Subscriptor
     get_intermediate_points_on_face (const typename Triangulation<dim,spacedim>::face_iterator &face,
                                     std::vector<Point<spacedim> > &points) const;
 
+                                    /**
+                                     * Return the normal vector to the surface
+                                     * at the point p. If p is not in fact
+                                     * on the surface, but only closeby,
+                                     * try to return something reasonable,
+                                     * for example the normal vector
+                                     * at the surface point closest to p.
+                                     * (The point p will in fact not normally
+                                     * lie on the actual surface, but rather
+                                     * be a quadrature point mapped by some
+                                     * polynomial mapping; the mapped surface,
+                                     * however, will not usually coincide with
+                                     * the actual surface.)
+                                     * 
+                                     * The face iterator gives an indication
+                                     * which face this function is supposed
+                                     * to compute the normal vector for.
+                                     * This is useful if the boundary of
+                                     * the domain is composed of different
+                                     * nondifferential pieces (for example
+                                     * when using the StraightBoundary class
+                                     * to approximate a geometry that is
+                                     * completely described by the coarse mesh,
+                                     * with piecewise (bi-)linear components
+                                     * between the vertices, but where the
+                                     * boundary may have a kink at the vertices
+                                     * itself).
+                                     * 
+                                     * @note Implementations of this function 
+                                     * should be able to assume that the point p
+                                     * lies within or close to the face described by the
+                                     * first argument. In turn, callers of this
+                                     * function should ensure that this is
+                                     * in fact the case.
+                                     */
+    virtual
+    Tensor<1,spacedim>
+    normal_vector (const typename Triangulation<dim,spacedim>::face_iterator &face,
+                  const Point<spacedim> &p) const;
+    
                                     /**
                                      * Compute the normal vectors to
                                      * the boundary at each vertex of
@@ -276,7 +320,8 @@ class Boundary : public Subscriptor
                                      * with respect to points inside
                                      * the given face.
                                      */
-    virtual void
+    virtual
+    void
     get_normals_at_vertices (const typename Triangulation<dim,spacedim>::face_iterator &face,
                             FaceVertexNormals &face_vertex_normals) const;
 
@@ -414,7 +459,8 @@ class StraightBoundary : public Boundary<dim,spacedim>
                                      * base class for more
                                      * information.
                                      */
-    virtual Point<spacedim>
+    virtual
+    Point<spacedim>
     get_new_point_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad) const;
 
                                     /**
@@ -429,7 +475,8 @@ class StraightBoundary : public Boundary<dim,spacedim>
                                      * and the documentation of the
                                      * base class.
                                      */
-    virtual void
+    virtual
+    void
     get_intermediate_points_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line,
                                     std::vector<Point<spacedim> > &points) const;
 
@@ -445,11 +492,26 @@ class StraightBoundary : public Boundary<dim,spacedim>
                                      * and the documentation of the
                                      * base class.
                                      */
-    virtual void
+    virtual
+    void
     get_intermediate_points_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad,
                                     std::vector<Point<spacedim> > &points) const;
 
                                     /**
+                                     * Implementation of the function
+                                     * declared in the base class.
+                                     *
+                                     * Refer to the general
+                                     * documentation of this class
+                                     * and the documentation of the
+                                     * base class.
+                                     */
+    virtual
+    Tensor<1,spacedim>
+    normal_vector (const typename Triangulation<dim,spacedim>::face_iterator &face,
+                  const Point<spacedim> &p) const;
+
+                                    /**
                                      * Compute the normals to the
                                      * boundary at the vertices of
                                      * the given face.
@@ -459,7 +521,8 @@ class StraightBoundary : public Boundary<dim,spacedim>
                                      * and the documentation of the
                                      * base class.
                                      */
-    virtual void
+    virtual
+    void
     get_normals_at_vertices (const typename Triangulation<dim,spacedim>::face_iterator &face,
                             typename Boundary<dim,spacedim>::FaceVertexNormals &face_vertex_normals) const;
 
index 178f5983c09d908bfa257bf26bfa5c6b730e9828..eb5f68f804063c5cfc06a84e0b10b1a2e0e87869 100644 (file)
@@ -388,7 +388,8 @@ class HyperBallBoundary : public StraightBoundary<dim,spacedim>
                                      * this class and the documentation of the
                                      * base class.
                                      */
-    virtual Point<spacedim>
+    virtual 
+    Point<spacedim>
     get_new_point_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line) const;
 
                                     /**
@@ -396,7 +397,8 @@ class HyperBallBoundary : public StraightBoundary<dim,spacedim>
                                      * this class and the documentation of the
                                      * base class.
                                      */
-    virtual Point<spacedim>
+    virtual 
+    Point<spacedim>
     get_new_point_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad) const;
 
                                     /**
@@ -408,7 +410,8 @@ class HyperBallBoundary : public StraightBoundary<dim,spacedim>
                                      * Calls
                                      * @p get_intermediate_points_between_points.
                                      */
-    virtual void
+    virtual 
+    void
     get_intermediate_points_on_line (const typename Triangulation<dim,spacedim>::line_iterator &line,
                                     std::vector<Point<spacedim> > &points) const;
 
@@ -421,11 +424,26 @@ class HyperBallBoundary : public StraightBoundary<dim,spacedim>
                                      * Only implemented for <tt>dim=3</tt>
                                      * and for <tt>points.size()==1</tt>.
                                      */
-    virtual void
+    virtual 
+    void
     get_intermediate_points_on_quad (const typename Triangulation<dim,spacedim>::quad_iterator &quad,
                                     std::vector<Point<spacedim> > &points) const;
 
                                     /**
+                                     * Implementation of the function
+                                     * declared in the base class.
+                                     *
+                                     * Refer to the general
+                                     * documentation of this class
+                                     * and the documentation of the
+                                     * base class.
+                                     */
+    virtual
+    Tensor<1,spacedim>
+    normal_vector (const typename Triangulation<dim,spacedim>::face_iterator &face,
+                  const Point<spacedim> &p) const;
+
+                                    /**
                                      * Compute the normals to the
                                      * boundary at the vertices of
                                      * the given face.
@@ -435,19 +453,22 @@ class HyperBallBoundary : public StraightBoundary<dim,spacedim>
                                      * and the documentation of the
                                      * base class.
                                      */
-    virtual void
+    virtual 
+    void
     get_normals_at_vertices (const typename Triangulation<dim,spacedim>::face_iterator &face,
                             typename Boundary<dim,spacedim>::FaceVertexNormals &face_vertex_normals) const;
 
                                     /**
                                      * Return the center of the ball.
                                      */
-    Point<spacedim> get_center () const;
+    Point<spacedim>
+    get_center () const;
 
                                     /**
                                      * Return the radius of the ball.
                                      */
-    double get_radius () const;
+    double 
+    get_radius () const;
 
                                     /**
                                      * Exception. Thrown by the
index 6b4a049a11298bca8817c06b8a18dddbec9b3312..91b7f1e21360d6e3db93c8956e5ab31d95e5fddd 100644 (file)
@@ -17,6 +17,7 @@
 #include <deal.II/grid/tria.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/tria_accessor.h>
+#include <deal.II/fe/fe_q.h>
 #include <cmath>
 
 DEAL_II_NAMESPACE_OPEN
@@ -150,6 +151,18 @@ get_intermediate_points_on_face (const Triangulation<1,2>::face_iterator &,
 
 
 
+template <int dim, int spacedim>
+Tensor<1,spacedim>
+Boundary<dim, spacedim>::
+normal_vector (const typename Triangulation<dim, spacedim>::face_iterator &,
+              const Point<spacedim> &) const
+{
+  Assert (false, ExcPureFunctionCalled());
+  return Tensor<1,spacedim>();
+}
+
+
+
 template <int dim, int spacedim>
 void
 Boundary<dim, spacedim>::
@@ -472,6 +485,157 @@ get_intermediate_points_on_quad (const Triangulation<2,3>::quad_iterator &quad,
 
 
 
+template <>
+Tensor<1,1>
+StraightBoundary<1,1>::
+normal_vector (const typename Triangulation<1,1>::face_iterator &,
+              const Point<1> &) const
+{
+  Assert (false, ExcNotImplemented());
+  return Tensor<1,1>();
+}
+
+
+template <>
+Tensor<1,2>
+StraightBoundary<1,2>::
+normal_vector (const typename Triangulation<1,2>::face_iterator &,
+              const Point<2> &) const
+{
+  Assert (false, ExcNotImplemented());
+  return Tensor<1,2>();
+}
+
+
+
+namespace internal
+{
+  namespace
+  {
+    /**
+     * Compute the normalized cross product of a set of dim-1 basis
+     * vectors.
+     */
+    Tensor<1,2>
+    normalized_alternating_product (const Tensor<1,2> (&basis_vectors)[1])
+    {
+      Tensor<1,2> tmp;
+      cross_product (tmp, basis_vectors[0]);
+      return tmp/tmp.norm();
+    }
+
+
+
+    Tensor<1,3>
+    normalized_alternating_product (const Tensor<1,3> (&basis_vectors)[1])
+    {
+      // we get here from StraightBoundary<2,3>::normal_vector, but
+      // the implementation below is bogus for this case anyway
+      // (see the assert at the beginning of that function).
+      Assert (false, ExcNotImplemented());
+      return Tensor<1,3>();
+    }
+
+
+
+    Tensor<1,3>
+    normalized_alternating_product (const Tensor<1,3> (&basis_vectors)[2])
+    {
+      Tensor<1,3> tmp;
+      cross_product (tmp, basis_vectors[0], basis_vectors[1]);
+      return tmp/tmp.norm();
+    }
+
+  }
+}
+
+
+template <int dim, int spacedim>
+Tensor<1,spacedim>
+StraightBoundary<dim,spacedim>::
+normal_vector (const typename Triangulation<dim,spacedim>::face_iterator &face,
+              const Point<spacedim> &p) const
+{
+  // I don't think the implementation below will work when dim!=spacedim;
+  // in fact, I believe that we don't even have enough information here,
+  // because we would need to know not only about the tangent vectors
+  // of the face, but also of the cell, to compute the normal vector.
+  // Someone will have to think about this some more.
+  Assert (dim == spacedim, ExcNotImplemented());
+
+  // in order to find out what the normal vector is, we first need to
+  // find the reference coordinates of the point p on the given face,
+  // or at least the reference coordinates of the closest point on the
+  // face
+  //
+  // in other words, we need to find a point xi so that f(xi)=||F(xi)-p||^2->min
+  // where F(xi) is the mapping. this algorithm is implemented in
+  // MappingQ1<dim,spacedim>::transform_real_to_unit_cell but only for cells,
+  // while we need it for faces here. it's also implemented in somewhat
+  // more generality there using the machinery of the MappingQ1 class
+  // while we really only need it for a specific case here
+  //
+  // in any case, the iteration we use here is a Gauss-Newton's iteration with
+  //   xi^{n+1} = xi^n - H(xi^n)^{-1} J(xi^n)
+  // where
+  //   J(xi) = (grad F(xi))^T (F(xi)-p)
+  // and
+  //   H(xi) = [grad F(xi)]^T [grad F(xi)]
+  // In all this,
+  //   F(xi) = sum_v vertex[v] phi_v(xi)
+  // We get the shape functions phi_v from an object of type FE_Q<dim-1>(1)
+
+  // we start with the point xi=1/2, xi=(1/2,1/2), ...
+  const unsigned int facedim = dim-1;
+
+  Point<facedim> xi;
+  for (unsigned int i=0; i<facedim; ++i)
+    xi[i] = 1./2;
+
+  FE_Q<facedim> linear_fe(1);
+
+  const double eps = 1e-12;
+  Tensor<1,spacedim> grad_F[facedim];
+  while (true)
+  {
+    Point<spacedim> F;
+    for (unsigned int v=0; v<GeometryInfo<facedim>::vertices_per_cell; ++v)
+      F += face->vertex(v) * linear_fe.shape_value(v, xi);
+
+    for (unsigned int i=0; i<facedim; ++i)
+    {
+      grad_F[i] = 0;
+      for (unsigned int v=0; v<GeometryInfo<facedim>::vertices_per_cell; ++v)
+       grad_F[i] += face->vertex(v) * linear_fe.shape_grad(v, xi)[i];
+    }
+
+    Tensor<1,facedim> J;
+    for (unsigned int i=0; i<facedim; ++i)
+      for (unsigned int j=0; j<spacedim; ++j)
+        J[i] += grad_F[i][j] * (F-p)[j];
+
+    Tensor<2,facedim> H;
+    for (unsigned int i=0; i<facedim; ++i)
+      for (unsigned int j=0; j<facedim; ++j)
+       for (unsigned int k=0; k<spacedim; ++k)
+      H[i][j] += grad_F[i][k] * grad_F[j][k];
+
+    const Point<facedim> delta_xi = -invert(H) * J;
+    xi += delta_xi;
+
+    if (delta_xi.norm() < eps)
+      break;
+  }
+
+  // so now we have the reference coordinates xi of the point p.
+  // we then have to compute the normal vector, which we can do
+  // by taking the (normalize) alternating product of all the tangent
+  // vectors given by grad_F
+  return internal::normalized_alternating_product(grad_F);
+}
+
+
+
 template <>
 void
 StraightBoundary<1>::
index 774cc39dc005ceca4986284333213f1262300159..2cf45eac6409940316b539c11dfa880a1b70780c 100644 (file)
@@ -515,6 +515,7 @@ get_new_point_on_quad (const Triangulation<1,1>::quad_iterator &) const
   return Point<1>();
 }
 
+
 template <>
 Point<2>
 HyperBallBoundary<1,2>::
@@ -576,6 +577,7 @@ HyperBallBoundary<dim,spacedim>::get_intermediate_points_on_line (
 }
 
 
+
 template <int dim, int spacedim>
 void
 HyperBallBoundary<dim,spacedim>::get_intermediate_points_between_points (
@@ -732,6 +734,18 @@ HyperBallBoundary<dim,spacedim>::get_intermediate_points_on_quad (
 
 
 
+template <int dim, int spacedim>
+Tensor<1,spacedim>
+HyperBallBoundary<dim,spacedim>::
+normal_vector (const typename Triangulation<dim,spacedim>::face_iterator &,
+              const Point<spacedim> &p) const
+{
+  const Tensor<1,spacedim> unnormalized_normal = p-center;
+  return unnormalized_normal/unnormalized_normal.norm();
+}
+
+
+
 template <>
 void
 HyperBallBoundary<1>::
diff --git a/tests/deal.II/normal_vector_01.cc b/tests/deal.II/normal_vector_01.cc
new file mode 100644 (file)
index 0000000..6d8925e
--- /dev/null
@@ -0,0 +1,92 @@
+//----------------------------  normal_vector_01.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2005, 2008, 2011 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.
+//
+//----------------------------  normal_vector_01.cc  ---------------------------
+
+
+// test that at the vertices, Boundary::normal_vector returns the same as
+// Boundary::get_normals_at_vertices once the latter vectors are normalized
+
+
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_boundary.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+
+void create_triangulation(const unsigned int case_no,
+                         Triangulation<3> &tria)
+{
+  switch (case_no)
+    {
+      case 0:
+           GridGenerator::hyper_cube(tria, 1., 3.);
+           break;
+      case 1:
+      {
+       GridGenerator::hyper_cube(tria, 1., 3.);
+       Point<3> &v0=tria.begin_active()->vertex(0);
+       v0 = Point<3> (0,-0.5,-1);
+       Point<3> &v1=tria.begin_active()->vertex(1);
+       v1 = Point<3> (1.25, 0.25, 0.25);
+       break;
+      }
+      default:
+           Assert(false, ExcNotImplemented());
+    };
+}
+
+
+
+int main ()
+{
+  std::ofstream logfile ("normal_vector_01/output");
+  deallog << std::setprecision (3);
+  deallog << std::fixed;
+  deallog.attach(logfile);
+  deallog.depth_console (0);
+
+  Triangulation<3> tria;
+  StraightBoundary<3> boundary;
+  Boundary<3>::FaceVertexNormals normals;
+  for (unsigned int case_no=0; case_no<2; ++case_no)
+    {
+      deallog << "Case" << case_no << std::endl;
+      create_triangulation(case_no, tria);
+      const Triangulation<3>::active_cell_iterator cell=tria.begin_active();
+      Triangulation<3>::face_iterator face;
+      for (unsigned int face_no=0; face_no<GeometryInfo<3>::faces_per_cell; ++face_no)
+       {
+         face=cell->face(face_no);
+         boundary.get_normals_at_vertices(face, normals);
+         for (unsigned int v=0; v<GeometryInfo<3>::vertices_per_face; ++v)
+           Assert ((boundary.normal_vector (face,
+                                            face->vertex(v))
+                    -
+                    normals[v] / normals[v].norm()).norm()
+                   <
+                   1e-12,
+                   ExcInternalError());
+       }
+      tria.clear();
+    }
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/deal.II/normal_vector_01/cmp/generic b/tests/deal.II/normal_vector_01/cmp/generic
new file mode 100644 (file)
index 0000000..dfea164
--- /dev/null
@@ -0,0 +1,4 @@
+
+DEAL::Case0
+DEAL::Case1
+DEAL::OK
diff --git a/tests/deal.II/normal_vector_01_2d.cc b/tests/deal.II/normal_vector_01_2d.cc
new file mode 100644 (file)
index 0000000..39b0a7f
--- /dev/null
@@ -0,0 +1,92 @@
+//----------------------------  normal_vector_01_2d.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2005, 2008, 2011 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.
+//
+//----------------------------  normal_vector_01_2d.cc  ---------------------------
+
+
+// test that at the vertices, Boundary::normal_vector returns the same as
+// Boundary::get_normals_at_vertices once the latter vectors are normalized
+
+
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_boundary.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+
+void create_triangulation(const unsigned int case_no,
+                         Triangulation<2> &tria)
+{
+  switch (case_no)
+    {
+      case 0:
+           GridGenerator::hyper_cube(tria, 1., 3.);
+           break;
+      case 1:
+      {
+       GridGenerator::hyper_cube(tria, 1., 3.);
+       Point<2> &v0=tria.begin_active()->vertex(0);
+       v0 = Point<2> (-0.5,-1);
+       Point<2> &v1=tria.begin_active()->vertex(1);
+       v1 = Point<2> (0.25, 0.25);
+       break;
+      }
+      default:
+           Assert(false, ExcNotImplemented());
+    };
+}
+
+
+
+int main ()
+{
+  std::ofstream logfile ("normal_vector_01_2d/output");
+  deallog << std::setprecision (3);
+  deallog << std::fixed;
+  deallog.attach(logfile);
+  deallog.depth_console (0);
+
+  Triangulation<2> tria;
+  StraightBoundary<2> boundary;
+  Boundary<2>::FaceVertexNormals normals;
+  for (unsigned int case_no=0; case_no<2; ++case_no)
+    {
+      deallog << "Case" << case_no << std::endl;
+      create_triangulation(case_no, tria);
+      const Triangulation<2>::active_cell_iterator cell=tria.begin_active();
+      Triangulation<2>::face_iterator face;
+      for (unsigned int face_no=0; face_no<GeometryInfo<2>::faces_per_cell; ++face_no)
+       {
+         face=cell->face(face_no);
+         boundary.get_normals_at_vertices(face, normals);
+         for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_face; ++v)
+           Assert ((boundary.normal_vector (face,
+                                            face->vertex(v))
+                    -
+                    normals[v] / normals[v].norm()).norm()
+                   <
+                   1e-12,
+                   ExcInternalError());
+       }
+      tria.clear();
+    }
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/deal.II/normal_vector_01_2d/cmp/generic b/tests/deal.II/normal_vector_01_2d/cmp/generic
new file mode 100644 (file)
index 0000000..dfea164
--- /dev/null
@@ -0,0 +1,4 @@
+
+DEAL::Case0
+DEAL::Case1
+DEAL::OK
diff --git a/tests/deal.II/normal_vector_02.cc b/tests/deal.II/normal_vector_02.cc
new file mode 100644 (file)
index 0000000..15c9f26
--- /dev/null
@@ -0,0 +1,65 @@
+//----------------------------  normal_vector_02.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2005, 2008, 2011 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.
+//
+//----------------------------  normal_vector_02.cc  ---------------------------
+
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_boundary.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+
+
+int main ()
+{
+  std::ofstream logfile ("normal_vector_02/output");
+  deallog << std::setprecision (3);
+  deallog << std::fixed;
+  deallog.attach(logfile);
+  deallog.depth_console (0);
+
+  HyperBallBoundary<3> boundary (Point<3>(1,0,0));
+
+  Triangulation<3> tria;
+  Boundary<3>::FaceVertexNormals normals;
+
+  GridGenerator::hyper_ball (tria, Point<3>(1,0,0), 3);
+
+  Triangulation<3>::active_cell_iterator cell=tria.begin_active();
+  for (; cell!=tria.end(); ++cell)
+    for (unsigned int face_no=0;
+        face_no<GeometryInfo<3>::faces_per_cell; ++face_no)
+      if (cell->at_boundary(face_no))
+       {
+         Triangulation<3>::face_iterator face = cell->face(face_no);
+         boundary.get_normals_at_vertices(face, normals);
+         for (unsigned int v=0; v<GeometryInfo<3>::vertices_per_face; ++v)
+           Assert ((boundary.normal_vector (face,
+                                            face->vertex(v))
+                    -
+                    normals[v] / normals[v].norm()).norm()
+                   <
+                   1e-12,
+                   ExcInternalError());
+       }
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/deal.II/normal_vector_02/cmp/generic b/tests/deal.II/normal_vector_02/cmp/generic
new file mode 100644 (file)
index 0000000..0fd8fc1
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::OK
diff --git a/tests/deal.II/normal_vector_02_2d.cc b/tests/deal.II/normal_vector_02_2d.cc
new file mode 100644 (file)
index 0000000..6366127
--- /dev/null
@@ -0,0 +1,65 @@
+//----------------------------  normal_vector_02_2d.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2005, 2008, 2011 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.
+//
+//----------------------------  normal_vector_02_2d.cc  ---------------------------
+
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_boundary.h>
+#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+
+
+int main ()
+{
+  std::ofstream logfile ("normal_vector_02_2d/output");
+  deallog << std::setprecision (3);
+  deallog << std::fixed;
+  deallog.attach(logfile);
+  deallog.depth_console (0);
+
+  HyperBallBoundary<2> boundary (Point<2>(1,0));
+
+  Triangulation<2> tria;
+  Boundary<2>::FaceVertexNormals normals;
+
+  GridGenerator::hyper_ball (tria, Point<2>(1,0), 3);
+
+  Triangulation<2>::active_cell_iterator cell=tria.begin_active();
+  for (; cell!=tria.end(); ++cell)
+    for (unsigned int face_no=0;
+        face_no<GeometryInfo<2>::faces_per_cell; ++face_no)
+      if (cell->at_boundary(face_no))
+       {
+         Triangulation<2>::face_iterator face = cell->face(face_no);
+         boundary.get_normals_at_vertices(face, normals);
+         for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_face; ++v)
+           Assert ((boundary.normal_vector (face,
+                                            face->vertex(v))
+                    -
+                    normals[v] / normals[v].norm()).norm()
+                   <
+                   1e-12,
+                   ExcInternalError());
+       }
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/deal.II/normal_vector_02_2d/cmp/generic b/tests/deal.II/normal_vector_02_2d/cmp/generic
new file mode 100644 (file)
index 0000000..0fd8fc1
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::OK
diff --git a/tests/deal.II/normal_vector_03.cc b/tests/deal.II/normal_vector_03.cc
new file mode 100644 (file)
index 0000000..00f04aa
--- /dev/null
@@ -0,0 +1,114 @@
+//----------------------------  normal_vector_03.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2005, 2008, 2011 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.
+//
+//----------------------------  normal_vector_03.cc  ---------------------------
+
+
+// like _01 but instead for random points on the face. note that the normal
+// vector on a bilinearly mapped face is also a bilinear function
+
+
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_boundary.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/fe/fe_q.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+
+void create_triangulation(const unsigned int case_no,
+                         Triangulation<3> &tria)
+{
+  switch (case_no)
+    {
+      case 0:
+           GridGenerator::hyper_cube(tria, 1., 3.);
+           break;
+      case 1:
+      {
+       GridGenerator::hyper_cube(tria, 1., 3.);
+       Point<3> &v0=tria.begin_active()->vertex(0);
+       v0 = Point<3> (0,-0.5,-1);
+       Point<3> &v1=tria.begin_active()->vertex(1);
+       v1 = Point<3> (1.25, 0.25, 0.25);
+       break;
+      }
+      default:
+           Assert(false, ExcNotImplemented());
+    };
+}
+
+
+
+int main ()
+{
+  std::ofstream logfile ("normal_vector_03/output");
+  deallog << std::setprecision (3);
+  deallog << std::fixed;
+  deallog.attach(logfile);
+  deallog.depth_console (0);
+
+  FE_Q<2> linear_interpolator(1);
+
+  Triangulation<3> tria;
+  StraightBoundary<3> boundary;
+  Boundary<3>::FaceVertexNormals normals;
+  for (unsigned int case_no=0; case_no<2; ++case_no)
+    {
+      deallog << "Case" << case_no << std::endl;
+      create_triangulation(case_no, tria);
+      const Triangulation<3>::active_cell_iterator cell=tria.begin_active();
+      Triangulation<3>::face_iterator face;
+      for (unsigned int face_no=0; face_no<GeometryInfo<3>::faces_per_cell; ++face_no)
+       {
+         face=cell->face(face_no);
+         boundary.get_normals_at_vertices(face, normals);
+
+         for (double xi=0; xi<=1; xi+=0.234)
+           for (double eta=0; eta<=1; eta+=0.234)
+             {
+               Point<3> p;
+               Tensor<1,3> normal;
+
+               for (unsigned int v=0; v<GeometryInfo<3>::vertices_per_face; ++v)
+                 {
+                   p += face->vertex(v) * linear_interpolator.shape_value(v,Point<2>(xi,eta));
+                   normal += normals[v] *
+                             linear_interpolator.shape_value(v,Point<2>(xi,eta));
+                 }
+               normal /= normal.norm();
+
+               deallog << "p=" << p
+                       << ", n=" << boundary.normal_vector (face, p)
+                       << std::endl;
+
+               Assert ((boundary.normal_vector (face, p)
+                        -
+                       normal).norm()
+                       <
+                       1e-10,
+                       ExcInternalError());
+             }
+       }
+
+      tria.clear();
+    }
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/deal.II/normal_vector_03/cmp/generic b/tests/deal.II/normal_vector_03/cmp/generic
new file mode 100644 (file)
index 0000000..37c91d2
--- /dev/null
@@ -0,0 +1,304 @@
+
+DEAL::Case0
+DEAL::p=1.000 1.000 1.000, n=1.000 0.000 0.000
+DEAL::p=1.000 1.000 1.468, n=1.000 0.000 0.000
+DEAL::p=1.000 1.000 1.936, n=1.000 0.000 0.000
+DEAL::p=1.000 1.000 2.404, n=1.000 0.000 0.000
+DEAL::p=1.000 1.000 2.872, n=1.000 0.000 0.000
+DEAL::p=1.000 1.468 1.000, n=1.000 0.000 0.000
+DEAL::p=1.000 1.468 1.468, n=1.000 0.000 0.000
+DEAL::p=1.000 1.468 1.936, n=1.000 0.000 0.000
+DEAL::p=1.000 1.468 2.404, n=1.000 0.000 0.000
+DEAL::p=1.000 1.468 2.872, n=1.000 0.000 0.000
+DEAL::p=1.000 1.936 1.000, n=1.000 0.000 0.000
+DEAL::p=1.000 1.936 1.468, n=1.000 0.000 0.000
+DEAL::p=1.000 1.936 1.936, n=1.000 0.000 0.000
+DEAL::p=1.000 1.936 2.404, n=1.000 0.000 0.000
+DEAL::p=1.000 1.936 2.872, n=1.000 0.000 0.000
+DEAL::p=1.000 2.404 1.000, n=1.000 0.000 0.000
+DEAL::p=1.000 2.404 1.468, n=1.000 0.000 0.000
+DEAL::p=1.000 2.404 1.936, n=1.000 0.000 0.000
+DEAL::p=1.000 2.404 2.404, n=1.000 0.000 0.000
+DEAL::p=1.000 2.404 2.872, n=1.000 0.000 0.000
+DEAL::p=1.000 2.872 1.000, n=1.000 0.000 0.000
+DEAL::p=1.000 2.872 1.468, n=1.000 0.000 0.000
+DEAL::p=1.000 2.872 1.936, n=1.000 0.000 0.000
+DEAL::p=1.000 2.872 2.404, n=1.000 0.000 0.000
+DEAL::p=1.000 2.872 2.872, n=1.000 0.000 0.000
+DEAL::p=3.000 1.000 1.000, n=1.000 0.000 0.000
+DEAL::p=3.000 1.000 1.468, n=1.000 0.000 0.000
+DEAL::p=3.000 1.000 1.936, n=1.000 0.000 0.000
+DEAL::p=3.000 1.000 2.404, n=1.000 0.000 0.000
+DEAL::p=3.000 1.000 2.872, n=1.000 0.000 0.000
+DEAL::p=3.000 1.468 1.000, n=1.000 0.000 0.000
+DEAL::p=3.000 1.468 1.468, n=1.000 0.000 0.000
+DEAL::p=3.000 1.468 1.936, n=1.000 0.000 0.000
+DEAL::p=3.000 1.468 2.404, n=1.000 0.000 0.000
+DEAL::p=3.000 1.468 2.872, n=1.000 0.000 0.000
+DEAL::p=3.000 1.936 1.000, n=1.000 0.000 0.000
+DEAL::p=3.000 1.936 1.468, n=1.000 0.000 0.000
+DEAL::p=3.000 1.936 1.936, n=1.000 0.000 0.000
+DEAL::p=3.000 1.936 2.404, n=1.000 0.000 0.000
+DEAL::p=3.000 1.936 2.872, n=1.000 0.000 0.000
+DEAL::p=3.000 2.404 1.000, n=1.000 0.000 0.000
+DEAL::p=3.000 2.404 1.468, n=1.000 0.000 0.000
+DEAL::p=3.000 2.404 1.936, n=1.000 0.000 0.000
+DEAL::p=3.000 2.404 2.404, n=1.000 0.000 0.000
+DEAL::p=3.000 2.404 2.872, n=1.000 0.000 0.000
+DEAL::p=3.000 2.872 1.000, n=1.000 0.000 0.000
+DEAL::p=3.000 2.872 1.468, n=1.000 0.000 0.000
+DEAL::p=3.000 2.872 1.936, n=1.000 0.000 0.000
+DEAL::p=3.000 2.872 2.404, n=1.000 0.000 0.000
+DEAL::p=3.000 2.872 2.872, n=1.000 0.000 0.000
+DEAL::p=1.000 1.000 1.000, n=0.000 1.000 0.000
+DEAL::p=1.468 1.000 1.000, n=0.000 1.000 0.000
+DEAL::p=1.936 1.000 1.000, n=0.000 1.000 0.000
+DEAL::p=2.404 1.000 1.000, n=0.000 1.000 0.000
+DEAL::p=2.872 1.000 1.000, n=0.000 1.000 0.000
+DEAL::p=1.000 1.000 1.468, n=0.000 1.000 0.000
+DEAL::p=1.468 1.000 1.468, n=0.000 1.000 0.000
+DEAL::p=1.936 1.000 1.468, n=0.000 1.000 0.000
+DEAL::p=2.404 1.000 1.468, n=0.000 1.000 0.000
+DEAL::p=2.872 1.000 1.468, n=0.000 1.000 0.000
+DEAL::p=1.000 1.000 1.936, n=0.000 1.000 0.000
+DEAL::p=1.468 1.000 1.936, n=0.000 1.000 0.000
+DEAL::p=1.936 1.000 1.936, n=0.000 1.000 0.000
+DEAL::p=2.404 1.000 1.936, n=0.000 1.000 0.000
+DEAL::p=2.872 1.000 1.936, n=0.000 1.000 0.000
+DEAL::p=1.000 1.000 2.404, n=0.000 1.000 0.000
+DEAL::p=1.468 1.000 2.404, n=0.000 1.000 0.000
+DEAL::p=1.936 1.000 2.404, n=0.000 1.000 0.000
+DEAL::p=2.404 1.000 2.404, n=0.000 1.000 0.000
+DEAL::p=2.872 1.000 2.404, n=0.000 1.000 0.000
+DEAL::p=1.000 1.000 2.872, n=0.000 1.000 0.000
+DEAL::p=1.468 1.000 2.872, n=0.000 1.000 0.000
+DEAL::p=1.936 1.000 2.872, n=0.000 1.000 0.000
+DEAL::p=2.404 1.000 2.872, n=0.000 1.000 0.000
+DEAL::p=2.872 1.000 2.872, n=0.000 1.000 0.000
+DEAL::p=1.000 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=1.000 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=1.000 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=1.000 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=1.000 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=1.000 1.000 1.000, n=0.000 0.000 1.000
+DEAL::p=1.000 1.468 1.000, n=0.000 0.000 1.000
+DEAL::p=1.000 1.936 1.000, n=0.000 0.000 1.000
+DEAL::p=1.000 2.404 1.000, n=0.000 0.000 1.000
+DEAL::p=1.000 2.872 1.000, n=0.000 0.000 1.000
+DEAL::p=1.468 1.000 1.000, n=0.000 0.000 1.000
+DEAL::p=1.468 1.468 1.000, n=0.000 0.000 1.000
+DEAL::p=1.468 1.936 1.000, n=0.000 0.000 1.000
+DEAL::p=1.468 2.404 1.000, n=0.000 0.000 1.000
+DEAL::p=1.468 2.872 1.000, n=0.000 0.000 1.000
+DEAL::p=1.936 1.000 1.000, n=0.000 0.000 1.000
+DEAL::p=1.936 1.468 1.000, n=0.000 0.000 1.000
+DEAL::p=1.936 1.936 1.000, n=0.000 0.000 1.000
+DEAL::p=1.936 2.404 1.000, n=0.000 0.000 1.000
+DEAL::p=1.936 2.872 1.000, n=0.000 0.000 1.000
+DEAL::p=2.404 1.000 1.000, n=0.000 0.000 1.000
+DEAL::p=2.404 1.468 1.000, n=0.000 0.000 1.000
+DEAL::p=2.404 1.936 1.000, n=0.000 0.000 1.000
+DEAL::p=2.404 2.404 1.000, n=0.000 0.000 1.000
+DEAL::p=2.404 2.872 1.000, n=0.000 0.000 1.000
+DEAL::p=2.872 1.000 1.000, n=0.000 0.000 1.000
+DEAL::p=2.872 1.468 1.000, n=0.000 0.000 1.000
+DEAL::p=2.872 1.936 1.000, n=0.000 0.000 1.000
+DEAL::p=2.872 2.404 1.000, n=0.000 0.000 1.000
+DEAL::p=2.872 2.872 1.000, n=0.000 0.000 1.000
+DEAL::p=1.000 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=1.000 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=1.000 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=1.000 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=1.000 2.872 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 2.872 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 2.872 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 2.872 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 2.872 3.000, n=0.000 0.000 1.000
+DEAL::Case1
+DEAL::p=0.000 -0.500 -1.000, n=0.968 -0.176 -0.176
+DEAL::p=0.234 -0.149 -0.064, n=0.971 -0.145 -0.189
+DEAL::p=0.468 0.202 0.872, n=0.973 -0.108 -0.203
+DEAL::p=0.702 0.553 1.808, n=0.974 -0.065 -0.219
+DEAL::p=0.936 0.904 2.744, n=0.971 -0.015 -0.237
+DEAL::p=0.234 0.319 -0.532, n=0.970 -0.193 -0.148
+DEAL::p=0.413 0.588 0.294, n=0.974 -0.159 -0.159
+DEAL::p=0.592 0.857 1.121, n=0.978 -0.120 -0.173
+DEAL::p=0.772 1.126 1.947, n=0.979 -0.073 -0.189
+DEAL::p=0.951 1.394 2.774, n=0.978 -0.017 -0.207
+DEAL::p=0.468 1.138 -0.064, n=0.971 -0.213 -0.113
+DEAL::p=0.592 1.325 0.653, n=0.976 -0.178 -0.123
+DEAL::p=0.717 1.511 1.370, n=0.982 -0.135 -0.135
+DEAL::p=0.841 1.698 2.087, n=0.985 -0.084 -0.149
+DEAL::p=0.966 1.885 2.804, n=0.986 -0.020 -0.166
+DEAL::p=0.702 1.957 0.404, n=0.969 -0.237 -0.071
+DEAL::p=0.772 2.062 1.011, n=0.977 -0.200 -0.078
+DEAL::p=0.841 2.166 1.619, n=0.984 -0.154 -0.086
+DEAL::p=0.911 2.271 2.226, n=0.991 -0.097 -0.097
+DEAL::p=0.981 2.375 2.834, n=0.994 -0.024 -0.110
+DEAL::p=0.936 2.776 0.872, n=0.964 -0.266 -0.017
+DEAL::p=0.951 2.798 1.370, n=0.974 -0.228 -0.019
+DEAL::p=0.966 2.821 1.868, n=0.984 -0.179 -0.022
+DEAL::p=0.981 2.843 2.366, n=0.993 -0.115 -0.025
+DEAL::p=0.996 2.866 2.864, n=0.999 -0.029 -0.029
+DEAL::p=1.250 0.250 0.250, n=0.816 -0.408 -0.408
+DEAL::p=1.659 0.425 0.894, n=0.833 -0.336 -0.439
+DEAL::p=2.069 0.601 1.537, n=0.846 -0.250 -0.470
+DEAL::p=2.479 0.777 2.181, n=0.852 -0.149 -0.502
+DEAL::p=2.888 0.952 2.824, n=0.847 -0.034 -0.530
+DEAL::p=1.659 0.894 0.425, n=0.833 -0.439 -0.336
+DEAL::p=1.973 1.028 1.028, n=0.857 -0.365 -0.365
+DEAL::p=2.287 1.162 1.630, n=0.877 -0.274 -0.395
+DEAL::p=2.601 1.297 2.233, n=0.889 -0.166 -0.426
+DEAL::p=2.914 1.431 2.835, n=0.890 -0.038 -0.455
+DEAL::p=2.069 1.537 0.601, n=0.846 -0.470 -0.250
+DEAL::p=2.287 1.630 1.162, n=0.877 -0.395 -0.274
+DEAL::p=2.505 1.724 1.724, n=0.905 -0.301 -0.301
+DEAL::p=2.723 1.817 2.285, n=0.926 -0.184 -0.329
+DEAL::p=2.940 1.910 2.846, n=0.934 -0.043 -0.355
+DEAL::p=2.479 2.181 0.777, n=0.852 -0.502 -0.149
+DEAL::p=2.601 2.233 1.297, n=0.889 -0.426 -0.166
+DEAL::p=2.723 2.285 1.817, n=0.926 -0.329 -0.184
+DEAL::p=2.845 2.337 2.337, n=0.957 -0.204 -0.204
+DEAL::p=2.967 2.390 2.858, n=0.974 -0.048 -0.224
+DEAL::p=2.888 2.824 0.952, n=0.847 -0.530 -0.034
+DEAL::p=2.914 2.835 1.431, n=0.890 -0.455 -0.038
+DEAL::p=2.940 2.846 1.910, n=0.934 -0.355 -0.043
+DEAL::p=2.967 2.858 2.390, n=0.974 -0.224 -0.048
+DEAL::p=2.993 2.869 2.869, n=0.997 -0.053 -0.053
+DEAL::p=0.000 -0.500 -1.000, n=-0.276 0.921 -0.276
+DEAL::p=0.293 -0.325 -0.708, n=-0.326 0.918 -0.225
+DEAL::p=0.585 -0.149 -0.415, n=-0.395 0.906 -0.149
+DEAL::p=0.878 0.027 -0.122, n=-0.491 0.871 -0.031
+DEAL::p=1.170 0.202 0.170, n=-0.616 0.772 0.153
+DEAL::p=0.234 -0.149 -0.064, n=-0.170 0.936 -0.308
+DEAL::p=0.568 -0.015 0.160, n=-0.195 0.942 -0.275
+DEAL::p=0.901 0.120 0.384, n=-0.228 0.947 -0.228
+DEAL::p=1.235 0.254 0.608, n=-0.273 0.948 -0.162
+DEAL::p=1.568 0.389 0.832, n=-0.337 0.940 -0.062
+DEAL::p=0.468 0.202 0.872, n=-0.098 0.940 -0.328
+DEAL::p=0.843 0.295 1.028, n=-0.110 0.947 -0.303
+DEAL::p=1.217 0.389 1.183, n=-0.125 0.954 -0.272
+DEAL::p=1.592 0.482 1.339, n=-0.145 0.963 -0.229
+DEAL::p=1.967 0.575 1.494, n=-0.171 0.970 -0.171
+DEAL::p=0.702 0.553 1.808, n=-0.047 0.939 -0.340
+DEAL::p=1.118 0.605 1.895, n=-0.052 0.946 -0.321
+DEAL::p=1.533 0.658 1.982, n=-0.057 0.953 -0.298
+DEAL::p=1.949 0.710 2.069, n=-0.065 0.961 -0.268
+DEAL::p=2.365 0.762 2.157, n=-0.074 0.971 -0.229
+DEAL::p=0.936 0.904 2.744, n=0.009 0.937 -0.349
+DEAL::p=1.393 0.915 2.763, n=-0.010 0.943 -0.334
+DEAL::p=1.850 0.926 2.781, n=-0.010 0.949 -0.315
+DEAL::p=2.306 0.938 2.800, n=-0.012 0.956 -0.292
+DEAL::p=2.763 0.949 2.819, n=-0.013 0.964 -0.264
+DEAL::p=1.000 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 1.000, n=0.000 1.000 0.000
+DEAL::p=1.000 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 1.468, n=0.000 1.000 0.000
+DEAL::p=1.000 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 1.936, n=0.000 1.000 0.000
+DEAL::p=1.000 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 2.404, n=0.000 1.000 0.000
+DEAL::p=1.000 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=1.468 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=1.936 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=2.404 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=2.872 3.000 2.872, n=0.000 1.000 0.000
+DEAL::p=0.000 -0.500 -1.000, n=-0.600 -0.261 0.756
+DEAL::p=0.234 0.319 -0.532, n=-0.417 -0.358 0.835
+DEAL::p=0.468 1.138 -0.064, n=-0.255 -0.424 0.869
+DEAL::p=0.702 1.957 0.404, n=-0.125 -0.465 0.876
+DEAL::p=0.936 2.776 0.872, n=-0.024 -0.491 0.871
+DEAL::p=0.293 -0.325 -0.708, n=-0.652 -0.151 0.743
+DEAL::p=0.568 0.453 -0.308, n=-0.458 -0.272 0.846
+DEAL::p=0.843 1.231 0.092, n=-0.281 -0.358 0.890
+DEAL::p=1.118 2.009 0.491, n=-0.137 -0.414 0.900
+DEAL::p=1.393 2.787 0.891, n=-0.026 -0.450 0.893
+DEAL::p=0.585 -0.149 -0.415, n=-0.701 -0.020 0.713
+DEAL::p=0.901 0.588 -0.084, n=-0.503 -0.165 0.848
+DEAL::p=1.217 1.325 0.247, n=-0.309 -0.276 0.910
+DEAL::p=1.533 2.062 0.578, n=-0.150 -0.351 0.924
+DEAL::p=1.850 2.798 0.909, n=-0.028 -0.400 0.916
+DEAL::p=0.878 0.027 -0.122, n=-0.739 0.130 0.661
+DEAL::p=1.235 0.722 0.140, n=-0.548 -0.034 0.836
+DEAL::p=1.592 1.418 0.403, n=-0.340 -0.174 0.924
+DEAL::p=1.949 2.114 0.665, n=-0.164 -0.273 0.948
+DEAL::p=2.306 2.810 0.928, n=-0.030 -0.339 0.940
+DEAL::p=1.170 0.202 0.170, n=-0.758 0.287 0.586
+DEAL::p=1.568 0.857 0.364, n=-0.586 0.119 0.801
+DEAL::p=1.967 1.511 0.558, n=-0.373 -0.048 0.927
+DEAL::p=2.365 2.166 0.753, n=-0.181 -0.177 0.967
+DEAL::p=2.763 2.821 0.947, n=-0.033 -0.266 0.963
+DEAL::p=1.000 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=1.000 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=1.000 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=1.000 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=1.000 2.872 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=1.468 2.872 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=1.936 2.872 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=2.404 2.872 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 1.000 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 1.468 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 1.936 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 2.404 3.000, n=0.000 0.000 1.000
+DEAL::p=2.872 2.872 3.000, n=0.000 0.000 1.000
+DEAL::OK
diff --git a/tests/deal.II/normal_vector_03_2d.cc b/tests/deal.II/normal_vector_03_2d.cc
new file mode 100644 (file)
index 0000000..e79dd1c
--- /dev/null
@@ -0,0 +1,114 @@
+//----------------------------  normal_vector_03_2d.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2005, 2008, 2011 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.
+//
+//----------------------------  normal_vector_03_2d.cc  ---------------------------
+
+
+// like _01 but instead for random points on the face. note that the normal
+// vector on a bilinearly mapped face is also a bilinear function
+
+
+
+#include "../tests.h"
+#include <deal.II/base/quadrature_lib.h>
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_boundary.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/fe/fe_q.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+
+void create_triangulation(const unsigned int case_no,
+                         Triangulation<2> &tria)
+{
+  switch (case_no)
+    {
+      case 0:
+           GridGenerator::hyper_cube(tria, 1., 3.);
+           break;
+      case 1:
+      {
+       GridGenerator::hyper_cube(tria, 1., 3.);
+       Point<2> &v0=tria.begin_active()->vertex(0);
+       v0 = Point<2> (-0.5,-1);
+       Point<2> &v1=tria.begin_active()->vertex(1);
+       v1 = Point<2> (1.25, 0.25);
+       break;
+      }
+      default:
+           Assert(false, ExcNotImplemented());
+    };
+}
+
+
+
+int main ()
+{
+  std::ofstream logfile ("normal_vector_03_2d/output");
+  deallog << std::setprecision (3);
+  deallog << std::fixed;
+  deallog.attach(logfile);
+  deallog.depth_console (0);
+
+  FE_Q<2> linear_interpolator(1);
+
+  Triangulation<2> tria;
+  StraightBoundary<2> boundary;
+  Boundary<2>::FaceVertexNormals normals;
+  for (unsigned int case_no=0; case_no<2; ++case_no)
+    {
+      deallog << "Case" << case_no << std::endl;
+      create_triangulation(case_no, tria);
+      const Triangulation<2>::active_cell_iterator cell=tria.begin_active();
+      Triangulation<2>::face_iterator face;
+      for (unsigned int face_no=0; face_no<GeometryInfo<2>::faces_per_cell; ++face_no)
+       {
+         face=cell->face(face_no);
+         boundary.get_normals_at_vertices(face, normals);
+
+         for (double xi=0; xi<=1; xi+=0.234)
+           for (double eta=0; eta<=1; eta+=0.234)
+             {
+               Point<2> p;
+               Tensor<1,2> normal;
+
+               for (unsigned int v=0; v<GeometryInfo<2>::vertices_per_face; ++v)
+                 {
+                   p += face->vertex(v) * linear_interpolator.shape_value(v,Point<2>(xi,eta));
+                   normal += normals[v] *
+                             linear_interpolator.shape_value(v,Point<2>(xi,eta));
+                 }
+               normal /= normal.norm();
+
+               deallog << "p=" << p
+                       << ", n=" << boundary.normal_vector (face, p)
+                       << std::endl;
+
+               Assert ((boundary.normal_vector (face, p)
+                        -
+                       normal).norm()
+                       <
+                       1e-10,
+                       ExcInternalError());
+             }
+       }
+
+      tria.clear();
+    }
+
+  deallog << "OK" << std::endl;
+}
diff --git a/tests/deal.II/normal_vector_03_2d/cmp/generic b/tests/deal.II/normal_vector_03_2d/cmp/generic
new file mode 100644 (file)
index 0000000..209f98d
--- /dev/null
@@ -0,0 +1,204 @@
+
+DEAL::Case0
+DEAL::p=1.000 1.000, n=1.000 0.000
+DEAL::p=0.766 0.766, n=1.000 0.000
+DEAL::p=0.532 0.532, n=1.000 0.000
+DEAL::p=0.298 0.298, n=1.000 0.000
+DEAL::p=0.064 0.064, n=1.000 0.000
+DEAL::p=1.000 1.468, n=1.000 0.000
+DEAL::p=0.766 1.124, n=1.000 0.000
+DEAL::p=0.532 0.781, n=1.000 0.000
+DEAL::p=0.298 0.437, n=1.000 0.000
+DEAL::p=0.064 0.094, n=1.000 0.000
+DEAL::p=1.000 1.936, n=1.000 0.000
+DEAL::p=0.766 1.483, n=1.000 0.000
+DEAL::p=0.532 1.030, n=1.000 0.000
+DEAL::p=0.298 0.577, n=1.000 0.000
+DEAL::p=0.064 0.124, n=1.000 0.000
+DEAL::p=1.000 2.404, n=1.000 0.000
+DEAL::p=0.766 1.841, n=1.000 0.000
+DEAL::p=0.532 1.279, n=1.000 0.000
+DEAL::p=0.298 0.716, n=1.000 0.000
+DEAL::p=0.064 0.154, n=1.000 0.000
+DEAL::p=1.000 2.872, n=1.000 0.000
+DEAL::p=0.766 2.200, n=1.000 0.000
+DEAL::p=0.532 1.528, n=1.000 0.000
+DEAL::p=0.298 0.856, n=1.000 0.000
+DEAL::p=0.064 0.184, n=1.000 0.000
+DEAL::p=3.000 1.000, n=1.000 0.000
+DEAL::p=2.298 0.766, n=1.000 0.000
+DEAL::p=1.596 0.532, n=1.000 0.000
+DEAL::p=0.894 0.298, n=1.000 0.000
+DEAL::p=0.192 0.064, n=1.000 0.000
+DEAL::p=3.000 1.468, n=1.000 0.000
+DEAL::p=2.298 1.124, n=1.000 0.000
+DEAL::p=1.596 0.781, n=1.000 0.000
+DEAL::p=0.894 0.437, n=1.000 0.000
+DEAL::p=0.192 0.094, n=1.000 0.000
+DEAL::p=3.000 1.936, n=1.000 0.000
+DEAL::p=2.298 1.483, n=1.000 0.000
+DEAL::p=1.596 1.030, n=1.000 0.000
+DEAL::p=0.894 0.577, n=1.000 0.000
+DEAL::p=0.192 0.124, n=1.000 0.000
+DEAL::p=3.000 2.404, n=1.000 0.000
+DEAL::p=2.298 1.841, n=1.000 0.000
+DEAL::p=1.596 1.279, n=1.000 0.000
+DEAL::p=0.894 0.716, n=1.000 0.000
+DEAL::p=0.192 0.154, n=1.000 0.000
+DEAL::p=3.000 2.872, n=1.000 0.000
+DEAL::p=2.298 2.200, n=1.000 0.000
+DEAL::p=1.596 1.528, n=1.000 0.000
+DEAL::p=0.894 0.856, n=1.000 0.000
+DEAL::p=0.192 0.184, n=1.000 0.000
+DEAL::p=1.000 1.000, n=0.000 -1.000
+DEAL::p=0.766 0.766, n=0.000 -1.000
+DEAL::p=0.532 0.532, n=0.000 -1.000
+DEAL::p=0.298 0.298, n=0.000 -1.000
+DEAL::p=0.064 0.064, n=0.000 -1.000
+DEAL::p=1.468 1.000, n=0.000 -1.000
+DEAL::p=1.124 0.766, n=0.000 -1.000
+DEAL::p=0.781 0.532, n=0.000 -1.000
+DEAL::p=0.437 0.298, n=0.000 -1.000
+DEAL::p=0.094 0.064, n=0.000 -1.000
+DEAL::p=1.936 1.000, n=0.000 -1.000
+DEAL::p=1.483 0.766, n=0.000 -1.000
+DEAL::p=1.030 0.532, n=0.000 -1.000
+DEAL::p=0.577 0.298, n=0.000 -1.000
+DEAL::p=0.124 0.064, n=0.000 -1.000
+DEAL::p=2.404 1.000, n=0.000 -1.000
+DEAL::p=1.841 0.766, n=0.000 -1.000
+DEAL::p=1.279 0.532, n=0.000 -1.000
+DEAL::p=0.716 0.298, n=0.000 -1.000
+DEAL::p=0.154 0.064, n=0.000 -1.000
+DEAL::p=2.872 1.000, n=0.000 -1.000
+DEAL::p=2.200 0.766, n=0.000 -1.000
+DEAL::p=1.528 0.532, n=0.000 -1.000
+DEAL::p=0.856 0.298, n=0.000 -1.000
+DEAL::p=0.184 0.064, n=0.000 -1.000
+DEAL::p=1.000 3.000, n=0.000 -1.000
+DEAL::p=0.766 2.298, n=0.000 -1.000
+DEAL::p=0.532 1.596, n=0.000 -1.000
+DEAL::p=0.298 0.894, n=0.000 -1.000
+DEAL::p=0.064 0.192, n=0.000 -1.000
+DEAL::p=1.468 3.000, n=0.000 -1.000
+DEAL::p=1.124 2.298, n=0.000 -1.000
+DEAL::p=0.781 1.596, n=0.000 -1.000
+DEAL::p=0.437 0.894, n=0.000 -1.000
+DEAL::p=0.094 0.192, n=0.000 -1.000
+DEAL::p=1.936 3.000, n=0.000 -1.000
+DEAL::p=1.483 2.298, n=0.000 -1.000
+DEAL::p=1.030 1.596, n=0.000 -1.000
+DEAL::p=0.577 0.894, n=0.000 -1.000
+DEAL::p=0.124 0.192, n=0.000 -1.000
+DEAL::p=2.404 3.000, n=0.000 -1.000
+DEAL::p=1.841 2.298, n=0.000 -1.000
+DEAL::p=1.279 1.596, n=0.000 -1.000
+DEAL::p=0.716 0.894, n=0.000 -1.000
+DEAL::p=0.154 0.192, n=0.000 -1.000
+DEAL::p=2.872 3.000, n=0.000 -1.000
+DEAL::p=2.200 2.298, n=0.000 -1.000
+DEAL::p=1.528 1.596, n=0.000 -1.000
+DEAL::p=0.856 0.894, n=0.000 -1.000
+DEAL::p=0.184 0.192, n=0.000 -1.000
+DEAL::Case1
+DEAL::p=-0.500 -1.000, n=0.936 -0.351
+DEAL::p=-0.383 -0.766, n=0.936 -0.351
+DEAL::p=-0.266 -0.532, n=0.936 -0.351
+DEAL::p=-0.149 -0.298, n=0.936 -0.351
+DEAL::p=-0.032 -0.064, n=0.936 -0.351
+DEAL::p=-0.149 -0.064, n=0.936 -0.351
+DEAL::p=-0.114 -0.049, n=0.936 -0.351
+DEAL::p=-0.079 -0.034, n=0.936 -0.351
+DEAL::p=-0.044 -0.019, n=0.936 -0.351
+DEAL::p=-0.010 0.004, n=0.936 -0.351
+DEAL::p=0.202 0.872, n=0.936 -0.351
+DEAL::p=0.155 0.668, n=0.936 -0.351
+DEAL::p=0.107 0.464, n=0.936 -0.351
+DEAL::p=0.060 0.260, n=0.936 -0.351
+DEAL::p=0.013 0.056, n=0.936 -0.351
+DEAL::p=0.553 1.808, n=0.936 -0.351
+DEAL::p=0.424 1.385, n=0.936 -0.351
+DEAL::p=0.294 0.962, n=0.936 -0.351
+DEAL::p=0.165 0.539, n=0.936 -0.351
+DEAL::p=0.035 0.116, n=0.936 -0.351
+DEAL::p=0.904 2.744, n=0.936 -0.351
+DEAL::p=0.692 2.102, n=0.936 -0.351
+DEAL::p=0.481 1.460, n=0.936 -0.351
+DEAL::p=0.269 0.818, n=0.936 -0.351
+DEAL::p=0.058 0.176, n=0.936 -0.351
+DEAL::p=1.250 0.250, n=0.844 -0.537
+DEAL::p=0.958 0.192, n=0.844 -0.537
+DEAL::p=0.665 0.133, n=0.844 -0.537
+DEAL::p=0.372 0.074, n=0.844 -0.537
+DEAL::p=0.080 0.016, n=0.844 -0.537
+DEAL::p=1.659 0.894, n=0.844 -0.537
+DEAL::p=1.271 0.684, n=0.844 -0.537
+DEAL::p=0.883 0.475, n=0.844 -0.537
+DEAL::p=0.495 0.266, n=0.844 -0.537
+DEAL::p=0.106 0.057, n=0.844 -0.537
+DEAL::p=2.069 1.537, n=0.844 -0.537
+DEAL::p=1.585 1.177, n=0.844 -0.537
+DEAL::p=1.101 0.818, n=0.844 -0.537
+DEAL::p=0.617 0.458, n=0.844 -0.537
+DEAL::p=0.132 0.098, n=0.844 -0.537
+DEAL::p=2.479 2.181, n=0.844 -0.537
+DEAL::p=1.899 1.670, n=0.844 -0.537
+DEAL::p=1.319 1.160, n=0.844 -0.537
+DEAL::p=0.739 0.650, n=0.844 -0.537
+DEAL::p=0.159 0.140, n=0.844 -0.537
+DEAL::p=2.888 2.824, n=0.844 -0.537
+DEAL::p=2.212 2.163, n=0.844 -0.537
+DEAL::p=1.536 1.502, n=0.844 -0.537
+DEAL::p=0.861 0.842, n=0.844 -0.537
+DEAL::p=0.185 0.181, n=0.844 -0.537
+DEAL::p=-0.500 -1.000, n=0.581 -0.814
+DEAL::p=-0.383 -0.766, n=0.581 -0.814
+DEAL::p=-0.266 -0.532, n=0.581 -0.814
+DEAL::p=-0.149 -0.298, n=0.581 -0.814
+DEAL::p=-0.032 -0.064, n=0.581 -0.814
+DEAL::p=-0.090 -0.708, n=0.581 -0.814
+DEAL::p=-0.069 -0.542, n=0.581 -0.814
+DEAL::p=-0.048 -0.376, n=0.581 -0.814
+DEAL::p=-0.027 -0.211, n=0.581 -0.814
+DEAL::p=0.006 -0.045, n=0.581 -0.814
+DEAL::p=0.319 -0.415, n=0.581 -0.814
+DEAL::p=0.244 -0.318, n=0.581 -0.814
+DEAL::p=0.170 -0.221, n=0.581 -0.814
+DEAL::p=0.095 -0.124, n=0.581 -0.814
+DEAL::p=0.020 -0.027, n=0.581 -0.814
+DEAL::p=0.729 -0.122, n=0.581 -0.814
+DEAL::p=0.558 -0.094, n=0.581 -0.814
+DEAL::p=0.388 -0.065, n=0.581 -0.814
+DEAL::p=0.217 -0.037, n=0.581 -0.814
+DEAL::p=0.047 0.008, n=0.581 -0.814
+DEAL::p=1.138 0.170, n=0.581 -0.814
+DEAL::p=0.872 0.130, n=0.581 -0.814
+DEAL::p=0.605 0.090, n=0.581 -0.814
+DEAL::p=0.339 0.051, n=0.581 -0.814
+DEAL::p=0.073 0.011, n=0.581 -0.814
+DEAL::p=1.000 3.000, n=0.000 -1.000
+DEAL::p=0.766 2.298, n=0.000 -1.000
+DEAL::p=0.532 1.596, n=0.000 -1.000
+DEAL::p=0.298 0.894, n=0.000 -1.000
+DEAL::p=0.064 0.192, n=0.000 -1.000
+DEAL::p=1.468 3.000, n=0.000 -1.000
+DEAL::p=1.124 2.298, n=0.000 -1.000
+DEAL::p=0.781 1.596, n=0.000 -1.000
+DEAL::p=0.437 0.894, n=0.000 -1.000
+DEAL::p=0.094 0.192, n=0.000 -1.000
+DEAL::p=1.936 3.000, n=0.000 -1.000
+DEAL::p=1.483 2.298, n=0.000 -1.000
+DEAL::p=1.030 1.596, n=0.000 -1.000
+DEAL::p=0.577 0.894, n=0.000 -1.000
+DEAL::p=0.124 0.192, n=0.000 -1.000
+DEAL::p=2.404 3.000, n=0.000 -1.000
+DEAL::p=1.841 2.298, n=0.000 -1.000
+DEAL::p=1.279 1.596, n=0.000 -1.000
+DEAL::p=0.716 0.894, n=0.000 -1.000
+DEAL::p=0.154 0.192, n=0.000 -1.000
+DEAL::p=2.872 3.000, n=0.000 -1.000
+DEAL::p=2.200 2.298, n=0.000 -1.000
+DEAL::p=1.528 1.596, n=0.000 -1.000
+DEAL::p=0.856 0.894, n=0.000 -1.000
+DEAL::p=0.184 0.192, n=0.000 -1.000
+DEAL::OK

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.