dune-localfunctions  2.2.0
pk1d.hh
Go to the documentation of this file.
00001 // -*- tab-width: 4; indent-tabs-mode: nil -*-
00002 // vi: set ts=4 sw=2 et sts=2:
00003 #ifndef DUNE_PK1DLOCALFINITEELEMENT_HH
00004 #define DUNE_PK1DLOCALFINITEELEMENT_HH
00005 
00006 #include <cstddef>
00007 
00008 #include <dune/geometry/type.hh>
00009 
00010 #include <dune/localfunctions/common/localfiniteelementtraits.hh>
00011 #include <dune/localfunctions/common/localtoglobaladaptors.hh>
00012 #include "pk1d/pk1dlocalbasis.hh"
00013 #include "pk1d/pk1dlocalcoefficients.hh"
00014 #include "pk1d/pk1dlocalinterpolation.hh"
00015 
00016 namespace Dune 
00017 {
00018 
00021   template<class D, class R, unsigned int k>
00022   class Pk1DLocalFiniteElement 
00023   {
00024   public:
00027         typedef LocalFiniteElementTraits<Pk1DLocalBasis<D,R,k>,
00028                                          Pk1DLocalCoefficients<k>,
00029                                          Pk1DLocalInterpolation<Pk1DLocalBasis<D,R,k> > > Traits;
00030 
00033         Pk1DLocalFiniteElement ()
00034         {
00035           gt.makeLine();
00036         }
00037 
00040         Pk1DLocalFiniteElement (int variant) : coefficients(variant)
00041         {
00042           gt.makeLine();
00043         }
00044 
00049         Pk1DLocalFiniteElement (const unsigned int vertexmap[3]) : coefficients(vertexmap)
00050         {
00051           gt.makeLine();
00052         }
00053 
00056         const typename Traits::LocalBasisType& localBasis () const
00057         {
00058           return basis;
00059         }
00060         
00063         const typename Traits::LocalCoefficientsType& localCoefficients () const
00064         {
00065           return coefficients;
00066         }
00067         
00070         const typename Traits::LocalInterpolationType& localInterpolation () const
00071         {
00072           return interpolation;
00073         }
00074         
00077         GeometryType type () const
00078         {
00079           return gt;
00080         }
00081 
00082     Pk1DLocalFiniteElement* clone () const
00083     {
00084       return new Pk1DLocalFiniteElement(*this);
00085     }
00086 
00087   private:
00088     Pk1DLocalBasis<D,R,k> basis;
00089     Pk1DLocalCoefficients<k> coefficients;
00090     Pk1DLocalInterpolation<Pk1DLocalBasis<D,R,k> > interpolation;
00091     GeometryType gt;
00092   };
00093 
00095 
00102   template<class Geometry, class RF, std::size_t k>
00103   class Pk1DFiniteElement {
00104     typedef typename Geometry::ctype DF;
00105     typedef Pk1DLocalBasis<DF,RF,k> LocalBasis;
00106     typedef Pk1DLocalInterpolation<LocalBasis> LocalInterpolation;
00107 
00108   public:
00112     struct Traits {
00113       typedef ScalarLocalToGlobalBasisAdaptor<LocalBasis, Geometry> Basis;
00114       typedef LocalToGlobalInterpolationAdaptor<
00115         LocalInterpolation,
00116         typename Basis::Traits
00117         > Interpolation;
00118       typedef Pk1DLocalCoefficients<k> Coefficients;
00119     };
00120 
00121   private:
00122     static const GeometryType gt;
00123     static const LocalBasis localBasis;
00124     static const LocalInterpolation localInterpolation;
00125 
00126     typename Traits::Basis basis_;
00127     typename Traits::Interpolation interpolation_;
00128     typename Traits::Coefficients coefficients_;
00129 
00130   public:
00132 
00145     template<class VertexOrder>
00146     Pk1DFiniteElement(const Geometry &geometry,
00147                       const VertexOrder& vertexOrder) :
00148       basis_(localBasis, geometry), interpolation_(localInterpolation),
00149       coefficients_(vertexOrder.begin(0, 0))
00150     { }
00151 
00152     const typename Traits::Basis& basis() const { return basis_; }
00153     const typename Traits::Interpolation& interpolation() const
00154     { return interpolation_; }
00155     const typename Traits::Coefficients& coefficients() const
00156     { return coefficients_; }
00157     const GeometryType &type() const { return gt; }
00158   };
00159 
00160   template<class Geometry, class RF, std::size_t k>
00161   const GeometryType
00162   Pk1DFiniteElement<Geometry, RF, k>::gt(GeometryType::simplex, 2);
00163 
00164   template<class Geometry, class RF, std::size_t k>
00165   const typename Pk1DFiniteElement<Geometry, RF, k>::LocalBasis
00166   Pk1DFiniteElement<Geometry, RF, k>::localBasis = LocalBasis();
00167 
00168   template<class Geometry, class RF, std::size_t k>
00169   const typename Pk1DFiniteElement<Geometry, RF, k>::LocalInterpolation
00170   Pk1DFiniteElement<Geometry, RF, k>::localInterpolation =
00171     LocalInterpolation();
00172 
00174 
00184   template<class Geometry, class RF, std::size_t k>
00185   struct Pk1DFiniteElementFactory {
00186     typedef Pk1DFiniteElement<Geometry, RF, k> FiniteElement;
00187 
00189 
00203     template<class VertexOrder>
00204     const FiniteElement make(const Geometry& geometry,
00205                              const VertexOrder& vertexOrder)
00206     { return FiniteElement(geometry, vertexOrder); }
00207   };
00208 }
00209 
00210 #endif