#ifndef __LocalContributions_hpp__
#define __LocalContributions_hpp__

#include "CommonTemplateTypes.hpp"
#include "InputData.hpp"

namespace ho
{
  struct LocalContributions
  {
    typedef CommonTemplateTypes::EigenGlobalVectorType EigenGlobalVectorType;
    typedef CommonTemplateTypes::ExactSolutionType ExactSolutionType;
    typedef CommonTemplateTypes::DiffusivityValueType DiffusivityValueType;
    typedef CommonTemplateTypes::DiffusivityType DiffusivityType;
    typedef CommonTemplateTypes::LoadType LoadType;
    typedef CommonTemplateTypes::MeshType MeshType;
    typedef CommonTemplateTypes::FEManagerType FEManagerType; 
  
    LocalContributions(){;}
    
    // Initialize and compute local operators
    void assemble(const MeshType& mesh,
                  const SpaMesh::ElementIterator& elemIt,
      		  const unsigned int p_order,
                  const ExactSolutionType & exact_solution,
                  const DiffusivityType & diffusivity,
                  const LoadType & load,
                  const bool & weak_bc = false,
                  const double & eta = 1.);
    // Interpolate exact solution
    Eigen::VectorXd interpolate(const MeshType& mesh, 
      		                const unsigned int p_order,
	                        const SpaMesh::ElementIterator& elemIt, 
				const ExactSolutionType & u);
    // Reconstruct local solution
    Eigen::VectorXd reconstruct(const Eigen::VectorXd & uFh);
                                
    // Reconstruct local solution
    Eigen::VectorXd reconstruct(const EigenGlobalVectorType & uFh,
                                const Eigen::VectorXi& idx_T);
    // High-order potential reconstruction
    Eigen::Matrix<double, InputData::n_dofs_ppo, 1>
    potentialReconstruction(const Eigen::VectorXd & uh_TF);
  
    enum {
      nb_local_face_dofs = InputData::n_face_dofs,
      nb_cell_dofs       = InputData::n_dofs,
      nb_local_cell_dofs = InputData::n_dofs,
      NG                 = InputData::n_dofs_ppo - 1,
      NP                 = InputData::n_dofs_ppo
    };
  
    FEManagerType FEManager;
  
    Eigen::MatrixXd STF; // Stabilization
    Eigen::MatrixXd JTF; // Penalty for weakly enforced BCs
  
    Eigen::MatrixXd ATF;
    Eigen::VectorXd bTF;
    Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic> AF;
    Eigen::VectorXd bF;
  
    Eigen::Matrix<double, InputData::n_dofs_ppo, InputData::n_dofs_ppo> MTT;
    std::vector<Eigen::Matrix<double, InputData::n_face_dofs, InputData::n_face_dofs> > MFF;
  
    Eigen::MatrixXd GT;
  
    std::size_t nb_face_dofs;
    std::size_t nb_tot_dofs;
  };
}

#endif
