// -*- tab-width: 4; indent-tabs-mode: nil -*-

/** \brief A function that defines analytic solution as Dirichlet boundary conditions AND
    its extension to the interior */
template<typename T>
class BCAnalytic
  : public Dune::PDELab::GridFunctionBase<Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
									   typename T::Traits::RangeFieldType,
									   1,
									   Dune::FieldVector<typename T::Traits::RangeFieldType,1> >,
					  BCAnalytic<T> > {
public:
 typedef Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
					  typename T::Traits::RangeFieldType,
					  1,
					  Dune::FieldVector<typename T::Traits::RangeFieldType,1> > Traits;

  //! construct from grid view
  BCAnalytic  (const typename Traits::GridViewType& gv_, T& t_) : gv(gv_), t(t_) {}

  //! evaluate extended function on element
  inline void evaluate (const typename Traits::ElementType& e,
                        const typename Traits::DomainType& xlocal,
                        typename Traits::RangeType& y) const
  {
    y = t.ua(e, xlocal);
  }

  //! get a reference to the grid view
  inline const typename Traits::GridViewType& getGridView () const
  {
    return gv;
  }

private:
  const typename Traits::GridViewType gv;
  T& t;
};


/** \brief A function that defines gradient of analytic solution as Dirichlet boundary conditions AND
    its extension to the interior */
template<typename T>
class BCAnalyticGrad
  : public Dune::PDELab::GridFunctionBase<Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
									   typename T::Traits::RangeFieldType,
									   1,
									   Dune::FieldVector<typename T::Traits::RangeFieldType,1> >,
					  BCAnalyticGrad<T> > {
public:
 typedef Dune::PDELab::GridFunctionTraits<typename T::Traits::GridViewType,
					  typename T::Traits::RangeFieldType,
					  1,
					  Dune::FieldVector<typename T::Traits::RangeFieldType,1> > Traits;

  //! construct from grid view
  BCAnalyticGrad  (const typename Traits::GridViewType& gv_, T& t_) : gv(gv_), t(t_) {}

  const static int dim = T::Traits::GridViewType::Grid::dimension;
  // typedef typename T::Traits::GridViewType::Grid::ctype ctype;

  //! evaluate extended function on element
  inline void evaluate (const typename Traits::ElementType& e,
                        const typename Traits::DomainType& xlocal,
                        Dune::FieldVector<typename T::Traits::RangeFieldType,dim>& y) const
  {

    y = t.ugrada(e, xlocal);
  }

  //! get a reference to the grid view
  inline const typename Traits::GridViewType& getGridView () const
  {
    return gv;
  }

private:
  const typename Traits::GridViewType gv;
  T& t;

};




/*! \brief Adapter returning ||f1(x)-f2(x)||^2 for two given grid functions

  \tparam T1  a grid function type
  \tparam T2  a grid function type
*/
template<typename T1, typename T2>
class DifferenceSquaredAdapter
  : public Dune::PDELab::GridFunctionBase<
  Dune::PDELab::GridFunctionTraits<typename T1::Traits::GridViewType,
				   typename T1::Traits::RangeFieldType,
				   1,Dune::FieldVector<typename T1::Traits::RangeFieldType,1> >
  ,DifferenceSquaredAdapter<T1,T2> >
{
public:
  typedef Dune::PDELab::GridFunctionTraits<typename T1::Traits::GridViewType,
                                           typename T1::Traits::RangeFieldType,
                                           1,Dune::FieldVector<typename T1::Traits::RangeFieldType,1> > Traits;

  //! constructor
  DifferenceSquaredAdapter (const T1& t1_, const T2& t2_) : t1(t1_), t2(t2_) {}

  //! \copydoc GridFunctionBase::evaluate()
  inline void evaluate (const typename Traits::ElementType& e,
                        const typename Traits::DomainType& x,
                        typename Traits::RangeType& y) const
  {
    typename T1::Traits::RangeType v1,v2;
    t1.evaluate( e, x, v1);
    t2.evaluate( e, x, v2);
    v1 -= v2;
    y = v1*v1;
  }

  inline const typename Traits::GridViewType& getGridView () const
  {
    return t1.getGridView();
  }

private:
  const T1& t1;
  const T2& t2;
};