Skip to content
Snippets Groups Projects
Select Git revision
  • 3f559fd17289f9753b03837f1844d10fc3460722
  • master default
  • dune-tkr-article
  • p/lh188/test-assembler-bug
  • mention-makefunction
  • temp/test-CI-with-subgrid-enabled
  • releases/2.5-1
  • releases/2.4-2
  • improve_grid
  • experimental/test-core-ci
  • feature/dune-functions-assemblers
  • feature/DG-Transferoperator-Assembler
  • debug/functionspacebasistest
  • releases/2.4-1
  • feature/dune-functions-assembler-with-skeleton
  • feature/elastictensor-interface
  • releases/2.3-2
  • maxka/conformingbasis
  • releases/2.3-1
  • subgridassembler-rework
  • releases/2.2-1
  • dune-tkr-article-patched
  • dune-tkr-article-base
  • subversion->git
24 results

h1functionalassembler.hh

  • Forked from agnumpde / dune-fufem
    1365 commits behind the upstream repository.
    Carsten Gräser's avatar
    graeser authored and graeser committed
    [[Imported from SVN: r8367]]
    3f559fd1
    History
    Code owners
    Assign users and groups as approvers for specific file changes. Learn more.
    h1functionalassembler.hh 4.54 KiB
    #ifndef H1_FUNCTIONAL_ASSEMBLER_HH
    #define H1_FUNCTIONAL_ASSEMBLER_HH
    
    
    #include <dune/common/fvector.hh>
    #include <dune/common/function.hh>
    #include <dune/istl/bvector.hh>
    
    
    #include <dune/fufem/quadraturerules/quadraturerulecache.hh>
    #include <dune/fufem/functions/virtualgridfunction.hh>
    
    #include <dune/fufem/assemblers/localfunctionalassembler.hh>
    
    /**  \brief Local assembler for finite element H^1-functionals
    */
    template <class GridType, class T=Dune::FieldVector<double,1> >
    class H1FunctionalAssembler :
        public LocalFunctionalAssembler<GridType, T>
    
    {
        private:
            typedef typename GridType::template Codim<0>::Geometry::GlobalCoordinate GlobalCoordinate;
            typedef typename GridType::template Codim<0>::Geometry::LocalCoordinate LocalCoordinate;
    
            static const int dim = GridType::dimension;
            static const int dimworld = GridType::dimensionworld;
    
        public:
            typedef typename LocalFunctionalAssembler<GridType,T>::Element Element;
            typedef typename LocalFunctionalAssembler<GridType,T>::Element::Geometry Geometry;
            typedef typename LocalFunctionalAssembler<GridType,T>::LocalVector LocalVector;
    
            typedef typename Dune::VirtualFunction<GlobalCoordinate, typename DerivativeTypefier<GlobalCoordinate,T>::DerivativeType > Function;
    
            /** \brief constructor
              *
              * \param f the function representing the functional (if assembling e.g. \f$(\nabla u,\nabla v)\f$ you need to provide \f$f=\nabla u\f$)
              * \param order the quadrature order (DEFAULT 2)
              */
            H1FunctionalAssembler(const Function& f, int order=2) :
                f_(f),
                order_(order)
            {}
    
            /** \copydoc LocalFunctionalAssembler::assemble
             */
            template <class TrialLocalFE>
            void assemble(const Element& element, LocalVector& localVector, const TrialLocalFE& tFE)
            {
                typedef typename GridType::template Codim<0>::Geometry::Jacobian GeoJacobianInvTransposed;
                typedef typename Function::RangeType RangeType;
                typedef typename TrialLocalFE::Traits::LocalBasisType::Traits::JacobianType JacobianType;
    
                // Make sure we got suitable shape functions
                assert(tFE.type() == element.type());
    
                // get geometry and store it
                const Geometry geometry = element.geometry();
    
                localVector = 0.0;
    
                // get quadrature rule
    //            const Dune::template QuadratureRule<double, dim>& quad = Dune::template QuadratureRules<double, dim>::rule(element.type(), order_);
                const Dune::template QuadratureRule<double, dim>& quad = QuadratureRuleCache<double, dim>::rule(element.type(), order_, IsRefinedLocalFiniteElement<TrialLocalFE>::value(tFE) );
    
                // store gradients of shape functions and base functions
                std::vector<JacobianType> referenceGradients(tFE.localBasis().size());
                std::vector<GlobalCoordinate> gradients(tFE.localBasis().size());
    
                // loop over quadrature points
                for (size_t pt=0; pt < quad.size(); ++pt)
                {
                    // get quadrature point
                    const LocalCoordinate& quadPos = quad[pt].position();
    
                    // get transposed inverse of Jacobian of transformation
                    const GeoJacobianInvTransposed& invJacobian = geometry.jacobianInverseTransposed(quadPos);
    
                    // get integration factor
                    const double integrationElement = geometry.integrationElement(quadPos);
    
                    // get gradients of shape functions
                    tFE.localBasis().evaluateJacobian(quadPos, referenceGradients);
    
                    // transform gradients
                    for (size_t i=0; i<gradients.size(); ++i)
                        invJacobian.mv(referenceGradients[i][0], gradients[i]);
    
                    // compute values of function
                    RangeType f_pos;
                    if (dynamic_cast<const VirtualGridFunction<GridType,RangeType>*>(&f_))
                        dynamic_cast<const VirtualGridFunction<GridType,RangeType>*>(&f_)->evaluateLocal(element, quadPos, f_pos);
                    else
                        f_.evaluate(geometry.global(quadPos), f_pos);
    
                    // and vector entries
                    for (size_t i=0; i<gradients.size(); ++i)
                    {
                        T dummy;
                        f_pos.mv(gradients[i],dummy);
                        localVector[i].axpy(quad[pt].weight()*integrationElement, dummy);
                    }
                }
                return;
            }
    
        private:
            const Function& f_;
            const int order_;
    };
    
    #endif