1 #ifndef DUNE_FEM_SPACE_BREZZIDOUGLASMARINI_HH
2 #define DUNE_FEM_SPACE_BREZZIDOUGLASMARINI_HH
3 
4 #if HAVE_DUNE_LOCALFUNCTIONS
5 
6 #include <dune/localfunctions/brezzidouglasmarini/brezzidouglasmarini1cube2d.hh>
7 #include <dune/localfunctions/brezzidouglasmarini/brezzidouglasmarini1cube3d.hh>
8 #include <dune/localfunctions/brezzidouglasmarini/brezzidouglasmarini1simplex2d.hh>
9 #include <dune/localfunctions/brezzidouglasmarini/brezzidouglasmarini2cube2d.hh>
10 #include <dune/localfunctions/brezzidouglasmarini/brezzidouglasmarini2simplex2d.hh>
11 
12 #include <dune/fem/space/common/uniquefacetorientation.hh>
13 #include <dune/fem/space/basisfunctionset/piolatransformation.hh>
14 #include <dune/fem/space/localfiniteelement/space.hh>
15 
16 namespace Dune
17 {
18 
19   namespace Fem
20   {
21 
22     namespace Impl
23     {
24 
25       // BDMLocalFiniteElement
26       // ---------------------
27 
28       template< unsigned int id, class DomainField, class RangeField, int dimension, int order >
29       struct BDMLocalFiniteElement
30       {
31         static_assert( AlwaysFalse< DomainField >::value, "BDMLocalFiniteElement not implemented for your choice." );
32       };
33 
34       // The following local finite elements are implemented
35 
36       // 2d, Cube, first order
37       template< class D >
38       struct BDMLocalFiniteElement< Dune::GeometryTypes::cube(2).id(), D, D, 2, 1 >
39         : public BDM1Cube2DLocalFiniteElement< D, D >
40       {
41         static const int numOrientations = 16;
42         template< class ... Args >
BDMLocalFiniteElementDune::Fem::Impl::BDMLocalFiniteElement43         BDMLocalFiniteElement ( Args && ... args )
44           : BDM1Cube2DLocalFiniteElement< D, D >( std::forward< Args >( args ) ... ) {}
45       };
46 
47       // 3d, Cube, first order
48       template< class D >
49       struct BDMLocalFiniteElement< Dune::GeometryTypes::cube(3).id(), D, D, 3, 1 >
50         : public BDM1Cube3DLocalFiniteElement< D, D >
51       {
52         static const int numOrientations = 64;
53         template< class ... Args >
BDMLocalFiniteElementDune::Fem::Impl::BDMLocalFiniteElement54         BDMLocalFiniteElement ( Args && ... args )
55           : BDM1Cube3DLocalFiniteElement< D, D >( std::forward< Args >( args ) ... ) {}
56       };
57 
58       // 2d, Cube, second order
59       template< class D >
60       struct BDMLocalFiniteElement< Dune::GeometryTypes::cube(2).id(), D, D, 2, 2 >
61         : public BDM2Cube2DLocalFiniteElement< D, D >
62       {
63         static const int numOrientations = 16;
64         template< class ... Args >
BDMLocalFiniteElementDune::Fem::Impl::BDMLocalFiniteElement65         BDMLocalFiniteElement ( Args && ... args )
66           : BDM2Cube2DLocalFiniteElement< D, D >( std::forward< Args >( args ) ... ) {}
67       };
68 
69 
70       // 2d, simplex, first order
71       template< class D >
72       struct BDMLocalFiniteElement< Dune::GeometryTypes::simplex(2).id(), D, D, 2, 1 >
73         : public BDM1Simplex2DLocalFiniteElement< D, D >
74       {
75         static const int numOrientations = 8;
76         template< class ... Args >
BDMLocalFiniteElementDune::Fem::Impl::BDMLocalFiniteElement77         BDMLocalFiniteElement ( Args && ... args )
78           : BDM1Simplex2DLocalFiniteElement< D, D >( std::forward< Args >( args ) ... ) {}
79       };
80 
81       // 2d, simplex, second order
82       template< class D >
83       struct BDMLocalFiniteElement< Dune::GeometryTypes::simplex(2).id(), D, D, 2, 2 >
84         : public BDM2Simplex2DLocalFiniteElement< D, D >
85       {
86         static const int numOrientations = 8;
87         template< class ... Args >
BDMLocalFiniteElementDune::Fem::Impl::BDMLocalFiniteElement88         BDMLocalFiniteElement ( Args && ... args )
89           : BDM2Simplex2DLocalFiniteElement< D, D >( std::forward< Args >( args ) ... ) {}
90       };
91 
92     }
93 
94     // BrezziDouglasMariniLocalFiniteElementMap
95     // ----------------------------------------
96 
97     template< class GridPart, class FunctionSpace, int polOrder >
98     class BrezziDouglasMariniLocalFiniteElementMap
99     {
100       typedef BrezziDouglasMariniLocalFiniteElementMap< GridPart, FunctionSpace, polOrder > ThisType;
101       static_assert( Dune::Fem::GridPartCapabilities::hasSingleGeometryType< GridPart >::v,
102                      "GridPart has more than one geometry type." );
103 
104       static const unsigned int topologyId = Dune::Fem::GridPartCapabilities::hasSingleGeometryType< GridPart >::topologyId;
105 
106     public:
107       typedef std::tuple< > KeyType;
108 
109       typedef GridPart GridPartType;
110 
111       typedef typename FunctionSpace::DomainFieldType DomainFieldType;
112       typedef typename FunctionSpace::RangeFieldType RangeFieldType;
113 
114       typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
115 
116       typedef PiolaTransformation< typename EntityType::Geometry, FunctionSpace::dimRange > TransformationType;
117 
118       static const int dimLocal = GridPart::dimension;
119 
120       typedef Impl::BDMLocalFiniteElement< topologyId, DomainFieldType, RangeFieldType, dimLocal, polOrder > LocalFiniteElementType;
121       typedef typename LocalFiniteElementType::Traits::LocalBasisType LocalBasisType;
122       typedef typename LocalFiniteElementType::Traits::LocalCoefficientsType   LocalCoefficientsType;
123       typedef typename LocalFiniteElementType::Traits::LocalInterpolationType  LocalInterpolationType;
124 
125       template< class ... Args >
BrezziDouglasMariniLocalFiniteElementMap(const GridPart & gridPart,Args...args)126       BrezziDouglasMariniLocalFiniteElementMap ( const GridPart &gridPart, Args ... args )
127         : orientation_( gridPart )
128       {
129         for( std::size_t i = 0; i < LocalFiniteElementType::numOrientations; ++i )
130           map_[ i ] = LocalFiniteElementType( i );
131       }
132 
size()133       static std::size_t size () { return LocalFiniteElementType::numOrientations; }
134 
order() const135       int order () const { return polOrder; }
136 
137       template< class Entity >
order(const Entity & entity) const138       int order ( const Entity &entity ) const { return order(); }
139 
140       template< class Entity >
operator ()(const Entity & e) const141       auto operator() ( const Entity &e ) const
142       {
143         unsigned int orient = orientation_( e );
144         return std::tuple< std::size_t, const LocalBasisType&, const LocalInterpolationType& >
145         { static_cast< std::size_t >( orient ),
146           map_[ orient ].localBasis(),
147           map_[ orient ].localInterpolation() };
148       }
149 
hasCoefficients(const GeometryType & t) const150       bool hasCoefficients ( const GeometryType &t ) const
151       {
152         Dune::GeometryType type( GridPartCapabilities::hasSingleGeometryType< GridPart >::topologyId, GridPart::dimension );
153         return (type == t);
154       }
155 
localCoefficients(const GeometryType & type) const156       const LocalCoefficientsType &localCoefficients ( const GeometryType &type ) const
157       {
158         return map_[ 0 ].localCoefficients();
159       }
160 
gridPart() const161       const GridPartType &gridPart () const { return orientation_.gridPart(); }
162 
163     private:
164       UniqueFacetOrientation< GridPartType > orientation_;
165       std::array< LocalFiniteElementType, LocalFiniteElementType::numOrientations > map_;
166     };
167 
168 
169     // BrezziDouglasMariniSpace
170     // ------------------------
171 
172     template< class FunctionSpace, class GridPart, int polOrder, class Storage = CachingStorage >
173     using BrezziDouglasMariniSpace
174             = LocalFiniteElementSpace< BrezziDouglasMariniLocalFiniteElementMap< GridPart, FunctionSpace, polOrder >,
175                                        FunctionSpace, Storage >;
176 
177   } // namespace Fem
178 
179 } // namespace Dune
180 
181 #endif // #if HAVE_DUNE_LOCALFUNCTIONS
182 
183 #endif // #ifndef DUNE_FEM_SPACE_BREZZIDOUGLASMARINI_HH
184