NodeCache.hpp

#pragma once
 
#include <memory>
#include <tuple>
#include <vector>
 
#include <dune/common/ftraits.hh>
#include <dune/common/hash.hh>
#include <dune/common/indices.hh>
 
#include <dune/geometry/type.hh>
 
#include <dune/typetree/leafnode.hh>
#include <dune/typetree/powernode.hh>
#include <dune/typetree/compositenode.hh>
 
#include <amdis/common/ConcurrentCache.hpp>
 
namespace AMDiS
{
  namespace Impl
  {
    template <class Node, class NodeTag = typename Node::NodeTag>
    struct NodeCacheFactory;
 
  } // end namespace Impl
 
 
  template <class Node>
  using NodeCache_t = typename Impl::NodeCacheFactory<Node>::type;
 
  template <class Node>
  auto makeNodeCache (Node const& node)
  {
    return NodeCache_t<Node>::create(node);
  }
 
 
  template <class NodeType>
  class NodeWrapper
  {
  public:
    using Node = NodeType;
    using Element = typename Node::Element;
 
  public:
    NodeWrapper (Node const& node)
      : node_(&node)
    {}
 
    Node const& node () const
    {
      assert(node_ != nullptr);
      return *node_;
    }
 
    Element const& element () const
    {
      return node_->element();
    }
 
    auto localIndex (std::size_t i) const
    {
      return node_->localIndex(i);
    }
 
    auto size () const
    {
      return node_->size();
    }
 
    auto treeIndex () const
    {
      return node_->treeIndex();
    }
 
  protected:
    Node const* node_ = nullptr;
  };
 
 
 
  template <class Node>
  class LeafNodeCache
      : public Dune::TypeTree::LeafNode
      , public NodeWrapper<Node>
  {
  public:
    using BasisNode = Node;
 
    using FiniteElement = typename Node::FiniteElement;
    using LocalBasis = typename FiniteElement::Traits::LocalBasisType;
 
    using ShapeValues = std::vector<typename LocalBasis::Traits::RangeType>;
    using ShapeGradients = std::vector<typename LocalBasis::Traits::JacobianType>;
 
  private:
    using DomainType = typename LocalBasis::Traits::DomainType;
 
    // Pair of GeometryType and local coordinates
    struct CoordKey
    {
      unsigned int id;    // topologyId
      DomainType local;   // local coordinate
 
      struct hasher
      {
        std::size_t operator()(CoordKey const& t) const
        {
          std::size_t seed = 0;
          Dune::hash_combine(seed, t.id);
          Dune::hash_range(seed, t.local.begin(), t.local.end());
          return seed;
        }
      };
 
      friend bool operator==(CoordKey const& lhs, CoordKey const& rhs)
      {
        return std::tie(lhs.id,lhs.local) == std::tie(rhs.id,rhs.local);
      }
    };
 
  private:
    // Constructor storing a reference to the passed basis-node
    LeafNodeCache (Node const& basisNode)
      : NodeWrapper<Node>(basisNode)
    {}
 
  public:
    static LeafNodeCache create (Node const& basisNode)
    {
      return {basisNode};
    }
 
    FiniteElement const& finiteElement () const
    {
      return this->node_->finiteElement();
    }
 
    ShapeValues const& localBasisValuesAt (DomainType const& local) const
    {
      CoordKey key{this->element().type().id(), local};
      return shapeValues_.get(key, [&](CoordKey const&)
      {
        ShapeValues data;
        this->localBasis().evaluateFunction(local, data);
        return data;
      });
    }
 
    ShapeGradients const& localBasisJacobiansAt (DomainType const& local) const
    {
      CoordKey key{this->element().type().id(), local};
      return shapeGradients_.get(key, [&](CoordKey const&)
      {
        ShapeGradients data;
        this->localBasis().evaluateJacobian(local, data);
        return data;
      });
    }
 
  private:
    LocalBasis const& localBasis () const
    {
      return this->node_->finiteElement().localBasis();
    }
 
  private:
    template <class Value>
    using CoordCache = ConcurrentCache<CoordKey, Value, ConsecutivePolicy,
      std::unordered_map<CoordKey, Value, typename CoordKey::hasher>>;
 
    CoordCache<ShapeValues> shapeValues_;
    CoordCache<ShapeGradients> shapeGradients_;
  };
 
 
  template <class Node>
  class PowerNodeCache
      : public Impl::NodeCacheFactory<Node>::Base
      , public NodeWrapper<Node>
  {
    using Self = PowerNodeCache;
    using Super = typename Impl::NodeCacheFactory<Node>::Base;
 
  private:
    PowerNodeCache (Node const& basisNode)
      : NodeWrapper<Node>(basisNode)
    {}
 
  public:
    static auto create (Node const& basisNode)
    {
      Self cache{basisNode};
      for (std::size_t i = 0; i < std::size_t(basisNode.degree()); ++i)
        cache.setChild(i, Super::ChildType::create(basisNode.child(i)));
      return cache;
    }
  };
 
 
  template <class Node>
  class CompositeNodeCache
      : public Impl::NodeCacheFactory<Node>::Base
      , public NodeWrapper<Node>
  {
    using Self = CompositeNodeCache;
    using Super = typename Impl::NodeCacheFactory<Node>::Base;
 
  private:
    CompositeNodeCache (Node const& basisNode)
      : NodeWrapper<Node>(basisNode)
    {}
 
  public:
    static auto create (Node const& basisNode)
    {
      using TT = typename Super::ChildTypes;
      Self cache{basisNode};
      Dune::Hybrid::forEach(std::make_index_sequence<std::size_t(Node::degree())>{}, [&](auto ii) {
        cache.setChild(std::tuple_element_t<ii,TT>::create(basisNode.child(ii)), ii);
      });
      return cache;
    }
  };
 
  namespace Impl
  {
    template <class Node>
    struct NodeCacheFactory<Node, Dune::TypeTree::LeafNodeTag>
    {
      using Base = Dune::TypeTree::LeafNode;
      using type = LeafNodeCache<Node>;
    };
 
    template <class Node>
    struct NodeCacheFactory<Node, Dune::TypeTree::PowerNodeTag>
    {
      using Child = typename NodeCacheFactory<typename Node::ChildType>::type;
      using Base = Dune::TypeTree::PowerNode<Child, std::size_t(Node::degree())>;
      using type = PowerNodeCache<Node>;
    };
 
    template <class Node>
    struct NodeCacheFactory<Node, Dune::TypeTree::CompositeNodeTag>
    {
      template <class Indices> struct Childs;
      template <std::size_t... i>
      struct Childs<std::index_sequence<i...>> {
        using type = Dune::TypeTree::CompositeNode<
          typename NodeCacheFactory<typename Node::template Child<i>::type>::type...
          >;
      };
 
      using Base = typename Childs<std::make_index_sequence<std::size_t(Node::degree())>>::type;
      using type = CompositeNodeCache<Node>;
    };
 
  } // end namespace Impl
 
} // end namespace AMDiS