dune-localfunctions  2.1.1
pk2d.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_Pk2DLOCALFINITEELEMENT_HH
00004 #define DUNE_Pk2DLOCALFINITEELEMENT_HH
00005 
00006 #include <cstddef>
00007 
00008 #include <dune/common/geometrytype.hh>
00009 
00010 #include <dune/localfunctions/common/localfiniteelementtraits.hh>
00011 #include <dune/localfunctions/common/localtoglobaladaptors.hh>
00012 #include "pk2d/pk2dlocalbasis.hh"
00013 #include "pk2d/pk2dlocalcoefficients.hh"
00014 #include "pk2d/pk2dlocalinterpolation.hh"
00015 
00016 namespace Dune 
00017 {
00018 
00021   template<class D, class R, unsigned int k>
00022   class Pk2DLocalFiniteElement 
00023   {
00024   public:
00027         typedef LocalFiniteElementTraits<Pk2DLocalBasis<D,R,k>,
00028                                          Pk2DLocalCoefficients<k>,
00029                                          Pk2DLocalInterpolation<Pk2DLocalBasis<D,R,k> > > Traits;
00030 
00033         Pk2DLocalFiniteElement ()
00034         {
00035           gt.makeTriangle();
00036         }
00037 
00040         Pk2DLocalFiniteElement (int variant) : coefficients(variant)
00041         {
00042           gt.makeTriangle();
00043         }
00044 
00051         Pk2DLocalFiniteElement (const unsigned int vertexmap[3]) : coefficients(vertexmap)
00052         {
00053           gt.makeTriangle();
00054         }
00055 
00058         const typename Traits::LocalBasisType& localBasis () const
00059         {
00060           return basis;
00061         }
00062         
00065         const typename Traits::LocalCoefficientsType& localCoefficients () const
00066         {
00067           return coefficients;
00068         }
00069         
00072         const typename Traits::LocalInterpolationType& localInterpolation () const
00073         {
00074           return interpolation;
00075         }
00076         
00079         GeometryType type () const
00080         {
00081           return gt;
00082         }
00083 
00084     Pk2DLocalFiniteElement* clone () const
00085     {
00086       return new Pk2DLocalFiniteElement(*this);
00087     }
00088 
00089   private:
00090         Pk2DLocalBasis<D,R,k> basis;
00091         Pk2DLocalCoefficients<k> coefficients;
00092         Pk2DLocalInterpolation<Pk2DLocalBasis<D,R,k> > interpolation;
00093         GeometryType gt;
00094   };
00095 
00097 
00104   template<class Geometry, class RF, std::size_t k>
00105   class Pk2DFiniteElement {
00106     typedef typename Geometry::ctype DF;
00107     typedef Pk2DLocalBasis<DF,RF,k> LocalBasis;
00108     typedef Pk2DLocalInterpolation<LocalBasis> LocalInterpolation;
00109 
00110   public:
00114     struct Traits {
00115       typedef ScalarLocalToGlobalBasisAdaptor<LocalBasis, Geometry> Basis;
00116       typedef LocalToGlobalInterpolationAdaptor<
00117         LocalInterpolation,
00118         typename Basis::Traits
00119         > Interpolation;
00120       typedef Pk2DLocalCoefficients<k> Coefficients;
00121     };
00122 
00123   private:
00124     static const GeometryType gt;
00125     static const LocalBasis localBasis;
00126     static const LocalInterpolation localInterpolation;
00127 
00128     typename Traits::Basis basis_;
00129     typename Traits::Interpolation interpolation_;
00130     typename Traits::Coefficients coefficients_;
00131 
00132   public:
00134 
00147     template<class VertexOrder>
00148     Pk2DFiniteElement(const Geometry &geometry,
00149                       const VertexOrder& vertexOrder) :
00150       basis_(localBasis, geometry), interpolation_(localInterpolation),
00151       coefficients_(vertexOrder.begin(0, 0))
00152     { }
00153 
00154     const typename Traits::Basis& basis() const { return basis_; }
00155     const typename Traits::Interpolation& interpolation() const
00156     { return interpolation_; }
00157     const typename Traits::Coefficients& coefficients() const
00158     { return coefficients_; }
00159     const GeometryType &type() const { return gt; }
00160   };
00161 
00162   template<class Geometry, class RF, std::size_t k>
00163   const GeometryType
00164   Pk2DFiniteElement<Geometry, RF, k>::gt(GeometryType::simplex, 2);
00165 
00166   template<class Geometry, class RF, std::size_t k>
00167   const typename Pk2DFiniteElement<Geometry, RF, k>::LocalBasis
00168   Pk2DFiniteElement<Geometry, RF, k>::localBasis = LocalBasis();
00169 
00170   template<class Geometry, class RF, std::size_t k>
00171   const typename Pk2DFiniteElement<Geometry, RF, k>::LocalInterpolation
00172   Pk2DFiniteElement<Geometry, RF, k>::localInterpolation =
00173     LocalInterpolation();
00174 
00176 
00186   template<class Geometry, class RF, std::size_t k>
00187   struct Pk2DFiniteElementFactory {
00188     typedef Pk2DFiniteElement<Geometry, RF, k> FiniteElement;
00189 
00191 
00205     template<class VertexOrder>
00206     const FiniteElement make(const Geometry& geometry,
00207                              const VertexOrder& vertexOrder)
00208     { return FiniteElement(geometry, vertexOrder); }
00209   };
00210 }
00211 
00212 #endif