]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added test.
authorLuca Heltai <luca.heltai@sissa.it>
Sun, 20 Mar 2016 18:11:51 +0000 (19:11 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Sun, 20 Mar 2016 18:14:34 +0000 (19:14 +0100)
source/grid/manifold.cc
tests/manifold/spherical_manifold_07.cc [new file with mode: 0644]
tests/manifold/spherical_manifold_07.output [new file with mode: 0644]

index 52cece6198804b601cf354e53e71b5b57db66673..e3fb14c7f7515b905f002e5ac88b8fd8aa812138 100644 (file)
@@ -73,7 +73,7 @@ normal_vector (const typename Triangulation<2, 2>::face_iterator &face,
 
   Tensor<1,spacedim> tangent = ((p-face->vertex(0)).norm_square() > (p-face->vertex(1)).norm_square() ?
                                 get_tangent_vector(p, face->vertex(0)) :
-                                get_tangent_vector(p, face->vertex(1)));
+                                -get_tangent_vector(p, face->vertex(1)));
   return cross_product_2d(tangent);
 }
 
@@ -175,7 +175,10 @@ get_normals_at_vertices (const typename Triangulation<dim, spacedim>::face_itera
                          FaceVertexNormals &n) const
 {
   for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_face; ++v)
-    n[v] = normal_vector(face, face->vertex(v));
+    {
+      n[v] = normal_vector(face, face->vertex(v));
+      n[v] /= n[v].norm();
+    }
 }
 
 
diff --git a/tests/manifold/spherical_manifold_07.cc b/tests/manifold/spherical_manifold_07.cc
new file mode 100644 (file)
index 0000000..f92a567
--- /dev/null
@@ -0,0 +1,82 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2016 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// test GridTools::torus() and TorusManifold, output visually checked
+
+#include "../tests.h"
+#include <deal.II/base/logstream.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/manifold.h>
+#include <deal.II/grid/manifold_lib.h>
+
+#include <fstream>
+#include <iomanip>
+
+
+template <int dim, int spacedim>
+void test ()
+{
+  Triangulation<dim, spacedim> triangulation;
+
+  GridGenerator::hyper_ball(triangulation);
+
+  static const SphericalManifold<dim, spacedim> manifold;
+  triangulation.set_all_manifold_ids_on_boundary(0);
+  triangulation.set_manifold (0, manifold);
+
+  triangulation.refine_global(1);
+
+  unsigned int c = 0;
+  deallog.get_file_stream()
+      << "set view equal xyz" << std::endl
+      << "set size ratio -1" << std::endl
+      << "set multiplot" << std::endl
+      << "splot '-' with vectors " << std::endl;
+
+  for (typename Triangulation<dim,spacedim>::active_cell_iterator
+       cell = triangulation.begin_active();
+       cell != triangulation.end(); ++cell)
+    {
+      typename Manifold<dim,spacedim>::FaceVertexNormals n;
+
+      for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+        if (cell->face(f)->at_boundary())
+          {
+            manifold.get_normals_at_vertices(cell->face(f), n);
+
+            for (unsigned int v=0; v<GeometryInfo<dim>::vertices_per_face; ++v)
+              {
+
+                deallog.get_file_stream() << cell->face(f)->vertex(v)
+                                          << " " << n[v]*.1 << std::endl;
+              }
+          }
+    }
+  deallog.get_file_stream() << "e" << std::endl;
+}
+
+
+int main ()
+{
+  initlog ();
+
+  test<2,2> ();
+  test<3,3> ();
+
+  return 0;
+}
diff --git a/tests/manifold/spherical_manifold_07.output b/tests/manifold/spherical_manifold_07.output
new file mode 100644 (file)
index 0000000..9bd8aaa
--- /dev/null
@@ -0,0 +1,123 @@
+
+set view equal xyz
+set size ratio -1
+set multiplot
+splot '-' with vectors 
+-0.707107 -0.707107 0.0707107 0.0707107
+-1.83697e-16 -1.00000 1.83697e-17 0.100000
+-1.83697e-16 -1.00000 1.83697e-17 0.100000
+0.707107 -0.707107 -0.0707107 0.0707107
+-0.707107 -0.707107 -0.0707107 -0.0707107
+-1.00000 1.22465e-16 -0.100000 1.22465e-17
+-1.00000 1.22465e-16 -0.100000 1.22465e-17
+-0.707107 0.707107 -0.0707107 0.0707107
+0.707107 -0.707107 -0.0707107 0.0707107
+1.00000 0.00000 -0.100000 0.00000
+1.00000 0.00000 -0.100000 0.00000
+0.707107 0.707107 -0.0707107 -0.0707107
+-0.707107 0.707107 -0.0707107 0.0707107
+6.12323e-17 1.00000 6.12323e-18 0.100000
+6.12323e-17 1.00000 6.12323e-18 0.100000
+0.707107 0.707107 0.0707107 0.0707107
+e
+set view equal xyz
+set size ratio -1
+set multiplot
+splot '-' with vectors 
+-0.577350 -0.577350 -0.577350 0.0577350 0.0577350 0.0577350
+0.00000 -0.707107 -0.707107 1.68252e-17 0.0707107 0.0707107
+-0.707107 0.00000 -0.707107 0.0707107 0.00000 0.0707107
+0.00000 0.00000 -1.00000 1.22465e-17 0.00000 -0.100000
+0.00000 -0.707107 -0.707107 1.11038e-17 0.0707107 0.0707107
+0.577350 -0.577350 -0.577350 -0.0577350 0.0577350 0.0577350
+0.00000 0.00000 -1.00000 -1.22465e-17 0.00000 0.100000
+0.707107 0.00000 -0.707107 -0.0707107 3.18167e-18 0.0707107
+-0.707107 0.00000 -0.707107 0.0707107 -8.52818e-18 0.0707107
+0.00000 0.00000 -1.00000 -1.22465e-17 0.00000 0.100000
+-0.577350 0.577350 -0.577350 0.0577350 -0.0577350 0.0577350
+0.00000 0.707107 -0.707107 5.60839e-18 -0.0707107 0.0707107
+0.00000 0.00000 -1.00000 -1.22465e-17 0.00000 0.100000
+0.707107 0.00000 -0.707107 -0.0707107 -3.18167e-18 0.0707107
+0.00000 0.707107 -0.707107 -6.66229e-18 -0.0707107 0.0707107
+0.577350 0.577350 -0.577350 -0.0577350 -0.0577350 0.0577350
+0.577350 -0.577350 -0.577350 -0.0577350 0.0577350 0.0577350
+0.707107 -0.707107 0.00000 -0.0707107 0.0707107 -6.12323e-18
+0.707107 0.00000 -0.707107 -0.0707107 3.18167e-18 0.0707107
+1.00000 0.00000 0.00000 -0.100000 0.00000 -6.12323e-18
+0.707107 0.00000 -0.707107 -0.0707107 -3.18167e-18 0.0707107
+1.00000 0.00000 0.00000 -0.100000 0.00000 -6.12323e-18
+0.577350 0.577350 -0.577350 -0.0577350 -0.0577350 0.0577350
+0.707107 0.707107 0.00000 -0.0707107 -0.0707107 -6.12323e-18
+0.707107 -0.707107 0.00000 -0.0707107 0.0707107 -6.12323e-18
+0.577350 -0.577350 0.577350 -0.0577350 0.0577350 -0.0577350
+1.00000 0.00000 0.00000 -0.100000 0.00000 -6.12323e-18
+0.707107 0.00000 0.707107 -0.0707107 -3.18167e-18 -0.0707107
+1.00000 0.00000 0.00000 -0.100000 0.00000 -6.12323e-18
+0.707107 0.00000 0.707107 -0.0707107 3.18167e-18 -0.0707107
+0.707107 0.707107 0.00000 -0.0707107 -0.0707107 -6.12323e-18
+0.577350 0.577350 0.577350 -0.0577350 -0.0577350 -0.0577350
+-0.577350 -0.577350 0.577350 0.0577350 0.0577350 -0.0577350
+-0.707107 0.00000 0.707107 0.0707107 -4.72716e-17 -0.0707107
+0.00000 -0.707107 0.707107 5.60839e-18 0.0707107 -0.0707107
+0.00000 0.00000 1.00000 nan nan nan
+0.00000 -0.707107 0.707107 1.55453e-17 0.0707107 -0.0707107
+0.00000 0.00000 1.00000 nan nan nan
+0.577350 -0.577350 0.577350 -0.0577350 0.0577350 -0.0577350
+0.707107 0.00000 0.707107 -0.0707107 -3.18167e-18 -0.0707107
+-0.707107 0.00000 0.707107 0.0707107 -6.82254e-18 -0.0707107
+-0.577350 0.577350 0.577350 0.0577350 -0.0577350 -0.0577350
+0.00000 0.00000 1.00000 nan nan nan
+0.00000 0.707107 0.707107 -5.60839e-18 -0.0707107 -0.0707107
+0.00000 0.00000 1.00000 nan nan nan
+0.00000 0.707107 0.707107 -2.22076e-18 -0.0707107 -0.0707107
+0.707107 0.00000 0.707107 -0.0707107 3.18167e-18 -0.0707107
+0.577350 0.577350 0.577350 -0.0577350 -0.0577350 -0.0577350
+-0.577350 -0.577350 -0.577350 0.0577350 0.0577350 0.0577350
+-0.707107 0.00000 -0.707107 0.0707107 -6.36335e-18 0.0707107
+-0.707107 -0.707107 0.00000 0.0707107 0.0707107 -6.12323e-18
+-1.00000 0.00000 0.00000 0.100000 -1.22465e-17 -6.12323e-18
+-0.707107 0.00000 -0.707107 0.0707107 -9.54502e-18 0.0707107
+-0.577350 0.577350 -0.577350 0.0577350 -0.0577350 0.0577350
+-1.00000 0.00000 0.00000 0.100000 -1.22465e-17 -6.12323e-18
+-0.707107 0.707107 0.00000 0.0707107 -0.0707107 -6.12323e-18
+-0.707107 -0.707107 0.00000 0.0707107 0.0707107 -6.12323e-18
+-1.00000 0.00000 0.00000 0.100000 -1.22465e-17 -6.12323e-18
+-0.577350 -0.577350 0.577350 0.0577350 0.0577350 -0.0577350
+-0.707107 0.00000 0.707107 0.0707107 -1.27267e-17 -0.0707107
+-1.00000 0.00000 0.00000 0.100000 -1.22465e-17 -6.12323e-18
+-0.707107 0.707107 0.00000 0.0707107 -0.0707107 -6.12323e-18
+-0.707107 0.00000 0.707107 0.0707107 -6.36335e-18 -0.0707107
+-0.577350 0.577350 0.577350 0.0577350 -0.0577350 -0.0577350
+-0.577350 -0.577350 -0.577350 0.0577350 0.0577350 0.0577350
+-0.707107 -0.707107 0.00000 0.0707107 0.0707107 -6.12323e-18
+0.00000 -0.707107 -0.707107 1.59084e-17 0.0707107 0.0707107
+0.00000 -1.00000 0.00000 1.83697e-17 0.100000 -6.12323e-18
+0.00000 -0.707107 -0.707107 9.54502e-18 0.0707107 0.0707107
+0.00000 -1.00000 0.00000 1.83697e-17 0.100000 -6.12323e-18
+0.577350 -0.577350 -0.577350 -0.0577350 0.0577350 0.0577350
+0.707107 -0.707107 0.00000 -0.0707107 0.0707107 -6.12323e-18
+-0.707107 -0.707107 0.00000 0.0707107 0.0707107 -6.12323e-18
+-0.577350 -0.577350 0.577350 0.0577350 0.0577350 -0.0577350
+0.00000 -1.00000 0.00000 1.83697e-17 0.100000 -6.12323e-18
+0.00000 -0.707107 0.707107 9.54502e-18 0.0707107 -0.0707107
+0.00000 -1.00000 0.00000 1.83697e-17 0.100000 -6.12323e-18
+0.00000 -0.707107 0.707107 1.59084e-17 0.0707107 -0.0707107
+0.707107 -0.707107 0.00000 -0.0707107 0.0707107 -6.12323e-18
+0.577350 -0.577350 0.577350 -0.0577350 0.0577350 -0.0577350
+-0.577350 0.577350 -0.577350 0.0577350 -0.0577350 0.0577350
+0.00000 0.707107 -0.707107 0.00000 -0.0707107 0.0707107
+-0.707107 0.707107 0.00000 0.0707107 -0.0707107 -6.12323e-18
+0.00000 1.00000 0.00000 -6.12323e-18 -0.100000 -6.12323e-18
+0.00000 0.707107 -0.707107 -6.36335e-18 -0.0707107 0.0707107
+0.577350 0.577350 -0.577350 -0.0577350 -0.0577350 0.0577350
+0.00000 1.00000 0.00000 -6.12323e-18 -0.100000 -6.12323e-18
+0.707107 0.707107 0.00000 -0.0707107 -0.0707107 -6.12323e-18
+-0.707107 0.707107 0.00000 0.0707107 -0.0707107 -6.12323e-18
+0.00000 1.00000 0.00000 -6.12323e-18 -0.100000 -6.12323e-18
+-0.577350 0.577350 0.577350 0.0577350 -0.0577350 -0.0577350
+0.00000 0.707107 0.707107 -6.36335e-18 -0.0707107 -0.0707107
+0.00000 1.00000 0.00000 -6.12323e-18 -0.100000 -6.12323e-18
+0.707107 0.707107 0.00000 -0.0707107 -0.0707107 -6.12323e-18
+0.00000 0.707107 0.707107 0.00000 -0.0707107 -0.0707107
+0.577350 0.577350 0.577350 -0.0577350 -0.0577350 -0.0577350
+e

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.