diff --git a/dune/fufem/discretizationerror.hh b/dune/fufem/discretizationerror.hh index d21b4d095a1383aec7c9d00d75e78514dc8bba34..0c3b7555cf670a9fac4de0ac7a313725ffd62065 100644 --- a/dune/fufem/discretizationerror.hh +++ b/dune/fufem/discretizationerror.hh @@ -19,24 +19,11 @@ public: /** \brief Compute L2 error between a grid function and an arbitrary function */ - static double computeL2Error(const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1> >* a, - const Dune::VirtualFunction<Dune::FieldVector<ctype,dim>, Dune::FieldVector<double,1> >* b, - QuadratureRuleKey quadKey); - - static double computeH1HalfNormDifferenceSquared(const GridView& gridView, - const VirtualDifferentiableFunction<Dune::FieldVector<ctype,dim>, Dune::FieldVector<double,1> >* u, - const VirtualDifferentiableFunction<Dune::FieldVector<ctype,dim>, Dune::FieldVector<double,1> >* v, - QuadratureRuleKey quadKey); - - -}; - -template <class GridView> -double DiscretizationError<GridView>:: -computeL2Error(const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1> >* a, - const Dune::VirtualFunction<Dune::FieldVector<ctype,dim>, Dune::FieldVector<double,1> >* b, + template <int blocksize> + static double computeL2Error(const VirtualGridViewFunction<GridView,Dune::FieldVector<double,blocksize> >* a, + const Dune::VirtualFunction<Dune::FieldVector<ctype,dim>, Dune::FieldVector<double,blocksize> >* b, QuadratureRuleKey quadKey) -{ + { // The error to be computed double error = 0; @@ -55,13 +42,13 @@ computeL2Error(const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1 for (size_t i=0; i<quad.size(); i++) { // Evaluate function a - Dune::FieldVector<double,1> aValue; + Dune::FieldVector<double,blocksize> aValue; a->evaluateLocal(*eIt, quad[i].position(),aValue); // Evaluate function b. If it is a grid function use that to speed up the evaluation - Dune::FieldVector<double,1> bValue; - if (dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1> >*>(b)) - dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1> >*>(b)->evaluateLocal(*eIt, + Dune::FieldVector<double,blocksize> bValue; + if (dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,blocksize> >*>(b)) + dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,blocksize> >*>(b)->evaluateLocal(*eIt, quad[i].position(), bValue ); @@ -77,15 +64,15 @@ computeL2Error(const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1 return std::sqrt(error); -} + } -template <class GridView> -double DiscretizationError<GridView>:: -computeH1HalfNormDifferenceSquared(const GridView& gridView, - const VirtualDifferentiableFunction<Dune::FieldVector<ctype,dim>, Dune::FieldVector<double,1> >* u, - const VirtualDifferentiableFunction<Dune::FieldVector<ctype,dim>, Dune::FieldVector<double,1> >* v, - QuadratureRuleKey quadKey) -{ + + template <int blocksize> + static double computeH1HalfNormDifferenceSquared(const GridView& gridView, + const VirtualDifferentiableFunction<Dune::FieldVector<ctype,dim>, Dune::FieldVector<double,blocksize> >* u, + const VirtualDifferentiableFunction<Dune::FieldVector<ctype,dim>, Dune::FieldVector<double,blocksize> >* v, + QuadratureRuleKey quadKey) + { double norm = 0.0; typedef typename GridView::template Codim<0>::Iterator ElementIterator; @@ -106,16 +93,16 @@ computeH1HalfNormDifferenceSquared(const GridView& gridView, const double integrationElement = eIt->geometry().integrationElement(quadPos); // Evaluate integral - Dune::FieldMatrix<double,1,dim> u_di; - Dune::FieldMatrix<double,1,dim> v_di; + Dune::FieldMatrix<double,blocksize,dim> u_di; + Dune::FieldMatrix<double,blocksize,dim> v_di; - if (dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1> >*>(u)) - dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1> >*>(u)->evaluateDerivativeLocal(*eIt, quadPos, u_di); + if (dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,blocksize> >*>(u)) + dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,blocksize> >*>(u)->evaluateDerivativeLocal(*eIt, quadPos, u_di); else u->evaluateDerivative(eIt->geometry().global(quadPos), u_di); - if (dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1> >*>(v)) - dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,1> >*>(v)->evaluateDerivativeLocal(*eIt, quadPos, v_di); + if (dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,blocksize> >*>(v)) + dynamic_cast<const VirtualGridViewFunction<GridView,Dune::FieldVector<double,blocksize> >*>(v)->evaluateDerivativeLocal(*eIt, quadPos, v_di); else v->evaluateDerivative(eIt->geometry().global(quadPos), v_di); @@ -124,8 +111,10 @@ computeH1HalfNormDifferenceSquared(const GridView& gridView, } } return norm; -} + } +}; + #endif