]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added test for periodicity in chart manifold.
authorLuca Heltai <luca.heltai@sissa.it>
Thu, 17 Mar 2016 11:29:08 +0000 (12:29 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Sun, 27 Mar 2016 21:44:35 +0000 (23:44 +0200)
tests/manifold/chart_manifold_09.cc [new file with mode: 0644]
tests/manifold/chart_manifold_09.output [new file with mode: 0644]

diff --git a/tests/manifold/chart_manifold_09.cc b/tests/manifold/chart_manifold_09.cc
new file mode 100644 (file)
index 0000000..66fdebd
--- /dev/null
@@ -0,0 +1,71 @@
+// ---------------------------------------------------------------------
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+// Check ChartManifold for periodicity issues. Check that the it finds
+// the right intermediate points independent on the number of
+// surrounding points
+
+#include "../tests.h"
+
+#include <deal.II/base/utilities.h>
+#include <deal.II/grid/manifold_lib.h>
+
+
+int
+main()
+{
+  initlog();
+
+  const FunctionManifold<2,2,2> manifold("x;y", "x;y", Point<2>(1.0, 0.0));
+
+  // Some points that would cross the periodicity boundary
+  std::vector<Point<2> > points;
+
+  points.push_back(Point<2>(0.1, 0.0));
+  points.push_back(Point<2>(0.9, 0.0));
+
+  // And the weights
+  std::vector<double> weights(2);
+
+  unsigned int n_intermediates=10;
+
+  // First, use only two surrounding points. The second time around,
+  // instead, use four, with additional zero weights
+  for (unsigned int test_no=0; test_no<2; ++test_no)
+    {
+      if (test_no == 1)
+        {
+          // Test that we get the same result with four surrounding points
+          points.push_back(Point<2>(0.3, 1.0));
+          points.push_back(Point<2>(0.7, 1.0));
+          weights.push_back(0.0);
+          weights.push_back(0.0);
+        }
+
+      deallog << "Test " << test_no << std::endl;
+      for (unsigned int i=0; i<n_intermediates; ++i)
+        {
+          weights[0] = (double)i/((double)n_intermediates-1);
+          weights[1] = 1.0-weights[0];
+
+          deallog << manifold.get_new_point(Quadrature<2>(points, weights)) << std::endl;
+        }
+    }
+
+  return 0;
+}
+
+
+
diff --git a/tests/manifold/chart_manifold_09.output b/tests/manifold/chart_manifold_09.output
new file mode 100644 (file)
index 0000000..2705afc
--- /dev/null
@@ -0,0 +1,23 @@
+
+DEAL::Test 0
+DEAL::0.900000 0.00000
+DEAL::0.922222 0.00000
+DEAL::0.944444 0.00000
+DEAL::0.966667 0.00000
+DEAL::0.988889 0.00000
+DEAL::0.0111111 0.00000
+DEAL::0.0333333 0.00000
+DEAL::0.0555556 0.00000
+DEAL::0.0777778 0.00000
+DEAL::0.100000 0.00000
+DEAL::Test 1
+DEAL::0.900000 0.00000
+DEAL::0.922222 0.00000
+DEAL::0.944444 0.00000
+DEAL::0.966667 0.00000
+DEAL::0.988889 0.00000
+DEAL::0.0111111 0.00000
+DEAL::0.0333333 0.00000
+DEAL::0.0555556 0.00000
+DEAL::0.0777778 0.00000
+DEAL::0.100000 0.00000

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.