From: Wolfgang Bangerth Date: Fri, 4 Mar 2016 21:00:35 +0000 (-0600) Subject: Add non-trivial tests for the tangent vector on a cylinder surface. X-Git-Tag: v8.5.0-rc1~1246^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F2305%2Fhead;p=dealii.git Add non-trivial tests for the tangent vector on a cylinder surface. --- diff --git a/tests/manifold/chart_manifold_07.cc b/tests/manifold/chart_manifold_07.cc new file mode 100644 index 0000000000..ba50b51b58 --- /dev/null +++ b/tests/manifold/chart_manifold_07.cc @@ -0,0 +1,164 @@ +//------------------------------------------------------------------- +// Copyright (C) 2016 by the deal.II authors. +// +// This file is subject to LGPL 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. +// +//------------------------------------------------------------------- + + +// Test direction vector on a cylinder surface + +#include "../tests.h" +#include +#include +#include + + +Point<3> periodicity (/*r=*/0, + /*phi=*/2*numbers::PI, + /*z=*/0); + +class MyCylinderManifold : public ChartManifold<2,3,3> +{ +public: + static const int dim = 2; + static const int spacedim = 3; + static const int chartdim = 3; + + MyCylinderManifold () + : + ChartManifold(periodicity) + {} + + + virtual + Point + pull_back(const Point &space_point) const + { + const double x = space_point[0]; + const double y = space_point[1]; + const double z = space_point[2]; + + const double r = std::sqrt(x*x + y*y); + const double phi = std::atan2(y,x); + + return Point<3>(r, + phi, + z); + } + + + virtual + Point + push_forward(const Point &chart_point) const + { + const double r = chart_point[0]; + const double phi = chart_point[1]; + const double z = chart_point[2]; + + return Point<3>(r*std::cos(phi), + r*std::sin(phi), + z); + } + + virtual + DerivativeForm<1,spacedim,spacedim> + push_forward_gradient(const Point &chart_point) const + { + DerivativeForm<1,spacedim,spacedim> g; + + const double r = chart_point[0]; + const double phi = chart_point[1]; + const double z = chart_point[2]; + + g[0][0] = std::cos(phi); + g[0][1] = -r*std::sin(phi); + g[0][2] = 0; + + g[1][0] = std::sin(phi); + g[1][1] = r*std::cos(phi); + g[1][2] = 0; + + g[2][0] = 0; + g[2][1] = 0; + g[2][2] = 1; + + return g; + } +}; + + + +void test_direction (const Point<3> &x1, + const Point<3> &x2) +{ + static MyCylinderManifold manifold; + + // check both the direction x1->x2 and x2->x1 + deallog << '[' << x1 << "] -> [" << x2 << "]: " + << manifold.get_tangent_vector (x1, x2) << std::endl; + deallog << '[' << x2 << "] -> [" << x1 << "]: " + << manifold.get_tangent_vector (x2, x1) << std::endl; +} + + +void test() +{ + MyCylinderManifold manifold; + + // check two points that are straight up and down from each other + test_direction (manifold.push_forward (Point<3>(/*r =*/4, + /*phi=*/0.1, + /*z =*/-1)), + manifold.push_forward (Point<3>(/*r =*/4, + /*phi=*/0.1, + /*z =*/+2))); + + // check two points that are radial + test_direction (manifold.push_forward (Point<3>(/*r =*/1, + /*phi=*/0.1, + /*z =*/-1)), + manifold.push_forward (Point<3>(/*r =*/4, + /*phi=*/0.1, + /*z =*/-1))); + + // check two points that are horizontal + test_direction (manifold.push_forward (Point<3>(/*r =*/4, + /*phi=*/0, + /*z =*/-1)), + manifold.push_forward (Point<3>(/*r =*/4, + /*phi=*/numbers::PI/4, + /*z =*/-1))); + + // same but rotated + test_direction (manifold.push_forward (Point<3>(/*r =*/4, + /*phi=*/numbers::PI/4, + /*z =*/-1)), + manifold.push_forward (Point<3>(/*r =*/4, + /*phi=*/numbers::PI/2, + /*z =*/-1))); + + + // check two points that are at the same radius but not horizontal + test_direction (manifold.push_forward (Point<3>(/*r =*/4, + /*phi=*/0, + /*z =*/-1)), + manifold.push_forward (Point<3>(/*r =*/4, + /*phi=*/numbers::PI/4, + /*z =*/1))); +} + +int main () +{ + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.threshold_double(1.e-10); + + test(); + + return 0; +} + diff --git a/tests/manifold/chart_manifold_07.output b/tests/manifold/chart_manifold_07.output new file mode 100644 index 0000000000..3012a86e52 --- /dev/null +++ b/tests/manifold/chart_manifold_07.output @@ -0,0 +1,11 @@ + +DEAL::[3.98002 0.399334 -1.00000] -> [3.98002 0.399334 2.00000]: 0.00000 0.00000 3.00000 +DEAL::[3.98002 0.399334 2.00000] -> [3.98002 0.399334 -1.00000]: 0.00000 0.00000 -3.00000 +DEAL::[0.995004 0.0998334 -1.00000] -> [3.98002 0.399334 -1.00000]: 2.98501 0.299500 0.00000 +DEAL::[3.98002 0.399334 -1.00000] -> [0.995004 0.0998334 -1.00000]: -2.98501 -0.299500 0.00000 +DEAL::[4.00000 0.00000 -1.00000] -> [2.82843 2.82843 -1.00000]: 0.00000 3.14159 0.00000 +DEAL::[2.82843 2.82843 -1.00000] -> [4.00000 0.00000 -1.00000]: 2.22144 -2.22144 0.00000 +DEAL::[2.82843 2.82843 -1.00000] -> [2.44929e-16 4.00000 -1.00000]: -2.22144 2.22144 0.00000 +DEAL::[2.44929e-16 4.00000 -1.00000] -> [2.82843 2.82843 -1.00000]: 3.14159 -1.92367e-16 0.00000 +DEAL::[4.00000 0.00000 -1.00000] -> [2.82843 2.82843 1.00000]: 0.00000 3.14159 2.00000 +DEAL::[2.82843 2.82843 1.00000] -> [4.00000 0.00000 -1.00000]: 2.22144 -2.22144 -2.00000 diff --git a/tests/manifold/chart_manifold_08.cc b/tests/manifold/chart_manifold_08.cc new file mode 100644 index 0000000000..b20297f2bc --- /dev/null +++ b/tests/manifold/chart_manifold_08.cc @@ -0,0 +1,141 @@ +//------------------------------------------------------------------- +// Copyright (C) 2016 by the deal.II authors. +// +// This file is subject to LGPL 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. +// +//------------------------------------------------------------------- + + +// Test direction vector on a cylinder surface. this is like the _07 +// test, but choose points in such a way that we wrap around the +// periodicity in phi + +#include "../tests.h" +#include +#include +#include + + +Point<3> periodicity (/*r=*/0, + /*phi=*/2*numbers::PI, + /*z=*/0); + +class MyCylinderManifold : public ChartManifold<2,3,3> +{ +public: + static const int dim = 2; + static const int spacedim = 3; + static const int chartdim = 3; + + MyCylinderManifold () + : + ChartManifold(periodicity) + {} + + + virtual + Point + pull_back(const Point &space_point) const + { + const double x = space_point[0]; + const double y = space_point[1]; + const double z = space_point[2]; + + const double r = std::sqrt(x*x + y*y); + const double phi = std::atan2(y,x); + + return Point<3>(r, + phi, + z); + } + + + virtual + Point + push_forward(const Point &chart_point) const + { + const double r = chart_point[0]; + const double phi = chart_point[1]; + const double z = chart_point[2]; + + return Point<3>(r*std::cos(phi), + r*std::sin(phi), + z); + } + + virtual + DerivativeForm<1,spacedim,spacedim> + push_forward_gradient(const Point &chart_point) const + { + DerivativeForm<1,spacedim,spacedim> g; + + const double r = chart_point[0]; + const double phi = chart_point[1]; + const double z = chart_point[2]; + + g[0][0] = std::cos(phi); + g[0][1] = -r*std::sin(phi); + g[0][2] = 0; + + g[1][0] = std::sin(phi); + g[1][1] = r*std::cos(phi); + g[1][2] = 0; + + g[2][0] = 0; + g[2][1] = 0; + g[2][2] = 1; + + return g; + } +}; + + + +void test_direction (const Point<3> &x1, + const Point<3> &x2) +{ + static MyCylinderManifold manifold; + + // check both the direction x1->x2 and x2->x1 + deallog << '[' << x1 << "] -> [" << x2 << "]: " + << manifold.get_tangent_vector (x1, x2) << std::endl; + deallog << '[' << x2 << "] -> [" << x1 << "]: " + << manifold.get_tangent_vector (x2, x1) << std::endl; +} + + +void test() +{ + MyCylinderManifold manifold; + + // check two points that are horizontal + test_direction (manifold.push_forward (Point<3>(/*r =*/2, + /*phi=*/3*numbers::PI/4, + /*z =*/-1)), + manifold.push_forward (Point<3>(/*r =*/2, + /*phi=*/-3*numbers::PI/4, + /*z =*/-1))); + + // same but rotated + test_direction (manifold.push_forward (Point<3>(/*r =*/2, + /*phi=*/-numbers::PI/4, + /*z =*/-1)), + manifold.push_forward (Point<3>(/*r =*/2, + /*phi=*/numbers::PI/4, + /*z =*/-1))); +} + +int main () +{ + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.threshold_double(1.e-10); + + test(); + + return 0; +} + diff --git a/tests/manifold/chart_manifold_08.output b/tests/manifold/chart_manifold_08.output new file mode 100644 index 0000000000..e54409188c --- /dev/null +++ b/tests/manifold/chart_manifold_08.output @@ -0,0 +1,5 @@ + +DEAL::[-1.41421 1.41421 -1.00000] -> [-1.41421 -1.41421 -1.00000]: -2.22144 -2.22144 0.00000 +DEAL::[-1.41421 -1.41421 -1.00000] -> [-1.41421 1.41421 -1.00000]: -2.22144 2.22144 0.00000 +DEAL::[1.41421 -1.41421 -1.00000] -> [1.41421 1.41421 -1.00000]: 2.22144 2.22144 0.00000 +DEAL::[1.41421 1.41421 -1.00000] -> [1.41421 -1.41421 -1.00000]: 2.22144 -2.22144 0.00000