From 8318bd6cb561d0dfcf5d6e450bf7fd6c5bbd9152 Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Fri, 18 Mar 2016 20:10:02 +0100 Subject: [PATCH] Added tangent vector for FunctionManifold --- doc/news/changes.h | 1 - include/deal.II/base/function_parser.h | 4 +- include/deal.II/grid/manifold_lib.h | 16 ++++ source/base/function_parser.cc | 2 +- source/grid/manifold_lib.cc | 14 +++ tests/manifold/function_manifold_03.cc | 96 +++++++++++++++++++ ...tion_manifold_03.with_muparser=true.output | 58 +++++++++++ 7 files changed, 187 insertions(+), 4 deletions(-) create mode 100644 tests/manifold/function_manifold_03.cc create mode 100644 tests/manifold/function_manifold_03.with_muparser=true.output diff --git a/doc/news/changes.h b/doc/news/changes.h index e103c7af5a..503c06e199 100644 --- a/doc/news/changes.h +++ b/doc/news/changes.h @@ -98,7 +98,6 @@ inconvenience this causes.

Specific improvements

-
  1. Improved: The distribution of degrees of freedom on multigrid levels, DoFHandler::distribute_mg_dofs(), contained a few steps that scaled diff --git a/include/deal.II/base/function_parser.h b/include/deal.II/base/function_parser.h index 5ae72919df..3fe91a0c4e 100644 --- a/include/deal.II/base/function_parser.h +++ b/include/deal.II/base/function_parser.h @@ -19,7 +19,7 @@ #include #include -#include +#include #include #include #include @@ -177,7 +177,7 @@ template class Vector; * @author Luca Heltai, Timo Heister 2005, 2014 */ template -class FunctionParser : public Function +class FunctionParser : public AutoDerivativeFunction { public: /** diff --git a/include/deal.II/grid/manifold_lib.h b/include/deal.II/grid/manifold_lib.h index 21de7b0555..15b3d7e4dc 100644 --- a/include/deal.II/grid/manifold_lib.h +++ b/include/deal.II/grid/manifold_lib.h @@ -265,6 +265,22 @@ public: virtual Point push_forward(const Point &chart_point) const; + /** + * Given a point in the chartdim dimensional Euclidean space, this + * method returns the derivatives of the function $F$ that maps from + * the sub_manifold coordinate system to the Euclidean coordinate + * system. In other words, it is a matrix of size + * $\text{spacedim}\times\text{chartdim}$. + * + * This function is used in the computations required by the + * get_tangent_vector() function. + * + * Refer to the general documentation of this class for more information. + */ + virtual + DerivativeForm<1,chartdim,spacedim> + push_forward_gradient(const Point &chart_point) const; + /** * Given a point in the spacedim coordinate system, uses the * pull_back_function to compute the pull_back of points in @p spacedim diff --git a/source/base/function_parser.cc b/source/base/function_parser.cc index 50703b27a6..094869b40f 100644 --- a/source/base/function_parser.cc +++ b/source/base/function_parser.cc @@ -48,7 +48,7 @@ template FunctionParser::FunctionParser(const unsigned int n_components, const double initial_time) : - Function(n_components, initial_time) + AutoDerivativeFunction(1e-8, n_components, initial_time) {} diff --git a/source/grid/manifold_lib.cc b/source/grid/manifold_lib.cc index 4740f4e65c..0d039f8c86 100644 --- a/source/grid/manifold_lib.cc +++ b/source/grid/manifold_lib.cc @@ -320,6 +320,20 @@ FunctionManifold::push_forward(const Point &cha } +template +DerivativeForm<1,chartdim, spacedim> +FunctionManifold::push_forward_gradient(const Point &chart_point) const +{ + DerivativeForm<1, chartdim, spacedim> DF; + std::vector > gradients(spacedim); + push_forward_function->vector_gradient(chart_point, gradients); + for (unsigned int i=0; i Point FunctionManifold::pull_back(const Point &space_point) const diff --git a/tests/manifold/function_manifold_03.cc b/tests/manifold/function_manifold_03.cc new file mode 100644 index 0000000000..a698f7ca1b --- /dev/null +++ b/tests/manifold/function_manifold_03.cc @@ -0,0 +1,96 @@ +//---------------------------- function_manifold_chart --------------------------- +// Copyright (C) 2011 - 2015 by the mathLab team. +// +// 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. +// +//---------------------------- function_manifold_chart --------------------------- + + +// Test a simple parabolic manifold, including gradients and tangent vector + +#include "../tests.h" +#include +#include + + +// all include files you need here +#include +#include +#include +#include +#include +#include +#include + +// Helper function +template +void test(unsigned int ref=1) +{ + deallog << "Testing dim " << dim + << ", spacedim " << spacedim << std::endl; + + // Here the only allowed axis is z. In cylinder the default is x. + std::string push_forward_expression; + std::string pull_back_expression; + + switch (spacedim) + { + case 2: + push_forward_expression = "x; x^2"; + pull_back_expression = "x"; + break; + case 3: + push_forward_expression = "x; x^2; 0"; + pull_back_expression = "x"; + break; + default: + Assert(false, ExcInternalError()); + } + + FunctionManifold manifold(push_forward_expression, + pull_back_expression); + + // Two points and two weights + std::vector > p(2); + p[1][0] = 1.0; + p[1][1] = 1.0; + std::vector w(2); + + unsigned int n_intermediates = 16; + + deallog << "P0: " << p[0] + << ", P1: " << p[1] << std::endl; + + for (unsigned int i=0; i ip = manifold.get_new_point(Quadrature(p, w)); + Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, p[0]); + Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, p[1]); + + deallog << "P: " << ip + << ", T(P, P0): " << t1 + << ", T(P, P1): " << t2 << std::endl; + + } +} + +int main () +{ + std::ofstream logfile("output"); + deallog.attach(logfile); + deallog.threshold_double(1.e-10); + + + test<2,2>(); + test<2,3>(); + test<3,3>(); + + return 0; +} + diff --git a/tests/manifold/function_manifold_03.with_muparser=true.output b/tests/manifold/function_manifold_03.with_muparser=true.output new file mode 100644 index 0000000000..c901d847a5 --- /dev/null +++ b/tests/manifold/function_manifold_03.with_muparser=true.output @@ -0,0 +1,58 @@ + +DEAL::Testing dim 2, spacedim 2 +DEAL::P0: 0.00000 0.00000, P1: 1.00000 1.00000 +DEAL::P: 0.00000 0.00000, T(P, P0): 0.00000 0.00000, T(P, P1): 1.00000 0.00000 +DEAL::P: 0.0625000 0.00390625, T(P, P0): -0.0625000 -0.00781250, T(P, P1): 0.937500 0.117187 +DEAL::P: 0.125000 0.0156250, T(P, P0): -0.125000 -0.0312500, T(P, P1): 0.875000 0.218750 +DEAL::P: 0.187500 0.0351562, T(P, P0): -0.187500 -0.0703125, T(P, P1): 0.812500 0.304687 +DEAL::P: 0.250000 0.0625000, T(P, P0): -0.250000 -0.125000, T(P, P1): 0.750000 0.375000 +DEAL::P: 0.312500 0.0976562, T(P, P0): -0.312500 -0.195312, T(P, P1): 0.687500 0.429687 +DEAL::P: 0.375000 0.140625, T(P, P0): -0.375000 -0.281250, T(P, P1): 0.625000 0.468750 +DEAL::P: 0.437500 0.191406, T(P, P0): -0.437500 -0.382812, T(P, P1): 0.562500 0.492187 +DEAL::P: 0.500000 0.250000, T(P, P0): -0.500000 -0.500000, T(P, P1): 0.500000 0.500000 +DEAL::P: 0.562500 0.316406, T(P, P0): -0.562500 -0.632813, T(P, P1): 0.437500 0.492188 +DEAL::P: 0.625000 0.390625, T(P, P0): -0.625000 -0.781250, T(P, P1): 0.375000 0.468750 +DEAL::P: 0.687500 0.472656, T(P, P0): -0.687500 -0.945313, T(P, P1): 0.312500 0.429688 +DEAL::P: 0.750000 0.562500, T(P, P0): -0.750000 -1.12500, T(P, P1): 0.250000 0.375000 +DEAL::P: 0.812500 0.660156, T(P, P0): -0.812500 -1.32031, T(P, P1): 0.187500 0.304688 +DEAL::P: 0.875000 0.765625, T(P, P0): -0.875000 -1.53125, T(P, P1): 0.125000 0.218750 +DEAL::P: 0.937500 0.878906, T(P, P0): -0.937500 -1.75781, T(P, P1): 0.0625000 0.117188 +DEAL::P: 1.00000 1.00000, T(P, P0): -1.00000 -2.00000, T(P, P1): 0.00000 0.00000 +DEAL::Testing dim 2, spacedim 3 +DEAL::P0: 0.00000 0.00000 0.00000, P1: 1.00000 1.00000 0.00000 +DEAL::P: 0.00000 0.00000 0.00000, T(P, P0): 0.00000 0.00000 0.00000, T(P, P1): 1.00000 0.00000 0.00000 +DEAL::P: 0.0625000 0.00390625 0.00000, T(P, P0): -0.0625000 -0.00781250 0.00000, T(P, P1): 0.937500 0.117187 0.00000 +DEAL::P: 0.125000 0.0156250 0.00000, T(P, P0): -0.125000 -0.0312500 0.00000, T(P, P1): 0.875000 0.218750 0.00000 +DEAL::P: 0.187500 0.0351562 0.00000, T(P, P0): -0.187500 -0.0703125 0.00000, T(P, P1): 0.812500 0.304687 0.00000 +DEAL::P: 0.250000 0.0625000 0.00000, T(P, P0): -0.250000 -0.125000 0.00000, T(P, P1): 0.750000 0.375000 0.00000 +DEAL::P: 0.312500 0.0976562 0.00000, T(P, P0): -0.312500 -0.195312 0.00000, T(P, P1): 0.687500 0.429687 0.00000 +DEAL::P: 0.375000 0.140625 0.00000, T(P, P0): -0.375000 -0.281250 0.00000, T(P, P1): 0.625000 0.468750 0.00000 +DEAL::P: 0.437500 0.191406 0.00000, T(P, P0): -0.437500 -0.382812 0.00000, T(P, P1): 0.562500 0.492187 0.00000 +DEAL::P: 0.500000 0.250000 0.00000, T(P, P0): -0.500000 -0.500000 0.00000, T(P, P1): 0.500000 0.500000 0.00000 +DEAL::P: 0.562500 0.316406 0.00000, T(P, P0): -0.562500 -0.632813 0.00000, T(P, P1): 0.437500 0.492188 0.00000 +DEAL::P: 0.625000 0.390625 0.00000, T(P, P0): -0.625000 -0.781250 0.00000, T(P, P1): 0.375000 0.468750 0.00000 +DEAL::P: 0.687500 0.472656 0.00000, T(P, P0): -0.687500 -0.945313 0.00000, T(P, P1): 0.312500 0.429688 0.00000 +DEAL::P: 0.750000 0.562500 0.00000, T(P, P0): -0.750000 -1.12500 0.00000, T(P, P1): 0.250000 0.375000 0.00000 +DEAL::P: 0.812500 0.660156 0.00000, T(P, P0): -0.812500 -1.32031 0.00000, T(P, P1): 0.187500 0.304688 0.00000 +DEAL::P: 0.875000 0.765625 0.00000, T(P, P0): -0.875000 -1.53125 0.00000, T(P, P1): 0.125000 0.218750 0.00000 +DEAL::P: 0.937500 0.878906 0.00000, T(P, P0): -0.937500 -1.75781 0.00000, T(P, P1): 0.0625000 0.117188 0.00000 +DEAL::P: 1.00000 1.00000 0.00000, T(P, P0): -1.00000 -2.00000 0.00000, T(P, P1): 0.00000 0.00000 0.00000 +DEAL::Testing dim 3, spacedim 3 +DEAL::P0: 0.00000 0.00000 0.00000, P1: 1.00000 1.00000 0.00000 +DEAL::P: 0.00000 0.00000 0.00000, T(P, P0): 0.00000 0.00000 0.00000, T(P, P1): 1.00000 0.00000 0.00000 +DEAL::P: 0.0625000 0.00390625 0.00000, T(P, P0): -0.0625000 -0.00781250 0.00000, T(P, P1): 0.937500 0.117187 0.00000 +DEAL::P: 0.125000 0.0156250 0.00000, T(P, P0): -0.125000 -0.0312500 0.00000, T(P, P1): 0.875000 0.218750 0.00000 +DEAL::P: 0.187500 0.0351562 0.00000, T(P, P0): -0.187500 -0.0703125 0.00000, T(P, P1): 0.812500 0.304687 0.00000 +DEAL::P: 0.250000 0.0625000 0.00000, T(P, P0): -0.250000 -0.125000 0.00000, T(P, P1): 0.750000 0.375000 0.00000 +DEAL::P: 0.312500 0.0976562 0.00000, T(P, P0): -0.312500 -0.195312 0.00000, T(P, P1): 0.687500 0.429687 0.00000 +DEAL::P: 0.375000 0.140625 0.00000, T(P, P0): -0.375000 -0.281250 0.00000, T(P, P1): 0.625000 0.468750 0.00000 +DEAL::P: 0.437500 0.191406 0.00000, T(P, P0): -0.437500 -0.382812 0.00000, T(P, P1): 0.562500 0.492187 0.00000 +DEAL::P: 0.500000 0.250000 0.00000, T(P, P0): -0.500000 -0.500000 0.00000, T(P, P1): 0.500000 0.500000 0.00000 +DEAL::P: 0.562500 0.316406 0.00000, T(P, P0): -0.562500 -0.632813 0.00000, T(P, P1): 0.437500 0.492188 0.00000 +DEAL::P: 0.625000 0.390625 0.00000, T(P, P0): -0.625000 -0.781250 0.00000, T(P, P1): 0.375000 0.468750 0.00000 +DEAL::P: 0.687500 0.472656 0.00000, T(P, P0): -0.687500 -0.945313 0.00000, T(P, P1): 0.312500 0.429688 0.00000 +DEAL::P: 0.750000 0.562500 0.00000, T(P, P0): -0.750000 -1.12500 0.00000, T(P, P1): 0.250000 0.375000 0.00000 +DEAL::P: 0.812500 0.660156 0.00000, T(P, P0): -0.812500 -1.32031 0.00000, T(P, P1): 0.187500 0.304688 0.00000 +DEAL::P: 0.875000 0.765625 0.00000, T(P, P0): -0.875000 -1.53125 0.00000, T(P, P1): 0.125000 0.218750 0.00000 +DEAL::P: 0.937500 0.878906 0.00000, T(P, P0): -0.937500 -1.75781 0.00000, T(P, P1): 0.0625000 0.117188 0.00000 +DEAL::P: 1.00000 1.00000 0.00000, T(P, P0): -1.00000 -2.00000 0.00000, T(P, P1): 0.00000 0.00000 0.00000 -- 2.39.5