GridFunctions

Summary

GridFunctions are objects that can be evaluated in global coordinates and can be restricted to grid elements and evaluated there in local coordinates.

ProblemStat<Traits> prob("name");
prob.initialize(INIT_ALL);

// create a grid function
auto gridFct = makeGridFunction(EXPRESSION, prob.gridView());

// eval GridFunction at global coordinates
auto global = Dune::FieldVector<double,2>{1.0, 2.0};
auto value = gridFct(global);

// create a local representation of the grid function
auto localFct = localFunction(gridFct);
for (auto const& element : elements(prob.gridView()))
{
  // bind the local function to a grid element
  localFct.bind(element);

  // eval LocalFunction at local coordinates
  auto local = element.geometry().center();
  value = localFct(local);

  localFct.unbind();
}

GridFunctions are build from expressions that are a composition of some elementary terms.

Examples of expressions

Before we give examples where and how to use GridFunctions, we demonstrate what an Expression could be, to create a GridFunction from. In the following examples, we assume that a ProblemStat named prob is already created and initialized.

1. Discrete Functions

Components of a solution vector are discrete functions, i.e. grid function build by a linear combination of local basis functions with a vector of coefficients.

auto expr1 = prob.solution(0);

2. Constants and functions

Constant values (or reference to values) and functions of global coordinates are expressions

auto expr2 = 1.0;

double var = 1.0;
auto expr3 = std::ref(var);

auto expr4 = [](auto const& x) { return x[0] + x[1]; };

Warning

Combining std::ref and constant expressions could result in unexpected behavior. Note that a reference wrapper implicitly converts to the underlying data type and thus std::ref<double> + double results in double and the reference is gone. Use a lambda function and a function wrapper instead.

3. Composition of expressions

Combining expressions using arithmetic operations or some mathematical functions build a new expressions:

auto expr5 = prob.solution(0) + 1.0;

auto expr6 = max(evalAtQP(expr4), expr3);

This works, if at least one of the expressions involved is already a grid functions, so one can not simply add two lambda functions

// ERROR:
auto no_expr = expr4 + 1.0;

Examples of usage

In the following examples, an Expression is anything a GridFunction can be created from, sometimes also called PreGridFunction. It includes constants, functors callable with global coordinates, and any combination of GridFunctions.

1. Usage of GridFunctions to build Operators:

ProblemStat<Traits> prob("name");
prob.initialize(INIT_ALL);

auto opB = makeOperator(BiLinearForm, Expression);
prob.addMatrixOperator(opB, Row, Col);

auto opL = makeOperator(LinearForm, Expression);
prob.addVectorOperator(opL, Row);
See also makeOperator().

2. Usage of GridFunctions in BoundaryConditions:

prob.addDirichletBC(Predicate, Row, Col, Expression);

3. Interpolate a GridFunction to a DOFVector:

prob.solution(_0).interpol(Expression);
prob.solution() << Expression;

4. Integrate a GridFunction on a GridView:

auto value = integrate(Expression, prob.gridView());

List of Grid Functions and Expressions

Function Descriptions
evalAtQP(f) Evaluates a functor in global coordinates.
gradientOf(gf) Differentiate a grid function w.r.t. global coords
invokeAtQP(f,gf_0,gf_1...) Apply a functor to the evaluated grid functions
X(), X(comp) Return (component of) the global coordinate
constant Any (constant) scalar value or FieldVector or FieldMatrix
std::ref(variable) A reference to an lvalue object
+, -, *, / Arithmetic expressions
get(gf,i), get<i>(gf) Retrieve the i-th component of a vector expression: g(x)i g(x)_i
cmath Function Descriptions
max(gf_0,gf_1) Maximum of two grid functions: max{g0(x),g1(x)} \max\{g_0(x),g_1(x)\}
min(gf_0,gf_1) Minimum of two grid functions: min{g0(x),g1(x)} \min\{g_0(x),g_1(x)\}
abs_max(gf_0,gf_1) Maximum of the absolute value of two grid functions: max{g0(x),g1(x)} \max\{\vert g_0(x)\vert,\vert g_1(x)\vert\}
abs_min(gf_0,gf_1) Minimum of the absolute value of two grid functions: min{g0(x),g1(x)} \min\{\vert g_0(x)\vert,\vert g_1(x)\vert\}
clamp(gf,a,b) A grid functions clamped btween two values a and b: max{a,min{b,g(x)}} \max\{a,\min\{b,g(x)\}\}
abs(gf) The absolute value of a grid functions: g(x) \vert g(x)\vert
sqr(gf) The square value of a grid functions: g(x)2 g(x)^2
pow(gf,p), pow<p>(gf) A grid functions taken to the power of p: g(x)p g(x)^p
vector/matrix Function Descriptions
sum(gf) The sum of the components of a vector-valued GridFunction: ig(x)i \sum_i g(x)_i
unary_dot(gf) The inner product of a vector-valued GridFunction with itself: ig(x)ig(x)i \sum_i g(x)_i\cdot g(x)_i
one_norm(gf) The 1-norm of a vector-valued GridFunction: ig(x)i \sum_i \vert g(x)_i\vert
two_norm(gf) The 2-norm of a vector-valued GridFunction: ig(x)i2 \sqrt{\sum_i \vert g(x)_i\vert^2}
p_norm<p>(gf) The p-norm of a vector-valued GridFunction: (ig(x)ip)1/p (\sum_i \vert g(x)_i\vert^p)^{1/p}
infty_norm(gf) The infty-norm of a vector-valued GridFunction: maxig(x)i \max_i \vert g(x)_i\vert
trans(gf) The transposed of a matrix-valued GridFunction: G(x)T G(x)^T
dot(gf_0,gf_1) The inner product of two vector-valued GridFunctions: ig0(x)ig1(x)i \sum_i g_0(x)_i\cdot g_1(x)_i
distance(gf_0,gf_1) The euclidean distance of two vector-valued GridFunctions: ig0(x)ig1(x)i2 \sqrt{\sum_i \vert g_0(x)_i - g_1(x)_i\vert^2}
outer(gf_0,gf_1) The outer product of two vector-valued GridFunctions: G(x)ij=g0(x)ig1(x)j G(x)_{ij} = g_0(x)_i\cdot g_1(x)_j

function evalAtQP()

Defined in header <amdis/gridfunctions/AnalyticGridFunction.hpp>

template <class Function>
auto evalAtQP(Function const& f);

Creates a GridFunction that evaluates a functor in global coordinates.

Arguments

Function f
A callable object that can be called with global coordinates, i.e. typically FieldVector<double, dow> where dow is the world dimension. The range-type of the function determines the range-type of the expression.

Requirements

  • Function models Concepts::CallableDomain

function X()

Defined in header <amdis/gridfunctions/CoordsGridFunction.hpp>

inline auto X();          // (1)

inline auto X(int comp);  // (2)

A GridFunction that evaluates to the global coordinates (1) or a component of the global coordinates (2).

Arguments

int comp
The component of the global coordinate vector. Should be 0 <= comp < dow where dow is the world dimension.

function gradientOf()

Defined in header <amdis/gridfunctions/DerivativeGridFunction.hpp>

template <class Expr>
auto gradientOf(Expr const& expr);

Creates a GridFunction representing the gradient w.r.t. global coordinates of the wrapped expression.

Arguments

Expr expr
An expression that can be transformed into a grid function

Requirements

  • The GridFunction of gf = makeGridFunction(expr, gridView)with some GridView models Concepts::GridFunction and its LocalFunction has a derivative, i.e. the expression derivativeOf(localFunction(gf),tag::gradient{}) does not fail.

Examples

We assume there is a ProblemStat with the name prob.

gradientOf(prob.solution(_0))
gradientOf(X(0) + X(1) + prob.solution(_0))

function invokeAtQP()

Defined in header <amdis/gridfunctions/FunctorGridFunction.hpp>

template <class Functor, class... Exprs>
auto invokeAtQP(Functor const& f, Exprs&&... g_i);

Creates a Gridfunction that applies a functor to the evaluated expressions. It creates a composition of GridFunctions g_i by applying a functor f locally, i.e. locally it is evaluated f(g0(x),g1(x),...) f(g_0(x), g_1(x), ...)

Arguments

Functor f
A Functor f(...)f(...) that accepts as many arguments as there are grid functions
Exprs g_i...
A number of expressions that can be transformed into grid functions gig_i

Requirements

  • arity(f) == sizeof...(Exprs)
  • arity(g_i) == arity(g_j) for i != j
  • g_i models concept Concepts::GridFunction

Note

The composition of grid functions can be differentiated using gradientOf() if all the grid functions are differentiable and the functor is differentiable, i.e. there exist valid functions derivativeOf(g_i,type) and partial(f, i).

Examples

invokeAtQP([](Dune::FieldVector<double, 2> const& x) { return two_norm(x); }, X());
invokeAtQP([](double u, auto const& x) { return u + x[0]; }, 1.0, X());
invokeAtQP(Operation::Plus{}, X(0), X(1));