function.hpp

// FEAT3: Finite Element Analysis Toolbox, Version 3
// Copyright (C) 2010 by Stefan Turek & the FEAT group
// FEAT3 is released under the GNU General Public License version 3,
// see the file 'copyright.txt' in the top level directory for details.
 
#pragma once
 
// includes, FEAT
#include <kernel/base_header.hpp>
#include <kernel/util/tiny_algebra.hpp>
 
namespace FEAT
{
  namespace Analytic
  {
    namespace Image
    {
      struct Scalar
      {
        static constexpr bool is_scalar = true;
        static constexpr bool is_vector = false;
        static constexpr int scalar_components = 1;
      };
 
      template<int m_>
      struct Vector
      {
        static_assert(m_ > 0, "invalid vector field dimension");
 
        static constexpr bool is_scalar = false;
        static constexpr bool is_vector = true;
        static constexpr int scalar_components = m_;
 
        static constexpr int image_dim = m_;
      };
    } // namespace Image
 
    template<typename DataType_, int domain_dim_, typename ImageType_>
    struct EvalTraitsBase;
 
    template<typename DataType_, int domain_dim_>
    struct EvalTraitsBase<DataType_, domain_dim_, Image::Scalar>
    {
      typedef DataType_ DataType;
      static constexpr int domain_dim = domain_dim_;
      static constexpr int image_dim = 1;
 
      typedef Tiny::Vector<DataType_, domain_dim_> PointType;
 
      typedef DataType_ ValueType;
      typedef Tiny::Vector<DataType_, domain_dim_> GradientType;
      typedef Tiny::Matrix<DataType_, domain_dim_, domain_dim_> HessianType;
    };
 
    template<typename DataType_, int domain_dim_, int image_dim_>
    struct EvalTraitsBase<DataType_, domain_dim_, Image::Vector<image_dim_>>
    {
      typedef DataType_ DataType;
      static constexpr int domain_dim = domain_dim_;
      static constexpr int image_dim = image_dim_;
 
      typedef Tiny::Vector<DataType_, domain_dim_> PointType;
 
      typedef Tiny::Vector<DataType_, image_dim_> ValueType;
      typedef Tiny::Matrix<DataType_, image_dim_, domain_dim_> GradientType;
      typedef Tiny::Tensor3<DataType_, image_dim_, domain_dim_, domain_dim_> HessianType;
    };
 
    template<typename DataType_, typename Function_>
    struct EvalTraits :
      public EvalTraitsBase<DataType_, Function_::domain_dim, typename Function_::ImageType>
    {
    };
 
    class Function
    {
    public:
      // Note:
      // The following block is only visible for doxygen. Its contents must be specified
      // by the derived classes.
#ifdef DOXYGEN
      static constexpr int domain_dim = ...;
 
      typedef ... ImageType;
 
      static constexpr bool can_value = true or false;
      static constexpr bool can_grad = true or false;
      static constexpr bool can_hess = true or false;
#else
      static constexpr bool can_value = false;
      static constexpr bool can_grad = false;
      static constexpr bool can_hess = false;
#endif // DOXYGEN
 
      template<typename Traits_>
      class Evaluator
      {
      public:
        typedef typename Traits_::DataType DataType;
 
        typedef typename Traits_::PointType PointType;
 
        typedef typename Traits_::ValueType ValueType;
 
        typedef typename Traits_::GradientType GradientType;
 
        typedef typename Traits_::HessianType HessianType;
 
      public:
#ifdef DOXYGEN
        explicit Evaluator(const Function& function)
        {
        }
 
        ValueType value(const PointType& point)
        {
        }
 
        GradientType gradient(const PointType& point)
        {
        }
 
        HessianType hessian(const PointType& point)
        {
        }
#endif // DOXYGEN
      }; // class Function::Evaluator<...>
    }; // class Function
 
    template<typename Function_, typename DT_, int dim_, int s_>
    typename Analytic::EvalTraits<DT_, Function_>::ValueType
    eval_value(const Function_& function, const Tiny::Vector<DT_, dim_, s_>& point)
    {
      typedef Analytic::EvalTraits<DT_, Function_> Traits;
      static_assert(Function_::can_value, "function does not support evaluation of values");
      static_assert(dim_ == Traits::domain_dim, "invalid point dimension");
      typename Function_::template Evaluator<Traits> evaluator(function);
      return evaluator.value(point);
    }
 
    // \see eval_value()
    template<typename Function_, typename DT_>
    typename Analytic::EvalTraits<DT_, Function_>::ValueType
    eval_value_x(const Function_& function, const DT_ x)
    {
      Tiny::Vector<DT_, 1> p;
      p.v[0] = x;
      return eval_value(function, p);
    }
 
    // \see eval_value()
    template<typename Function_, typename DT_>
    typename Analytic::EvalTraits<DT_, Function_>::ValueType
    eval_value_x(const Function_& function, const DT_ x, const DT_ y)
    {
      Tiny::Vector<DT_, 2> p;
      p.v[0] = x;
      p.v[1] = y;
      return eval_value(function, p);
    }
 
    // \see eval_value()
    template<typename Function_, typename DT_>
    typename Analytic::EvalTraits<DT_, Function_>::ValueType
    eval_value_x(const Function_& function, const DT_ x, const DT_ y, const DT_ z)
    {
      Tiny::Vector<DT_, 3> p;
      p.v[0] = x;
      p.v[1] = y;
      p.v[2] = z;
      return eval_value(function, p);
    }
 
    template<typename Function_, typename DT_, int dim_, int s_>
    typename Analytic::EvalTraits<DT_, Function_>::GradientType
    eval_gradient(const Function_& function, const Tiny::Vector<DT_, dim_, s_>& point)
    {
      typedef Analytic::EvalTraits<DT_, Function_> Traits;
      static_assert(Function_::can_grad, "function does not support evaluation of gradients");
      static_assert(dim_ == Traits::domain_dim, "invalid point dimension");
      typename Function_::template Evaluator<Traits> evaluator(function);
      return evaluator.gradient(point);
    }
 
    // \see eval_gradient()
    template<typename Function_, typename DT_>
    typename Analytic::EvalTraits<DT_, Function_>::GradientType
    eval_gradient_x(const Function_& function, const DT_ x)
    {
      Tiny::Vector<DT_, 1> p;
      p.v[0] = x;
      return eval_gradient(function, p);
    }
 
    // \see eval_gradient()
    template<typename Function_, typename DT_>
    typename Analytic::EvalTraits<DT_, Function_>::GradientType
    eval_gradient_x(const Function_& function, const DT_ x, const DT_ y)
    {
      Tiny::Vector<DT_, 2> p;
      p.v[0] = x;
      p.v[1] = y;
      return eval_gradient(function, p);
    }
 
    // \see eval_gradient()
    template<typename Function_, typename DT_>
    typename Analytic::EvalTraits<DT_, Function_>::GradientType
    eval_gradient_x(const Function_& function, const DT_ x, const DT_ y, const DT_ z)
    {
      Tiny::Vector<DT_, 3> p;
      p.v[0] = x;
      p.v[1] = y;
      p.v[2] = z;
      return eval_gradient(function, p);
    }
 
    template<typename Function_, typename DT_, int dim_, int s_>
    typename Analytic::EvalTraits<DT_, Function_>::HessianType
    eval_hessian(const Function_& function, const Tiny::Vector<DT_, dim_, s_>& point)
    {
      typedef Analytic::EvalTraits<DT_, Function_> Traits;
      static_assert(Function_::can_hess, "function does not support evaluation of hessians");
      static_assert(dim_ == Traits::domain_dim, "invalid point dimension");
      typename Function_::template Evaluator<Traits> evaluator(function);
      return evaluator.hessian(point);
    }
 
    // \see eval_hessian()
    template<typename Function_, typename DT_>
    typename Analytic::EvalTraits<DT_, Function_>::HessianType
    eval_hessian_x(const Function_& function, const DT_ x)
    {
      Tiny::Vector<DT_, 1> p;
      p.v[0] = x;
      return eval_hessian(function, p);
    }
 
    // \see eval_hessian()
    template<typename Function_, typename DT_>
    typename Analytic::EvalTraits<DT_, Function_>::HessianType
    eval_hessian_x(const Function_& function, const DT_ x, const DT_ y)
    {
      Tiny::Vector<DT_, 2> p;
      p.v[0] = x;
      p.v[1] = y;
      return eval_hessian(function, p);
    }
 
    // \see eval_hessian()
    template<typename Function_, typename DT_>
    typename Analytic::EvalTraits<DT_, Function_>::HessianType
    eval_hessian_x(const Function_& function, const DT_ x, const DT_ y, const DT_ z)
    {
      Tiny::Vector<DT_, 3> p;
      p.v[0] = x;
      p.v[1] = y;
      p.v[2] = z;
      return eval_hessian(function, p);
    }
 
  } // namespace Analytic
} // namespace FEAT