evaluator.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/trafo/evaluator_base.hpp>
#include <kernel/trafo/standard/details.hpp>
 
namespace FEAT
{
  namespace Trafo
  {
    namespace Standard
    {
      template<
        typename Trafo_,
        typename EvalPolicy_,
        typename Shape_ = typename EvalPolicy_::ShapeType>
      class Evaluator:
        public EvaluatorBase<Trafo_, Evaluator<Trafo_, EvalPolicy_, Shape_ >, EvalPolicy_>
      {
      public:
        typedef EvaluatorBase<Trafo_, Evaluator, EvalPolicy_> BaseClass;
        typedef Shape_ ShapeType;
        typedef Trafo_ TrafoType;
        typedef EvalPolicy_ EvalPolicy;
 
        typedef typename TrafoType::MeshType MeshType;
 
        typedef typename EvalPolicy::DataType DataType;
        typedef typename EvalPolicy::DomainPointType DomainPointType;
        typedef typename EvalPolicy::ImagePointType ImagePointType;
        typedef typename EvalPolicy::JacobianMatrixType JacobianMatrixType;
        typedef typename EvalPolicy::JacobianInverseType JacobianInverseType;
        typedef typename EvalPolicy::JacobianDeterminantType JacobianDeterminantType;
        typedef typename EvalPolicy::HessianTensorType HessianTensorType;
        typedef typename EvalPolicy::HessianInverseType HessianInverseType;
 
        static constexpr int domain_dim = EvalPolicy::domain_dim;
        static constexpr int image_dim = EvalPolicy::image_dim;
        static constexpr int num_verts = Shape::FaceTraits<ShapeType, 0>::count;
 
        typedef EvalHelper<DataType, DomainPointType, ImagePointType, JacobianMatrixType, JacobianInverseType, JacobianDeterminantType,
                          HessianTensorType, HessianInverseType ,ShapeType, image_dim> EvalHp;
 
        static constexpr TrafoTags eval_caps =
          TrafoTags::dom_point | TrafoTags::img_point |
          TrafoTags::jac_mat | TrafoTags::jac_det | TrafoTags::hess_ten |
          (domain_dim == image_dim ? (TrafoTags::jac_inv | TrafoTags::hess_inv) : TrafoTags::none);
 
      protected:
        Tiny::Matrix<DataType, image_dim, num_verts> _coeff;
 
      public:
        explicit Evaluator(const TrafoType& trafo) :
          BaseClass(trafo)
        {
        }
 
        void prepare(Index cell_index)
        {
          // prepare base-class
          BaseClass::prepare(cell_index);
 
          // fetch the mesh from the trafo
          const MeshType& mesh = this->_trafo.get_mesh();
 
          // fetch the vertex set from the mesh
          typedef typename MeshType::VertexSetType VertexSetType;
          const VertexSetType& vertex_set = mesh.get_vertex_set();
 
          // fetch the index set
          typedef typename MeshType::template IndexSet<domain_dim, 0>::Type IndexSetType;
          const IndexSetType& index_set = mesh.template get_index_set<domain_dim, 0>();
 
          EvalHp::set_coefficients(_coeff, vertex_set, index_set, cell_index);
        }
 
        void map_point(ImagePointType& img_point, const DomainPointType& dom_point) const
        {
          EvalHp::map_point(img_point, dom_point, _coeff);
        }
 
        void calc_jac_mat(JacobianMatrixType& jac_mat, const DomainPointType& dom_point) const
        {
          EvalHp::calc_jac_mat(jac_mat, dom_point, _coeff);
        }
 
        void calc_hess_ten(HessianTensorType& hess_ten, const DomainPointType& dom_point) const
        {
          EvalHp::calc_hess_ten(hess_ten, dom_point, _coeff);
        }
 
        DataType volume() const
        {
          return EvalHp::volume(_coeff);
        }
 
        DataType width_directed(const ImagePointType& ray) const
        {
          return EvalHp::width_directed(ray, _coeff);
        }
      }; // class Evaluator<Hypercube<1>,...>
 
      // /* ************************************************************************************* */
 
      template<
        typename Trafo_,
        typename EvalPolicy_>
      class Evaluator<Trafo_, EvalPolicy_, Shape::Vertex > :
        public EvaluatorBase<Trafo_, Evaluator<Trafo_, EvalPolicy_, Shape::Vertex >, EvalPolicy_>
      {
      public:
        typedef EvaluatorBase<Trafo_, Evaluator, EvalPolicy_> BaseClass;
        typedef Shape::Vertex ShapeType;
        typedef Trafo_ TrafoType;
        typedef EvalPolicy_ EvalPolicy;
 
        typedef typename TrafoType::MeshType MeshType;
 
        typedef typename EvalPolicy::DataType DataType;
        typedef typename EvalPolicy::DomainPointType DomainPointType;
        typedef typename EvalPolicy::ImagePointType ImagePointType;
        typedef typename EvalPolicy::JacobianMatrixType JacobianMatrixType;
        typedef typename EvalPolicy::JacobianInverseType JacobianInverseType;
        typedef typename EvalPolicy::JacobianDeterminantType JacobianDeterminantType;
        typedef typename EvalPolicy::HessianTensorType HessianTensorType;
        typedef typename EvalPolicy::HessianInverseType HessianInverseType;
 
        static constexpr int domain_dim = EvalPolicy::domain_dim;
        static constexpr int image_dim = EvalPolicy::image_dim;
 
        static constexpr TrafoTags eval_caps = TrafoTags::dom_point | TrafoTags::img_point;
 
      protected:
        Tiny::Vector<DataType, image_dim> _coeff;
 
      public:
        explicit Evaluator(const TrafoType& trafo) :
          BaseClass(trafo)
        {
        }
 
        void prepare(Index cell_index)
        {
          // prepare base-class
          BaseClass::prepare(cell_index);
 
          // fetch the mesh from the trafo
          const MeshType& mesh = this->_trafo.get_mesh();
 
          // fetch the vertex set from the mesh
          typedef typename MeshType::VertexSetType VertexSetType;
          const VertexSetType& vertex_set = mesh.get_vertex_set();
 
          // fetch the vertex
          typedef typename VertexSetType::VertexType VertexType;
          const VertexType& vtx = vertex_set[cell_index];
 
          // calculate transformation coefficients
          for(int i(0); i < image_dim; ++i)
          {
            _coeff[i] = DataType(vtx[i]);
          }
        }
 
        void map_point(ImagePointType& img_point, const DomainPointType& DOXY(dom_point)) const
        {
          for(int i(0); i < image_dim; ++i)
          {
            img_point[i] = _coeff[i];
          }
        }
      }; // class Evaluator<Vertex,...>
    } // namespace Standard
  } // namespace Trafo
} // namespace FEAT