Add a Mapping of a Geometry

Summary

Essentially, a geometry is a local functions that is differentiable and invertible, plus some extra utility functions for volume and integration element. It would be nice to have a way to chain these "functions", or to combine it with a classical differentiable function.

This MR adds a transformation/mapping of a geometry.

Examples

In the following three examples, functions defined on different domains with different interpretations of derivatives are shown.

1. In Example 1, a function that can be evaluated in global coordinates with a derivative w.r.t. global coordinates is used. This requires that we map first from the reference element to the grid element and also requires that the jacobian is transformed from global coordinates to local coordinates
2. In Example 2, the function is defined purely local, in terms of a local basis and a linear-combination of the basis functions. Thus, the evaluation is in local coordinates and the derivatives are already w.r.t. local coordinates. No additional transformation into a real element is required and we can use the reference-element geometry instead.
3. The last example, Example 3, is in between the first two examples. Thereby, the function can be evaluated in local coordinates, but the derivatives are defined w.r.t. global coordinates. Thus, we need to transform only the derivative, using a special wrapper around the reference-element geometry, returning the jacobianTransposed() from the real geometry.

Example 1

An analytic function for the parametrization:

// construct a base grid e.g. with piecewise flat elements
YaspGrid<2> grid({1.0,1.0}, {10,10});

using Domain = FieldVector<double,2>;
using Range = FieldVector<double,3>;

auto sig = Functions::SignatureTag<Range(Domain)>{};
auto f = Functions::makeDifferentiableFunctionFromCallables(sig, 
  [](const Domain& x) -> FieldVector<double,3> {
    return {x[0], x[1], std::sin(x[0] * x[1])};
  },
  [](const Domain& x) -> FieldMatrix<double,3,2> {
    return {
      {1.0, 0.0},
      {0.0, 1.0},
      {x[1] * std::cos(x[0] * x[1]), x[0] * std::cos(x[0] * x[1])}
    };
  });

for (const auto& e : elements(grid.leafGridView()))
{
  // construct the mapped geometry
  MappedGeometry geometry{f, e.geometry()};
}

Examples 2

A discrete function build from dune-localfunction local finite-element basis functions:

template <class LocalBasis, class Coefficients>
class LocalFunction
{
  // implementation of a discrete function based on linear-combination of basis functions
};

using LocalFiniteElement = LagrangeCubeLocalFiniteElement<double,double,2,order>;
using Range = FieldVector<double,3>;

// construct a base grid e.g. with piecewise flat elements
YaspGrid<2> grid({1.0,1.0}, {10,10});
LocalFiniteElement localFE{};

for (const auto& e : elements(grid.leafGridView()))
{
  auto X = [geo=e.geometry()](const auto& local) -> Range
  {
    auto x = geo.global(local);
    return {x[0], x[1], std::sin(x[0] * x[1])};
  };

  std::vector<Range> vertices;
  localFE.localInterpolation().interpolate(X, vertices);

  // construct the mapped geometry
  LocalFunction lf{localFE.localBasis(), vertices};
  MappedGeometry geometry{lf, ReferenceElementGeometry{referenceElement(e)}};
  // or
  // auto refElem = referenceElement(e);
  // MappendGeometry geometry{lf, refElem.template geometry<0>(0)};
}

Example 3

In case you get a local function created from a grid-function with a derivative in global coordinates, but the argument is already transformed by a geometry, you can use the additional wrapper LocalDerivativeGeometry:

// construct a base grid e.g. with piecewise flat elements
YaspGrid<2> grid({1.0,1.0}, {10,10});

using Domain = FieldVector<double,2>;
using Range = FieldVector<double,3>;

auto sig = Functions::SignatureTag<Range(Domain)>{};
auto f = Functions::makeDifferentiableFunctionFromCallables(sig, 
  [](const Domain& x) -> FieldVector<double,3> {
    return {x[0], x[1], std::sin(x[0] * x[1])};
  },
  [](const Domain& x) -> FieldMatrix<double,3,2> {
    return {
      {1.0, 0.0},
      {0.0, 1.0},
      {x[1] * std::cos(x[0] * x[1]), x[0] * std::cos(x[0] * x[1])}
    };
  });

// create a differentiable analytic grid function
auto fGridFunction = Functions::makeAnalyticGridViewFunction(f, grid.leafGridView());

// build the corresponding local function
auto fLocalFunction = localFunction(fGridFunction);

for (const auto& e : elements(grid.leafGridView()))
{
  // bind the local function to the grid element
  fLocalFunction.bind(e);

  // construct the mapped geometry
  MappedGeometry geometry{fLocalFunction, LocalDerivativeGeometry{e.geometry()}};

  // unbind the local function from the element
  fLocalFunction.unbind();
}