diff --git a/NumLib/NamedFunction.cpp b/NumLib/NamedFunction.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..fbb5124beec0fea0c632a86e0b6541b2bd03a29f
--- /dev/null
+++ b/NumLib/NamedFunction.cpp
@@ -0,0 +1,198 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "NamedFunction.h"
+#include "BaseLib/TMPUtil.h"
+
+//! Helper struct used in conjunction with BaseLib::IntegerSequence to get a
+//! sequence of types, where each type is double.
+template <int>
+struct Double
+{
+    using type = double;
+};
+
+/*! Calls the given \c function with the given  \c arguments.
+ *
+ * \tparam Indices sequence of integers used to expand the \c arguments vector
+ * to individual arguments.
+ */
+template <int... Indices>
+double call_(void* function, std::vector<double> const& arguments)
+{
+    assert(arguments.size() == sizeof...(Indices));
+    auto* fct = reinterpret_cast<
+        std::function<double(typename Double<Indices>::type...)>*>(function);
+    auto* args = arguments.data();
+    return (*fct)(args[Indices]...);
+}
+
+typedef double (*CallerFunction) (void*, std::vector<double> const&);
+
+//! Helps instantiating the call_() function.
+template <int... Indices>
+CallerFunction generateCallerFromIntegerSequence(
+    BaseLib::IntegerSequence<Indices...>)
+{
+    return call_<Indices...>;
+}
+
+//! Instantiates the call_() function for the provided number of arguments.
+template <int NumArguments>
+CallerFunction generateCaller()
+{
+    return generateCallerFromIntegerSequence(
+        typename BaseLib::GenerateIntegerSequence<NumArguments>::type{});
+}
+
+//! Holds instantiations of the call_() function for various numbers of
+//! arguments.
+static const CallerFunction callers[] = {
+    generateCaller<0>(),  generateCaller<1>(),  generateCaller<2>(),
+    generateCaller<3>(),  generateCaller<4>(),  generateCaller<5>(),
+    generateCaller<6>(),  generateCaller<7>(),  generateCaller<8>(),
+    generateCaller<9>(),  generateCaller<10>(), generateCaller<11>(),
+    generateCaller<12>(), generateCaller<13>(), generateCaller<14>(),
+    generateCaller<15>(), generateCaller<16>(), generateCaller<17>(),
+    generateCaller<18>(), generateCaller<19>(), generateCaller<20>(),
+    generateCaller<21>(), generateCaller<22>(), generateCaller<23>(),
+    generateCaller<24>(), generateCaller<25>(), generateCaller<26>(),
+    generateCaller<27>(), generateCaller<28>(), generateCaller<29>(),
+    generateCaller<30>(), generateCaller<31>(), generateCaller<32>()};
+static_assert(sizeof(callers) / sizeof(CallerFunction) ==
+                  NumLib::NamedFunction::MAX_FUNCTION_ARGS + 1,
+              "You did not instantiate the right number of callers.");
+
+/*! Deletes the given \c function.
+ *
+ * \tparam Indices sequence of integers used to cast the \c function to the
+ * correct type.
+ */
+template <int... Indices>
+void delete_(void* function)
+{
+    auto* fct = reinterpret_cast<
+        std::function<double(typename Double<Indices>::type...)>*>(function);
+    delete fct;
+}
+
+typedef void (*DeleterFunction) (void*);
+
+//! Helps instantiating the delete_() function.
+template <int... Indices>
+DeleterFunction generateDeleterFromIntegerSequence(
+    BaseLib::IntegerSequence<Indices...>)
+{
+    return delete_<Indices...>;
+}
+
+//! Instantiates the delete_() function for the provided number of arguments.
+template <int NumArguments>
+DeleterFunction generateDeleter()
+{
+    return generateDeleterFromIntegerSequence(
+        typename BaseLib::GenerateIntegerSequence<NumArguments>::type{});
+}
+
+//! Holds instantiations of the delete_() function for various numbers of
+//! arguments.
+static const DeleterFunction deleters[] = {
+    generateDeleter<0>(),  generateDeleter<1>(),  generateDeleter<2>(),
+    generateDeleter<3>(),  generateDeleter<4>(),  generateDeleter<5>(),
+    generateDeleter<6>(),  generateDeleter<7>(),  generateDeleter<8>(),
+    generateDeleter<9>(),  generateDeleter<10>(), generateDeleter<11>(),
+    generateDeleter<12>(), generateDeleter<13>(), generateDeleter<14>(),
+    generateDeleter<15>(), generateDeleter<16>(), generateDeleter<17>(),
+    generateDeleter<18>(), generateDeleter<19>(), generateDeleter<20>(),
+    generateDeleter<21>(), generateDeleter<22>(), generateDeleter<23>(),
+    generateDeleter<24>(), generateDeleter<25>(), generateDeleter<26>(),
+    generateDeleter<27>(), generateDeleter<28>(), generateDeleter<29>(),
+    generateDeleter<30>(), generateDeleter<31>(), generateDeleter<32>()};
+static_assert(sizeof(deleters) / sizeof(DeleterFunction) ==
+                  NumLib::NamedFunction::MAX_FUNCTION_ARGS + 1,
+              "You did not instantiate the right number of deleters.");
+
+/*! Copies the given \c function.
+ *
+ * \tparam Indices sequence of integers used to cast the \c function to the
+ * correct type.
+ */
+template <int... Indices>
+void* copy_(void* function)
+{
+    auto* fct = reinterpret_cast<
+        std::function<double(typename Double<Indices>::type...)>*>(function);
+    return new std::function<double(typename Double<Indices>::type...)>(*fct);
+}
+
+typedef void* (*CopierFunction) (void*);
+
+//! Helps instantiating the copy_() function.
+template <int... Indices>
+CopierFunction generateCopierFromIntegerSequence(
+    BaseLib::IntegerSequence<Indices...>)
+{
+    return copy_<Indices...>;
+}
+
+//! Instantiates the copy_() function for the provided number of arguments.
+template <int NumArguments>
+CopierFunction generateCopier()
+{
+    return generateCopierFromIntegerSequence(
+        typename BaseLib::GenerateIntegerSequence<NumArguments>::type{});
+}
+
+//! Holds instantiations of the copy_() function for various numbers of
+//! arguments.
+static const CopierFunction copiers[] = {
+    generateCopier<0>(),  generateCopier<1>(),  generateCopier<2>(),
+    generateCopier<3>(),  generateCopier<4>(),  generateCopier<5>(),
+    generateCopier<6>(),  generateCopier<7>(),  generateCopier<8>(),
+    generateCopier<9>(),  generateCopier<10>(), generateCopier<11>(),
+    generateCopier<12>(), generateCopier<13>(), generateCopier<14>(),
+    generateCopier<15>(), generateCopier<16>(), generateCopier<17>(),
+    generateCopier<18>(), generateCopier<19>(), generateCopier<20>(),
+    generateCopier<21>(), generateCopier<22>(), generateCopier<23>(),
+    generateCopier<24>(), generateCopier<25>(), generateCopier<26>(),
+    generateCopier<27>(), generateCopier<28>(), generateCopier<29>(),
+    generateCopier<30>(), generateCopier<31>(), generateCopier<32>()};
+static_assert(sizeof(copiers) / sizeof(CopierFunction) ==
+                  NumLib::NamedFunction::MAX_FUNCTION_ARGS + 1,
+              "You did not instantiate the right number of deleters.");
+
+namespace NumLib
+{
+NamedFunction::NamedFunction(NamedFunction&& other)
+    : _name(std::move(other._name)),
+      _argument_names(std::move(other._argument_names)),
+      _function(other._function)
+{
+    other._function = nullptr;
+}
+
+NamedFunction::NamedFunction(NamedFunction const& other)
+    : _name(other._name),
+      _argument_names(other._argument_names),
+      _function(copiers[_argument_names.size()](other._function))
+{
+}
+
+NamedFunction::~NamedFunction()
+{
+    deleters[_argument_names.size()](_function);
+}
+
+double NamedFunction::call(const std::vector<double>& arguments) const
+{
+    assert(arguments.size() == _argument_names.size());
+    return callers[_argument_names.size()](_function, arguments);
+}
+
+}  // namespace NumLib
diff --git a/NumLib/NamedFunction.h b/NumLib/NamedFunction.h
new file mode 100644
index 0000000000000000000000000000000000000000..5cdea7ec71f38046b75bfd37cbd8675558372d2e
--- /dev/null
+++ b/NumLib/NamedFunction.h
@@ -0,0 +1,116 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef NUMLIB_NAMED_FUNCTION
+#define NUMLIB_NAMED_FUNCTION
+
+#include <cassert>
+#include <functional>
+#include <memory>
+#include <string>
+#include <type_traits>
+#include <vector>
+
+namespace detail
+{
+//! Checks if all types \c Ts are the same as \c TRef.
+template<typename TRef, typename... Ts>
+struct AllTypesSameAs;
+
+template<typename TRef, typename T1, typename... Ts>
+struct AllTypesSameAs<TRef, T1, Ts...>
+{
+    static const bool value = std::is_same<TRef, T1>::value
+        && AllTypesSameAs<TRef, Ts...>::value;
+};
+
+template<typename TRef>
+struct AllTypesSameAs<TRef>
+{
+    static const bool value = true;
+};
+
+template<typename TRef, typename T1, typename... Ts>
+const bool AllTypesSameAs<TRef, T1, Ts...>::value;
+
+template<typename TRef>
+const bool AllTypesSameAs<TRef>::value;
+
+} // namespace detail
+
+namespace NumLib
+{
+//! Stores a function object along with a name for it and information about its
+//! arguments.
+class NamedFunction final
+{
+public:
+    /*! Constructs a new named function.
+     *
+     * \param name the function's name
+     * \param argument_names names  of arguments of the function
+     * \param function the actual function object
+     */
+    template <typename ReturnType, typename... Arguments>
+    NamedFunction(std::string const& name,
+                  std::vector<std::string>&& argument_names,
+                  std::function<ReturnType(Arguments...)>&& function);
+
+    NamedFunction(NamedFunction&& other);
+    NamedFunction(NamedFunction const&);
+
+    ~NamedFunction();
+
+    //! Returns the function's name.
+    std::string const& getName() const { return _name; }
+    //! Returns the names of the function's arguments.
+    std::vector<std::string> const& getArgumentNames() const
+    {
+        return _argument_names;
+    }
+
+    //! Call the function with the supplied arguments.
+    double call(std::vector<double> const& arguments) const;
+
+    //! Maximum number of function arguments supported by NamedFunction.
+    static const int MAX_FUNCTION_ARGS = 32;
+
+private:
+    //! The function's name.
+    std::string _name;
+
+    //! Information about the function's arguments.
+    std::vector<std::string> _argument_names;
+
+    //! The function handle.
+    void* _function;
+};
+
+template <typename ReturnType, typename... Arguments>
+NamedFunction::NamedFunction(std::string const& name,
+                             std::vector<std::string>&& argument_names,
+                             std::function<ReturnType(Arguments...)>&& function)
+    : _name(name),
+      _argument_names(std::move(argument_names)),
+      _function(
+          new std::function<ReturnType(Arguments...)>(std::move(function)))
+{
+    static_assert(
+        ::detail::AllTypesSameAs<double, ReturnType, Arguments...>::value,
+        "Some of the arguments or the return type of the passed function are "
+        "not of the type double.");
+    static_assert(sizeof...(Arguments) <= MAX_FUNCTION_ARGS,
+                  "The function you passed has too many arguments.");
+
+    assert(sizeof...(Arguments) == _argument_names.size());
+}
+
+}  // namespace NumLib
+
+#endif // NUMLIB_NAMED_FUNCTION
diff --git a/NumLib/NamedFunctionCaller.cpp b/NumLib/NamedFunctionCaller.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..b7111612e0f45bf3bee11137e8ef745d16dd414c
--- /dev/null
+++ b/NumLib/NamedFunctionCaller.cpp
@@ -0,0 +1,310 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#include "NamedFunctionCaller.h"
+
+#include <algorithm>
+
+#include "BaseLib/uniqueInsert.h"
+#include "BaseLib/Error.h"
+
+/// To understand the working of the algorithm the following two lemmata are
+/// useful:
+/// If a directed graph has a topological order, then it is a Directed Acyclic
+/// Graph (DAG).
+/// If a graph is a DAG, then it has a node with no incoming edges.
+///
+/// \note Negative indices in the graph adjacency list are leaves of the graph.
+///
+/// \param graph represents directed graph given as an adjacency list.
+/// \return true if the given directed graph is acyclic.
+bool hasTopologicalOrdering(std::vector<std::vector<int>> const& graph)
+{
+    // Number of incoming edges for each node of the graph
+    std::vector<unsigned> number_incoming_edges(graph.size());
+
+    // Count all incoming edges (i,j) for each node (i).
+    for (auto const& node_i_adjacencies : graph)
+    {
+        for (int node_j : node_i_adjacencies)
+        {
+            if (node_j >= 0)    // ignore negative indices.
+                ++number_incoming_edges[node_j];
+        }
+    }
+
+    // Working queue: a set of nodes with no incoming edges.
+    std::vector<std::size_t> q;
+    for (std::size_t node_i = 0; node_i < number_incoming_edges.size();
+         ++node_i)
+    {
+        if (number_incoming_edges[node_i] == 0)
+        {
+            q.push_back(node_i);
+        }
+    }
+
+
+    auto num_dependent = number_incoming_edges.size() - q.size();
+
+    while (!q.empty())
+    {
+        auto const node_i = q.back();
+        q.pop_back();
+
+        // Decrement counts for all edges (i,j).
+        for (int node_j : graph[node_i])
+        {
+            if (node_j < 0)
+                continue;  // ignore negative indices
+            if (--number_incoming_edges[node_j] == 0)
+            {  // Add a node without incoming edges to the queue.
+                q.push_back(node_j);
+                --num_dependent;
+            }
+        }
+    }
+
+    return num_dependent == 0;
+}
+
+enum class TraversePosition { StartNode, BetweenChildren, EndNode };
+
+/*! Traverses the graph given by the adjacency list \c map_sink_source in a
+ * depth-first manner starting at the node \c sink_fct.
+ *
+ * At the beginning and end of each node as well as between every two child
+ * nodes the given \c callback is called.
+ */
+template <typename Callback>
+void traverse(std::vector<std::vector<int>> const& map_sink_source,
+              int sink_fct, Callback&& callback)
+{
+    assert(sink_fct < static_cast<int>(map_sink_source.size()));
+    callback(sink_fct, TraversePosition::StartNode);
+
+    if (sink_fct < 0)
+        return;
+
+    auto const& si_so = map_sink_source[sink_fct];
+    std::size_t const num_args = si_so.size();
+    for (std::size_t sink_arg = 0; sink_arg != num_args; ++sink_arg) {
+        if (sink_arg != 0)
+            callback(sink_fct, TraversePosition::BetweenChildren);
+        traverse(map_sink_source, si_so[sink_arg], callback);
+    }
+
+    callback(sink_fct, TraversePosition::EndNode);
+}
+
+namespace NumLib
+{
+NamedFunctionCaller::NamedFunctionCaller(
+    std::initializer_list<std::string> unbound_argument_names)
+    : _uninitialized(-1 - unbound_argument_names.size())
+{
+    int idx = -1;
+    for (auto arg : unbound_argument_names) {
+        BaseLib::insertIfKeyUniqueElseError(
+            _map_name_idx, arg, idx,
+            "The name of the unbound argument is not unique.");
+        --idx;
+    }
+}
+
+void NamedFunctionCaller::addNamedFunction(NamedFunction&& fct)
+{
+    DBUG("Adding named function `%s'", fct.getName().c_str());
+
+    BaseLib::insertIfKeyUniqueElseError(
+        _map_name_idx, fct.getName(), _named_functions.size(),
+        "The name of the function is not unique.");
+
+    _map_sink_source.emplace_back(fct.getArgumentNames().size(),
+                                  _uninitialized);
+    _named_functions.push_back(std::move(fct));
+}
+
+void NamedFunctionCaller::plug(const std::string& sink_fct,
+                               const std::string& sink_arg,
+                               const std::string& source_fct)
+{
+    _deferred_plugs.push_back({sink_fct, sink_arg, source_fct});
+}
+
+void NamedFunctionCaller::applyPlugs()
+{
+    while (!_deferred_plugs.empty())
+    {
+        auto const& plug = _deferred_plugs.back();
+        auto const& sink_fct = plug.sink_fct;
+        auto const& sink_arg = plug.sink_arg;
+        auto const& source = plug.source;
+
+        auto const source_it = _map_name_idx.find(source);
+        if (source_it == _map_name_idx.end())
+        {
+            OGS_FATAL("A function with the name `%s' has not been found.",
+                      source.c_str());
+        }
+        auto const source_idx = source_it->second;
+
+        auto const sink_it = _map_name_idx.find(sink_fct);
+        if (sink_it == _map_name_idx.end())
+        {
+            OGS_FATAL("A function with the name `%s' has not been found.",
+                      sink_fct.c_str());
+        }
+        auto const sink_fct_idx = sink_it->second;
+
+        auto const& sink_args =
+            _named_functions[sink_it->second].getArgumentNames();
+        auto const sink_arg_it =
+            std::find(sink_args.begin(), sink_args.end(), sink_arg);
+        if (sink_arg_it == sink_args.end())
+        {
+            OGS_FATAL(
+                "An argument with the name `%s' has not been found for the "
+                "function `%s'.",
+                sink_arg.c_str(), sink_fct.c_str());
+        }
+        std::size_t const sink_arg_idx =
+            std::distance(sink_args.begin(), sink_arg_it);
+
+        auto& sis_sos = _map_sink_source[sink_fct_idx];
+        if (sis_sos[sink_arg_idx] != _uninitialized)
+        {
+            OGS_FATAL("A dependency for `%s'.`%s' has already been introduced.",
+                      sink_fct.c_str(), sink_arg.c_str());
+        }
+        sis_sos[sink_arg_idx] = source_idx;
+        if (!hasTopologicalOrdering(_map_sink_source))
+        {
+            OGS_FATAL(
+                "The call graph being plugged together must be an acyclic "
+                "graph. The added dependency for `%s'.`%s' introduces a cycle "
+                "into the graph.",
+                sink_fct.c_str(), sink_arg.c_str());
+        }
+
+        _deferred_plugs.pop_back();
+    }
+}
+
+double NamedFunctionCaller::call(
+    std::size_t function_idx,
+    const std::vector<double>& unbound_arguments) const
+{
+    assert(_deferred_plugs.empty() &&
+           "You must call applyPlugs() before this method!");
+
+    DBUG("Preparing call of fct #%lu %s()", function_idx,
+         _named_functions[function_idx].getName().c_str());
+    auto const& sis_sos = _map_sink_source[function_idx];
+    assert(sis_sos.size() ==
+           _named_functions[function_idx].getArgumentNames().size());
+    std::vector<double> fct_args(sis_sos.size());
+
+    for (std::size_t sink=0; sink<sis_sos.size(); ++sink)
+    {
+        auto const source = sis_sos[sink];
+
+        if (source >= 0) {
+            fct_args[sink] = call(source, unbound_arguments);
+            DBUG("setting %luth argument to %g", sink, fct_args[sink]);
+        } else {
+            assert(source != _uninitialized);
+            fct_args[sink] = unbound_arguments[-source-1];
+            DBUG("setting %luth argument to %g", sink, fct_args[sink]);
+        }
+    }
+
+    DBUG("Finished preparing call of fct #%lu %s()", function_idx,
+         _named_functions[function_idx].getName().c_str());
+
+    return _named_functions[function_idx].call(fct_args);
+}
+
+std::string NamedFunctionCaller::getCallExpression(
+    std::string const& function_name) const
+{
+    auto const fct_it = _map_name_idx.find(function_name);
+    if (fct_it == _map_name_idx.end()) {
+        OGS_FATAL("A function with the name `%s' has not been found.",
+                  function_name.c_str());
+    }
+
+    std::string expr;
+    auto callback = [&](int fct_idx, TraversePosition pos)
+    {
+        switch (pos) {
+        case TraversePosition::StartNode:
+        {
+            if (fct_idx < 0) {
+                auto it = std::find_if(
+                    _map_name_idx.begin(), _map_name_idx.end(),
+                    [fct_idx](std::pair<std::string, int> const& e) {
+                        return e.second == fct_idx;
+                    });
+                if (it == _map_name_idx.end()) {
+                    OGS_FATAL("The function index %i has not been found.", fct_idx);
+                }
+                expr += it->first;
+            } else {
+                expr += _named_functions[fct_idx].getName() + "(";
+            }
+            break;
+        }
+        case TraversePosition::BetweenChildren:
+            expr += ", ";
+            break;
+        case TraversePosition::EndNode:
+            expr += ")";
+        }
+    };
+
+    traverse(_map_sink_source, fct_it->second, callback);
+    DBUG("expression: %s", expr.c_str());
+    return expr;
+}
+
+SpecificFunctionCaller
+NamedFunctionCaller::getSpecificFunctionCaller(const std::string &function_name)
+{
+    auto const fct_it = _map_name_idx.find(function_name);
+    if (fct_it == _map_name_idx.end()) {
+        OGS_FATAL("A function with the name `%s' has not been found.",
+                  function_name.c_str());
+    }
+    return SpecificFunctionCaller(fct_it->second, *this);
+}
+
+std::size_t NamedFunctionCaller::getNumberOfUnboundArguments() const
+{
+    return -_uninitialized - 1;
+}
+
+SpecificFunctionCaller::SpecificFunctionCaller(const std::size_t function_idx,
+                                             const NamedFunctionCaller& caller)
+    : _function_idx(function_idx), _caller(caller)
+{
+}
+
+double SpecificFunctionCaller::call(
+    const std::vector<double>& unbound_arguments) const
+{
+    return _caller.call(_function_idx, unbound_arguments);
+}
+
+std::size_t SpecificFunctionCaller::getNumberOfUnboundArguments() const
+{
+    return _caller.getNumberOfUnboundArguments();
+}
+
+} // namespace NumLib
diff --git a/NumLib/NamedFunctionCaller.h b/NumLib/NamedFunctionCaller.h
new file mode 100644
index 0000000000000000000000000000000000000000..f1851b3006b280e5c3a91db59f863622af810ade
--- /dev/null
+++ b/NumLib/NamedFunctionCaller.h
@@ -0,0 +1,130 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef NUMLIB_NAMEDFUNCTIONCALLER_H
+#define NUMLIB_NAMEDFUNCTIONCALLER_H
+
+#include <map>
+#include <vector>
+#include "NamedFunction.h"
+
+namespace NumLib
+{
+class SpecificFunctionCaller;
+
+//! Builds expression trees of named functions dynamically at runtime.
+class NamedFunctionCaller final
+{
+public:
+    //! Constructs an instance whose unbound arguments have the given names.
+    explicit NamedFunctionCaller(
+        std::initializer_list<std::string> unbound_argument_names);
+
+    //! Adds the given named function
+    void addNamedFunction(NamedFunction&& fct);
+
+    //! Returns all named functions associated with the caller instance.
+    std::vector<NamedFunction> const& getNamedFunctions() const
+    {
+        return _named_functions;
+    }
+
+    //! Declares that the argument with name \c sink_arg of the function \c
+    //! sink_fct is being computed by the function \c source_fct.
+    //!
+    //! The functions involved need not already be known to the
+    //! NamedFunctionCaller.
+    void plug(std::string const& sink_fct, std::string const& sink_arg,
+              std::string const& source_fct);
+
+    //! Actually plug all the plugs previously declared.
+    //!
+    //! \pre All functions involved must have been added.
+    void applyPlugs();
+
+    //! Creates a function caller that is able to call the function with the
+    //! given name.
+    SpecificFunctionCaller getSpecificFunctionCaller(
+        std::string const& function_name);
+
+    //! Returns a string representing the expression graph of the given
+    //! function.
+    //!
+    //! \pre applyPlugs() must have been called before.
+    std::string getCallExpression(std::string const& function_name) const;
+
+    //! Returns the number of unbound arguments.
+    std::size_t getNumberOfUnboundArguments() const;
+
+private:
+    //! Calls the function with the given index with the given unbound
+    //! arguments.
+    double call(std::size_t function_idx,
+                std::vector<double> const& unbound_arguments) const;
+
+    //! Maps function names to indices.
+    //! Negative indices refer to unbound arguments.
+    std::map<std::string, int> _map_name_idx;
+
+    //! Contains all named functions.
+    std::vector<NamedFunction> _named_functions;
+
+    //! The expression graph.
+    //! Contains for each named function (outer vector) a vector which maps each
+    //! function argument to the source function index that computes this
+    //! argument.
+    std::vector<std::vector<int>> _map_sink_source;
+
+    //! Magic number used to mark function arguments in \c _map_sink_source
+    //! whose source functions have not yet been set up.
+    const int _uninitialized;
+
+    struct SinkSource
+    {
+        std::string const sink_fct;
+        std::string const sink_arg;
+        std::string const source;
+    };
+
+    //! Saves plugs declared by plug().
+    std::vector<SinkSource> _deferred_plugs;
+
+    friend class SpecificFunctionCaller;
+};
+
+//! A function caller that can call one specific function.
+//!
+//! \todo Use this class to provide some optimizations of the expression
+//! evaluation.
+class SpecificFunctionCaller final
+{
+public:
+    //! Constructs a new instance.
+    SpecificFunctionCaller(std::size_t const function_idx,
+                          NamedFunctionCaller const& caller);
+
+    //! Call the function set up with the given unbound arguments.
+    double call(std::vector<double> const& unbound_arguments) const;
+
+    //! Returns the number of unbound arguments.
+    std::size_t getNumberOfUnboundArguments() const;
+
+private:
+    //! Index of the referenced function.
+    std::size_t const _function_idx;
+
+    //! The named function caller used for the evaluation.
+    NamedFunctionCaller const& _caller;
+};
+
+} // namespace NumLib
+
+
+
+#endif // NUMLIB_NAMEDFUNCTIONCALLER_H
diff --git a/NumLib/NamedFunctionProvider.h b/NumLib/NamedFunctionProvider.h
new file mode 100644
index 0000000000000000000000000000000000000000..a814ea5772628b8d4799bf06e3187c703c92f12a
--- /dev/null
+++ b/NumLib/NamedFunctionProvider.h
@@ -0,0 +1,32 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef NUMLIB_NAMEDFUNCTIONPROVIDER_H
+#define NUMLIB_NAMEDFUNCTIONPROVIDER_H
+
+#include "NamedFunction.h"
+
+namespace NumLib
+{
+
+//! Interface used for providing named functions.
+class NamedFunctionProvider
+{
+public:
+    virtual std::vector<NamedFunction> getNamedFunctions() const
+    {
+        return std::vector<NamedFunction>{};
+    }
+
+    virtual ~NamedFunctionProvider() = default;
+};
+
+} // namespace NumLib
+
+#endif // NUMLIB_NAMEDFUNCTIONPROVIDER_H
diff --git a/Tests/BaseLib/TestFunctional.cpp b/Tests/BaseLib/TestFunctional.cpp
index 763359c1eca54bdcf36d3ca48a244e03bd17b80a..4e7da46359a355a0c0125739ec071206026a4673 100644
--- a/Tests/BaseLib/TestFunctional.cpp
+++ b/Tests/BaseLib/TestFunctional.cpp
@@ -11,52 +11,9 @@
 #include <logog/include/logog.hpp>
 
 #include "BaseLib/Functional.h"
+#include "Tests/InstanceCounter.h"
 
-class InstanceCounter
-{
-public:
-    InstanceCounter() {
-        ++_num_constructed;
-    }
-    InstanceCounter(InstanceCounter const&) {
-        ++_num_copied;
-    }
-    InstanceCounter(InstanceCounter&&) {
-        ++_num_moved;
-    }
-    virtual ~InstanceCounter() {
-        ++_num_destroyed;
-    }
-
-    static int getNumberOfConstructions() { return _num_constructed; }
-    static int getNumberOfCopies() { return _num_copied; }
-    static int getNumberOfMoves() { return _num_moved; }
-    static int getNumberOfInstances()
-    {
-        return _num_constructed + _num_moved + _num_copied - _num_destroyed;
-    }
-
-    static void update(int& num_const, int& num_move, int& num_copy, int& num_inst)
-    {
-        num_const = getNumberOfConstructions();
-        num_move = getNumberOfMoves();
-        num_copy = getNumberOfCopies();
-        num_inst = getNumberOfInstances();
-    }
-
-private:
-    static int _num_constructed;
-    static int _num_copied;
-    static int _num_moved;
-    static int _num_destroyed;
-};
-
-int InstanceCounter::_num_constructed = 0;
-int InstanceCounter::_num_copied = 0;
-int InstanceCounter::_num_moved = 0;
-int InstanceCounter::_num_destroyed = 0;
-
-class A : public InstanceCounter
+class A : private InstanceCounter<A>
 {
 public:
     A(const double value) : _value(value) {}
@@ -79,52 +36,45 @@ private:
     double _value;
 };
 
-#define EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst) \
-    EXPECT_EQ((num_const), InstanceCounter::getNumberOfConstructions());    \
-    EXPECT_EQ((num_move), InstanceCounter::getNumberOfMoves());             \
-    EXPECT_EQ((num_copy), InstanceCounter::getNumberOfCopies());            \
-    EXPECT_EQ((num_inst), InstanceCounter::getNumberOfInstances())
-
-
 TEST(BaseLib, Functional)
 {
-    auto num_const = InstanceCounter::getNumberOfConstructions();
-    auto num_move = InstanceCounter::getNumberOfMoves();
-    auto num_copy = InstanceCounter::getNumberOfCopies();
-    auto num_inst = InstanceCounter::getNumberOfInstances();
+    auto num_const = InstanceCounter<A>::getNumberOfConstructions();
+    auto num_move = InstanceCounter<A>::getNumberOfMoves();
+    auto num_copy = InstanceCounter<A>::getNumberOfCopies();
+    auto num_inst = InstanceCounter<A>::getNumberOfInstances();
 
     // Base line: measure how many copies and moves
     // std::function<>(std::bind(...)) needs.
     A a_base(0.0);
 
     // move the object to std::bind()
-    InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+    InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
     std::function<void(A)> fct_mult(
         std::bind(&A::multiply, std::move(a_base), std::placeholders::_1));
     auto const num_copy_base_move =
-        InstanceCounter::getNumberOfCopies() - num_copy;
+        InstanceCounter<A>::getNumberOfCopies() - num_copy;
     auto const num_move_base_move =
-        InstanceCounter::getNumberOfMoves() - num_move;
+        InstanceCounter<A>::getNumberOfMoves() - num_move;
 
     // call std::function using pass-by-value
     A a_base2(0.0);
-    InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+    InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
     fct_mult(a_base2);
     auto const num_copy_base_pass =
-        InstanceCounter::getNumberOfCopies() - num_copy;
+        InstanceCounter<A>::getNumberOfCopies() - num_copy;
     auto const num_move_base_pass =
-        InstanceCounter::getNumberOfMoves() - num_move;
+        InstanceCounter<A>::getNumberOfMoves() - num_move;
     // end base line
 
     // self test
-    InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+    InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
     A a1(3.0);
     A a2(a1);
-    EXPECT_INSTANCES(num_const+1, num_move, num_copy+1, num_inst+2);
-    InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+    EXPECT_INSTANCES(A, num_const+1, num_move, num_copy+1, num_inst+2);
+    InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
 
     auto f1_get = BaseLib::easyBind(&A::getValue, a1);
-    EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+    EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
     EXPECT_EQ(3.0, f1_get());
 
     // check that really a reference is returned
@@ -135,76 +85,76 @@ TEST(BaseLib, Functional)
         value_ref = 4.0;
         EXPECT_EQ(4.0, a2.getValue());
     }
-    EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+    EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
 
     // test binding to pointers
     {
         A* ap = &a1;
         auto fp_get = BaseLib::easyBind(&A::getValue, ap);
-        EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+        EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
         EXPECT_EQ(3.0, fp_get());
 
         A const* apc = &a1;
         auto fpc_get = BaseLib::easyBind(&A::getValue, apc);
-        EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+        EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
         EXPECT_EQ(3.0, fpc_get());
     }
-    EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+    EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
 
     // check that referenced objects are not copied
     {
         A& a3 = a2;
         auto f3_get = BaseLib::easyBind(&A::getValue, a3);
-        EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+        EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
         EXPECT_EQ(4.0, f3_get());
     }
     {
         A const& a3 = a2;
         auto f3_get = BaseLib::easyBind(&A::getValue, a3);
-        EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+        EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
         EXPECT_EQ(4.0, f3_get());
     }
-    EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+    EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
 
     // temporaries must be moved
     {
         auto ftemp_get = BaseLib::easyBind(&A::getValue, A(5.0));
 
-        EXPECT_INSTANCES(num_const + 1, num_move + num_move_base_move,
+        EXPECT_INSTANCES(A, num_const + 1, num_move + num_move_base_move,
                          num_copy + num_copy_base_move, num_inst + 1);
-        InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+        InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
 
         EXPECT_EQ(5.0, ftemp_get());
     }
     // ftemp_get destroyed
-    EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst-1);
-    InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+    EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst-1);
+    InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
 
     // testing explicit move
     {
         A a_move(5.0);
-        EXPECT_INSTANCES(num_const+1, num_move, num_copy, num_inst+1);
-        InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+        EXPECT_INSTANCES(A, num_const+1, num_move, num_copy, num_inst+1);
+        InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
 
         auto ftemp_get = BaseLib::easyBind(&A::getValue, std::move(a_move));
 
-        EXPECT_INSTANCES(num_const, num_move + num_move_base_move,
+        EXPECT_INSTANCES(A, num_const, num_move + num_move_base_move,
                          num_copy + num_copy_base_move, num_inst+1);
-        InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+        InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
 
         EXPECT_EQ(5.0, ftemp_get());
     }
     // ftemp_get destroyed and a_move
-    EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst-2);
-    InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+    EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst-2);
+    InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
 
     // test binding a callable object
     {
         auto f1_op = BaseLib::easyBind(a1);
-        EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+        EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
         EXPECT_EQ(21.0, f1_op(7.0));
     }
-    EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+    EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
 
     // test binding a lambda
     {
@@ -220,21 +170,21 @@ TEST(BaseLib, Functional)
     // check that parameters passed by reference are not copied
     {
         auto f1_add = BaseLib::easyBind(&A::add, a1);
-        EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+        EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
         f1_add(a2);
-        EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+        EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
         EXPECT_EQ(7.0, f1_get());
     }
 
     // check that parameters passed by value are copied
     {
         auto f1_mult = BaseLib::easyBind(&A::multiply, a1);
-        EXPECT_INSTANCES(num_const, num_move, num_copy, num_inst);
+        EXPECT_INSTANCES(A, num_const, num_move, num_copy, num_inst);
         f1_mult(a2);
 
-        EXPECT_INSTANCES(num_const, num_move + num_move_base_pass,
+        EXPECT_INSTANCES(A, num_const, num_move + num_move_base_pass,
                          num_copy + num_copy_base_pass, num_inst);
-        InstanceCounter::update(num_const, num_move, num_copy, num_inst);
+        InstanceCounter<A>::update(num_const, num_move, num_copy, num_inst);
 
         EXPECT_EQ(28.0, f1_get());
         EXPECT_EQ(4.0, a2.getValue());
diff --git a/Tests/InstanceCounter.h b/Tests/InstanceCounter.h
new file mode 100644
index 0000000000000000000000000000000000000000..4d1699459db7e7a51b2af762ba95484d9c35eba0
--- /dev/null
+++ b/Tests/InstanceCounter.h
@@ -0,0 +1,75 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#ifndef TESTS_INSTANCECOUNTER_H
+#define TESTS_INSTANCECOUNTER_H
+
+template <class T>
+class InstanceCounter
+{
+public:
+    InstanceCounter() {
+        ++_num_constructed;
+    }
+    InstanceCounter(InstanceCounter<T> const&) {
+        ++_num_copied;
+    }
+    InstanceCounter(InstanceCounter<T>&&) {
+        ++_num_moved;
+    }
+    virtual ~InstanceCounter() {
+        ++_num_destroyed;
+    }
+
+    static int getNumberOfConstructions() { return _num_constructed; }
+    static int getNumberOfCopies() { return _num_copied; }
+    static int getNumberOfMoves() { return _num_moved; }
+    static int getNumberOfDestructions() { return _num_destroyed; }
+    static int getNumberOfInstances()
+    {
+        return _num_constructed + _num_moved + _num_copied - _num_destroyed;
+    }
+
+    static void update(int& num_const, int& num_move, int& num_copy, int& num_inst)
+    {
+        num_const = getNumberOfConstructions();
+        num_move = getNumberOfMoves();
+        num_copy = getNumberOfCopies();
+        num_inst = getNumberOfInstances();
+    }
+
+private:
+    static int _num_constructed;
+    static int _num_copied;
+    static int _num_moved;
+    static int _num_destroyed;
+};
+
+template <class T>
+int InstanceCounter<T>::_num_constructed = 0;
+template <class T>
+int InstanceCounter<T>::_num_copied = 0;
+template <class T>
+int InstanceCounter<T>::_num_moved = 0;
+template <class T>
+int InstanceCounter<T>::_num_destroyed = 0;
+
+#define EXPECT_INSTANCES(type, num_const, num_move, num_copy, num_inst)        \
+    EXPECT_EQ((num_const), InstanceCounter<type>::getNumberOfConstructions()); \
+    EXPECT_EQ((num_move), InstanceCounter<type>::getNumberOfMoves());          \
+    EXPECT_EQ((num_copy), InstanceCounter<type>::getNumberOfCopies());         \
+    EXPECT_EQ((num_inst), InstanceCounter<type>::getNumberOfInstances())
+
+#define UPDATE_INSTANCES(type, num_const, num_move, num_copy, num_inst) \
+    (num_const) = InstanceCounter<type>::getNumberOfConstructions();    \
+    (num_move) = InstanceCounter<type>::getNumberOfMoves();             \
+    (num_copy) = InstanceCounter<type>::getNumberOfCopies();            \
+    (num_inst) = InstanceCounter<type>::getNumberOfInstances()
+
+#endif // TESTS_INSTANCECOUNTER_H
diff --git a/Tests/NumLib/TestNamedFunction.cpp b/Tests/NumLib/TestNamedFunction.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..bcb3b8b179707962e9a3a183e5b0a4ffeaa43af0
--- /dev/null
+++ b/Tests/NumLib/TestNamedFunction.cpp
@@ -0,0 +1,196 @@
+/**
+ * \copyright
+ * Copyright (c) 2012-2016, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ */
+
+#include <gtest/gtest.h>
+#include <logog/include/logog.hpp>
+
+#include "BaseLib/Functional.h"
+#include "NumLib/NamedFunctionProvider.h"
+#include "NumLib/NamedFunctionCaller.h"
+#include "Tests/InstanceCounter.h"
+
+class F : public NumLib::NamedFunctionProvider
+{
+public:
+    double f(double arg_g, double arg_y) const
+    {
+        return arg_g + arg_y;
+    }
+
+    std::vector<NumLib::NamedFunction>
+    getNamedFunctions() const override
+    {
+        return {{"f", {"g_arg", "y"}, BaseLib::easyBind(&F::f, this)}};
+    }
+};
+
+class G : public NumLib::NamedFunctionProvider
+{
+public:
+    double g(double arg_x) const
+    {
+        return -arg_x;
+    }
+
+    std::vector<NumLib::NamedFunction>
+    getNamedFunctions() const override
+    {
+        return {{"g", {"x"}, BaseLib::easyBind(&G::g, this)}};
+    }
+};
+
+class H : public InstanceCounter<H>
+{
+public:
+    H(double const z) : _z(z) {}
+    double h(double arg_x, double arg_y) { return arg_x * arg_y - _z; }
+    double setZ(double const z)
+    {
+        _z = z;
+        return z;
+    }
+
+private:
+    double _z;
+};
+
+class I : public NumLib::NamedFunctionProvider
+{
+public:
+    double i(double arg_x) const
+    {
+        return -arg_x;
+    }
+
+    std::vector<NumLib::NamedFunction>
+    getNamedFunctions() const override
+    {
+        return {{"i", {"x"}, BaseLib::easyBind(&I::i, this)}};
+    }
+};
+
+TEST(NumLib, NamedFunctionCaller)
+{
+    F f_inst;
+    G g_inst;
+
+    NumLib::NamedFunctionCaller caller{ "x", "y" };
+
+    for (auto&& f_named : f_inst.getNamedFunctions()) {
+        caller.addNamedFunction(std::move(f_named));
+    }
+
+    caller.plug("f", "g_arg", "g");
+    caller.plug("f", "y", "y");
+    caller.plug("g", "x", "x");
+
+    // test if adding function after plug works
+    for (auto&& g_named : g_inst.getNamedFunctions()) {
+        caller.addNamedFunction(std::move(g_named));
+    }
+
+    caller.applyPlugs();
+
+    double x = 1.0;
+    double y = 2.0;
+
+    auto const g_caller = caller.getSpecificFunctionCaller("g");
+    DBUG("calling %s", caller.getCallExpression("g").c_str());
+    EXPECT_EQ(g_inst.g(x), g_caller.call({x, y}));
+
+    auto const f_caller = caller.getSpecificFunctionCaller("f");
+    DBUG("calling %s", caller.getCallExpression("f").c_str());
+    EXPECT_EQ(f_inst.f(g_inst.g(x), y), f_caller.call({x, y}));
+}
+
+TEST(NumLib, NamedFunctionCallerCyclicGraph)
+{
+    // Construct a cyclic case with f(g(i(f(...), y)))
+    F f_inst;
+    G g_inst;
+    I i_inst;
+
+    NumLib::NamedFunctionCaller caller{ "x", "y" };
+
+    for (auto&& f_named : f_inst.getNamedFunctions()) {
+        caller.addNamedFunction(std::move(f_named));
+    }
+    for (auto&& g_named : g_inst.getNamedFunctions()) {
+        caller.addNamedFunction(std::move(g_named));
+    }
+    for (auto&& i_named : i_inst.getNamedFunctions()) {
+        caller.addNamedFunction(std::move(i_named));
+    }
+
+    caller.plug("f", "g_arg", "g");
+    caller.plug("f", "y", "y");
+    caller.plug("g", "x", "i");
+    caller.plug("i", "x", "f");
+
+    ASSERT_ANY_THROW(caller.applyPlugs());
+}
+
+TEST(NumLib, NamedFunctionNoLeaks)
+{
+    auto num_const = InstanceCounter<H>::getNumberOfConstructions();
+    auto num_move = InstanceCounter<H>::getNumberOfMoves();
+    auto num_copy = InstanceCounter<H>::getNumberOfCopies();
+    auto num_inst = InstanceCounter<H>::getNumberOfInstances();
+
+    {
+        H h_inst(1.0);
+        EXPECT_EQ(num_const+1, InstanceCounter<H>::getNumberOfConstructions());
+        EXPECT_EQ(num_inst+1, InstanceCounter<H>::getNumberOfInstances());
+        InstanceCounter<H>::update(num_const, num_move, num_copy, num_inst);
+
+        auto h_fct = NumLib::NamedFunction("h", {"x", "y"},
+                                           BaseLib::easyBind(&H::h, h_inst));
+        EXPECT_EQ(num_const, InstanceCounter<H>::getNumberOfConstructions());
+        EXPECT_EQ(num_inst, InstanceCounter<H>::getNumberOfInstances());
+
+        EXPECT_EQ(5.0, h_fct.call({ 2.0, 3.0 }));
+        h_inst.setZ(2.0);
+        EXPECT_EQ(4.0, h_fct.call({ 2.0, 3.0 }));
+
+        auto h_bind = BaseLib::easyBind(&H::h, H{3.0});
+        EXPECT_EQ(num_const+1, InstanceCounter<H>::getNumberOfConstructions());
+        EXPECT_EQ(num_inst+1, InstanceCounter<H>::getNumberOfInstances());
+        InstanceCounter<H>::update(num_const, num_move, num_copy, num_inst);
+
+        // Move an object. NamedFunction will implicitly do memory management.
+        auto h_fct2 = NumLib::NamedFunction("h", {"x", "y"}, std::move(h_bind));
+        EXPECT_EQ(num_copy, InstanceCounter<H>::getNumberOfCopies());
+        EXPECT_EQ(num_const, InstanceCounter<H>::getNumberOfConstructions());
+        EXPECT_EQ(num_inst, InstanceCounter<H>::getNumberOfInstances());
+
+        EXPECT_EQ(5.0, h_fct2.call({ 2.0, 4.0 }));
+
+        // copy
+        NumLib::NamedFunction h_fct3 = h_fct2;
+        EXPECT_EQ(num_const, InstanceCounter<H>::getNumberOfConstructions());
+        EXPECT_EQ(num_copy+1, InstanceCounter<H>::getNumberOfCopies());
+        EXPECT_EQ(num_inst+1, InstanceCounter<H>::getNumberOfInstances());
+        InstanceCounter<H>::update(num_const, num_move, num_copy, num_inst);
+
+        EXPECT_EQ(-1.0, h_fct3.call({ 2.0, 1.0 }));
+
+        // move
+        NumLib::NamedFunction h_fct4 = std::move(h_fct3);
+        EXPECT_INSTANCES(H, num_const, num_move, num_copy, num_inst);
+        InstanceCounter<H>::update(num_const, num_move, num_copy, num_inst);
+
+        EXPECT_EQ(3.0, h_fct4.call({ 3.0, 2.0 }));
+    }
+    EXPECT_EQ(num_const, InstanceCounter<H>::getNumberOfConstructions());
+    EXPECT_EQ(num_copy, InstanceCounter<H>::getNumberOfCopies());
+    // If zero instances are left, the destructor has been called the right
+    // number of times, i.e., all internal casts in NamedFunction have been
+    // successful.
+    EXPECT_EQ(0, InstanceCounter<H>::getNumberOfInstances());
+}