From e868d83048964f24cbb95ec46ab32ef975d5cd42 Mon Sep 17 00:00:00 2001 From: Luca Heltai Date: Thu, 24 Mar 2016 19:10:52 +0100 Subject: [PATCH] Initial CompositionManifold. --- include/deal.II/grid/manifold_tools.h | 115 ++++++++++++++++++ source/grid/CMakeLists.txt | 2 + source/grid/manifold_tools.cc | 83 +++++++++++++ source/grid/manifold_tools.inst.in | 24 ++++ tests/manifold/composition_manifold_01.cc | 75 ++++++++++++ ...tion_manifold_01.with_muparser=true.output | 19 +++ tests/manifold/composition_manifold_02.cc | 81 ++++++++++++ ...tion_manifold_02.with_muparser=true.output | 21 ++++ 8 files changed, 420 insertions(+) create mode 100644 include/deal.II/grid/manifold_tools.h create mode 100644 source/grid/manifold_tools.cc create mode 100644 source/grid/manifold_tools.inst.in create mode 100644 tests/manifold/composition_manifold_01.cc create mode 100644 tests/manifold/composition_manifold_01.with_muparser=true.output create mode 100644 tests/manifold/composition_manifold_02.cc create mode 100644 tests/manifold/composition_manifold_02.with_muparser=true.output diff --git a/include/deal.II/grid/manifold_tools.h b/include/deal.II/grid/manifold_tools.h new file mode 100644 index 0000000000..21341d3f06 --- /dev/null +++ b/include/deal.II/grid/manifold_tools.h @@ -0,0 +1,115 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + +#ifndef dealii__manifold_tools_h +#define dealii__manifold_tools_h + + +/*---------------------------- manifold_tools.h -----------------*/ + +#include +#include +#include +#include +#include +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + +/** + * CompositionManifold. Take two ChartManifold objects, and make + * their composition. The CompositionManifold object is a + * ChartManifold going from the chart of the first ChartManifold to + * the embedding space of the second ChartManifold. + * + * This class only works for dim <= chartdim <= intermediate_spacedim + * <= spacedim. If you try to instantiate anything different, an + * Exception will be thrown in one of the ChartManifold classes that + * violates this condition. + * + * Given the ChartManifold F and the ChartManifold G, this class + * represents the composition of G after F. + * + * The template parameters have the following meaning: + * + * @tparam dim The dimension of the resulting ChartManifold + * @tparam spacedim The space dimension of the resulting ChartManifold + * @tparam chartdim The chart dimension of the resulting ChartManifold + * @tparam intermediate_dim The space dimension of the first ChartManifold + * @tparam dim1 The dimension of the first ChartManifold, which coincides also + * with the chart dimension of the second ChartManifold + * @tparam dim2 The dimension of the second ChartManifold + * + * @ingroup manifold + * @author Luca Heltai, Timo Heister, 2016 + */ +template +class CompositionManifold : public ChartManifold +{ +public: + + /** + * Construct the composition of the two given manifolds. + */ + CompositionManifold(const ChartManifold &F, + const ChartManifold &G); + + /** + * Pull back the given point in spacedim to the Euclidean chartdim + * dimensional space. This function calls the pull_back() function + * of G, and then the pull_back() function of F. + */ + virtual + Point + pull_back(const Point &space_point) const; + + + /** + * Push forward the chartdim dimensional point to a spacedim + * Euclidean point. The function calls first the push_forward() of + * F, and then the push_foward() of G. + */ + virtual + Point + push_forward(const Point &chart_point) const; + + /** + * Return the derivative of the composition of G after F. + */ + virtual + DerivativeForm<1,chartdim,spacedim> + push_forward_gradient(const Point &chart_point) const; + +private: + /** + * The first ChartManifold. + */ + SmartPointer, + CompositionManifold > F; + + + /** + * The second ChartManifold. + */ + SmartPointer, + CompositionManifold > G; +}; + +DEAL_II_NAMESPACE_CLOSE + +#endif diff --git a/source/grid/CMakeLists.txt b/source/grid/CMakeLists.txt index 8279350bf2..18abd37211 100644 --- a/source/grid/CMakeLists.txt +++ b/source/grid/CMakeLists.txt @@ -26,6 +26,7 @@ SET(_src intergrid_map.cc manifold.cc manifold_lib.cc + manifold_tools.cc persistent_tria.cc tria_accessor.cc tria_boundary.cc @@ -45,6 +46,7 @@ SET(_inst intergrid_map.inst.in manifold.inst.in manifold_lib.inst.in + manifold_tools.inst.in tria_accessor.inst.in tria_boundary.inst.in tria_boundary_lib.inst.in diff --git a/source/grid/manifold_tools.cc b/source/grid/manifold_tools.cc new file mode 100644 index 0000000000..5cb2219fb8 --- /dev/null +++ b/source/grid/manifold_tools.cc @@ -0,0 +1,83 @@ +// --------------------------------------------------------------------- +// +// Copyright (C) 1998 - 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. +// +// --------------------------------------------------------------------- + +#include +#include +#include +#include +#include +#include +#include + +DEAL_II_NAMESPACE_OPEN + + +/* -------------------------- CompositionManifold --------------------- */ + + + +template +CompositionManifold::CompositionManifold +(const ChartManifold &F, + const ChartManifold &G) : + F(&F), + G(&G) +{} + +template +Point +CompositionManifold::pull_back(const Point &space_point) const +{ + return F->pull_back(G->pull_back(space_point)); +} + +template +Point +CompositionManifold::push_forward(const Point &chart_point) const +{ + return G->push_forward(F->push_forward(chart_point)); +} + + +template +DerivativeForm<1,chartdim,spacedim> +CompositionManifold::push_forward_gradient(const Point &chart_point) const +{ + DerivativeForm<1,chartdim,intermediate_dim> DF = + F->push_forward_gradient(chart_point); + + DerivativeForm<1,intermediate_dim,spacedim> DG = + G->push_forward_gradient(F->push_forward(chart_point)); + + DerivativeForm<1,chartdim,spacedim> DF_DG; + + for (unsigned int d=0; d; +#endif + } + diff --git a/tests/manifold/composition_manifold_01.cc b/tests/manifold/composition_manifold_01.cc new file mode 100644 index 0000000000..19af4633b4 --- /dev/null +++ b/tests/manifold/composition_manifold_01.cc @@ -0,0 +1,75 @@ +//---------------------------- 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. +// +//---------------------------- composition_manifold --------------------------- + + +// Test the combination of simple ChartManifolds: parabolic + translation + +#include "../tests.h" +#include +#include + + +// all include files you need here +#include +#include + + +int main () +{ + initlog(); + + const int dim=2, spacedim=2; + + FunctionManifold F("x;x^2", "x"); + FunctionManifold G("x;y+1", "x;y-1"); + + CompositionManifold manifold(F,G); + + // Chart points. + Point<1> cp[2]; + cp[1][0] = 1.0; + + // Spacedim points + std::vector > sp(2); + + // Weights + std::vector w(2); + + sp[0] = manifold.push_forward(cp[0]); + sp[1] = manifold.push_forward(cp[1]); + + for (unsigned int d=0; d<2; ++d) + if (cp[d].distance(manifold.pull_back(sp[d])) > 1e-10) + deallog << "Error!" << std::endl; + + unsigned int n_intermediates = 16; + + deallog << "P0: " << sp[0] + << ", P1: " << sp[1] << std::endl; + + for (unsigned int i=0; i ip = manifold.get_new_point(Quadrature(sp, w)); + Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); + Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); + + deallog << "P: " << ip + << ", T(P, P0): " << t1 + << ", T(P, P1): " << t2 << std::endl; + + } + + + return 0; +} + diff --git a/tests/manifold/composition_manifold_01.with_muparser=true.output b/tests/manifold/composition_manifold_01.with_muparser=true.output new file mode 100644 index 0000000000..3b65170f50 --- /dev/null +++ b/tests/manifold/composition_manifold_01.with_muparser=true.output @@ -0,0 +1,19 @@ + +DEAL::P0: 0.00000 1.00000, P1: 1.00000 2.00000 +DEAL::P: 0.00000 1.00000, T(P, P0): 0.00000 0.00000, T(P, P1): 1.00000 0.00000 +DEAL::P: 0.0625000 1.00391, T(P, P0): -0.0625000 -0.00781250, T(P, P1): 0.937500 0.117187 +DEAL::P: 0.125000 1.01562, T(P, P0): -0.125000 -0.0312500, T(P, P1): 0.875000 0.218750 +DEAL::P: 0.187500 1.03516, T(P, P0): -0.187500 -0.0703125, T(P, P1): 0.812500 0.304687 +DEAL::P: 0.250000 1.06250, T(P, P0): -0.250000 -0.125000, T(P, P1): 0.750000 0.375000 +DEAL::P: 0.312500 1.09766, T(P, P0): -0.312500 -0.195312, T(P, P1): 0.687500 0.429687 +DEAL::P: 0.375000 1.14062, T(P, P0): -0.375000 -0.281250, T(P, P1): 0.625000 0.468750 +DEAL::P: 0.437500 1.19141, T(P, P0): -0.437500 -0.382812, T(P, P1): 0.562500 0.492187 +DEAL::P: 0.500000 1.25000, T(P, P0): -0.500000 -0.500000, T(P, P1): 0.500000 0.500000 +DEAL::P: 0.562500 1.31641, T(P, P0): -0.562500 -0.632812, T(P, P1): 0.437500 0.492187 +DEAL::P: 0.625000 1.39062, T(P, P0): -0.625000 -0.781250, T(P, P1): 0.375000 0.468750 +DEAL::P: 0.687500 1.47266, T(P, P0): -0.687500 -0.945312, T(P, P1): 0.312500 0.429687 +DEAL::P: 0.750000 1.56250, T(P, P0): -0.750000 -1.12500, T(P, P1): 0.250000 0.375000 +DEAL::P: 0.812500 1.66016, T(P, P0): -0.812500 -1.32031, T(P, P1): 0.187500 0.304688 +DEAL::P: 0.875000 1.76562, T(P, P0): -0.875000 -1.53125, T(P, P1): 0.125000 0.218750 +DEAL::P: 0.937500 1.87891, T(P, P0): -0.937500 -1.75781, T(P, P1): 0.0625000 0.117187 +DEAL::P: 1.00000 2.00000, T(P, P0): -1.00000 -2.00000, T(P, P1): 0.00000 0.00000 diff --git a/tests/manifold/composition_manifold_02.cc b/tests/manifold/composition_manifold_02.cc new file mode 100644 index 0000000000..095f0a7480 --- /dev/null +++ b/tests/manifold/composition_manifold_02.cc @@ -0,0 +1,81 @@ +//---------------------------- 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. +// +//---------------------------- composition_manifold --------------------------- + + +// Test the combination of simple ChartManifolds: SphericalManifold + +// Translation + +#include "../tests.h" +#include +#include + + +// all include files you need here +#include +#include + + +int main () +{ + initlog(); + std::ostream &out = deallog.get_file_stream(); + + const int dim=2, spacedim=2; + + SphericalManifold<1,2> F; + FunctionManifold<2,2,2> G("x;y+1", "x;y-1"); + + CompositionManifold<2,2,2,2,1> manifold(F,G); + + // Chart points. + Point<2> cp[2]; + cp[0][0] = 1.0; + cp[1][0] = 1.0; + cp[1][1] = numbers::PI/2; + + // Spacedim points + std::vector > sp(2); + + // Weights + std::vector w(2); + + sp[0] = manifold.push_forward(cp[0]); + sp[1] = manifold.push_forward(cp[1]); + + for (unsigned int d=0; d<2; ++d) + if (cp[d].distance(manifold.pull_back(sp[d])) > 1e-10) + deallog << "Error!" << std::endl; + + unsigned int n_intermediates = 16; + + out << "set size ratio -1" << std::endl + << "plot '-' with vectors " << std::endl; + + for (unsigned int i=0; i ip = manifold.get_new_point(Quadrature(sp, w)); + Tensor<1,spacedim> t1 = manifold.get_tangent_vector(ip, sp[0]); + Tensor<1,spacedim> t2 = manifold.get_tangent_vector(ip, sp[1]); + + out << ip << " " + << t2 << std::endl + // << ip << " " + // << t2 << std::endl + ; + } + + out << "e" << std::endl; + + return 0; +} + diff --git a/tests/manifold/composition_manifold_02.with_muparser=true.output b/tests/manifold/composition_manifold_02.with_muparser=true.output new file mode 100644 index 0000000000..c54d0d7e10 --- /dev/null +++ b/tests/manifold/composition_manifold_02.with_muparser=true.output @@ -0,0 +1,21 @@ + +set size ratio -1 +plot '-' with vectors +1.00000 1.00000 0.00000 1.57080 +0.995185 1.09802 -0.144342 1.46553 +0.980785 1.19509 -0.268141 1.34804 +0.956940 1.29028 -0.370482 1.22132 +0.923880 1.38268 -0.450838 1.08842 +0.881921 1.47140 -0.509072 0.952407 +0.831470 1.55557 -0.545430 0.816293 +0.773010 1.63439 -0.560533 0.683011 +0.707107 1.70711 -0.555360 0.555360 +0.634393 1.77301 -0.531231 0.435970 +0.555570 1.83147 -0.489776 0.327258 +0.471397 1.88192 -0.432912 0.231396 +0.382683 1.92388 -0.362807 0.150279 +0.290285 1.95694 -0.281842 0.0854959 +0.195090 1.98079 -0.192577 0.0383059 +0.0980171 1.99518 -0.0977020 0.00962281 +6.12323e-17 2.00000 0.00000 0.00000 +e -- 2.39.5