#pragma once

#include <type_traits>

#include <dune/functions/functionspacebases/subentitydofs.hh>

#include <amdis/functions/Interpolate.hpp>
#include <amdis/LinearAlgebra.hpp>

namespace AMDiS {
namespace Impl {

/// Traverse all interior boundary DOFs and apply the functor f.
/// The Functor must have the signature
/// `void(int, typename Basis::LocalView, typename Basis::GridView::Intersection)`
template <class Basis, class F>
void forEachInteriorBoundaryDOF(Basis const& basis, F&& f)
{
  auto localView = basis.localView();
  auto seDOFs = Dune::Functions::subEntityDOFs(basis);
  auto const& gridView = basis.gridView();
  for (auto const& element : elements(gridView, typename BackendTraits<Basis>::PartitionSet{})) {
    if (element.hasBoundaryIntersections()) {
      localView.bind(element);
      for(auto const& intersection: intersections(gridView, element))
        if (intersection.boundary())
          for(auto localIndex: seDOFs.bind(localView,intersection))
            f(localIndex, localView, intersection);
    }
  }
}

} // end namespace Impl


template <class Mat, class Sol, class Rhs, class RB, class CB, class V>
void DirichletBC<Mat, Sol, Rhs, RB, CB, V>::
init()
{
  using LV = typename CB::LocalView;
  dirichletNodes_.assign(colBasis_.dimension(), false);
  Impl::forEachInteriorBoundaryDOF(colBasis_,
    [&](int localIndex, LV const& localView, Intersection const& intersection) {
      dirichletNodes_[localView.index(localIndex)] = dirichletNodes_[localView.index(localIndex)]
                                                     || boundarySubset_(intersection);
  });
}


template <class Mat, class Sol, class Rhs, class RB, class CB, class V>
void DirichletBC<Mat, Sol, Rhs, RB, CB, V>::
apply(Mat& matrix, Sol& solution, Rhs& rhs)
{
  AMDiS::interpolate(colBasis_, solution, values_, tag::defaulted{}, tag::defaulted{}, dirichletNodes_);
  dirichletBC(matrix, solution, rhs, dirichletNodes_);
}

} // end namespace AMDiS