1 #ifndef DUNE_FEM_SPACE_LAGRANGE_SPACE_HH 2 #define DUNE_FEM_SPACE_LAGRANGE_SPACE_HH 3 4 // C++ includes 5 #include <vector> 6 7 // dune-common includes 8 #include <dune/common/exceptions.hh> 9 10 // dune-geometry includes 11 #include <dune/geometry/type.hh> 12 13 // dune-fem includes 14 #include <dune/fem/common/hybrid.hh> 15 16 #include <dune/fem/space/basisfunctionset/default.hh> 17 #include <dune/fem/space/common/localinterpolation.hh> 18 #include <dune/fem/space/common/allgeomtypes.hh> 19 #include <dune/fem/space/common/basesetlocalkeystorage.hh> 20 #include <dune/fem/space/common/defaultcommhandler.hh> 21 #include <dune/fem/space/common/discretefunctionspace.hh> 22 #include <dune/fem/space/common/functionspace.hh> 23 #include <dune/fem/space/mapper/indexsetdofmapper.hh> 24 #include <dune/fem/space/mapper/nonblockmapper.hh> 25 #include <dune/fem/space/shapefunctionset/proxy.hh> 26 #include <dune/fem/space/shapefunctionset/selectcaching.hh> 27 #include <dune/fem/space/shapefunctionset/vectorial.hh> 28 #include <dune/fem/storage/singletonlist.hh> 29 30 // local includes 31 #include "adaptmanager.hh" 32 #include "capabilities.hh" 33 #include "interpolation.hh" 34 #include "lagrangepoints.hh" 35 #include "shapefunctionset.hh" 36 #include "storage.hh" 37 38 39 namespace Dune 40 { 41 42 namespace Fem 43 { 44 45 // Forward declaration 46 // ------------------- 47 48 template< class FunctionSpace, class GridPart, int maxPolOrder, class Storage = CachingStorage > 49 class LagrangeDiscreteFunctionSpace; 50 51 52 53 // LagrangeDiscreteFunctionSpaceTraits 54 // ----------------------------------- 55 56 template< class FunctionSpace, class GridPart, int maxPolOrder, class Storage > 57 struct LagrangeDiscreteFunctionSpaceTraits 58 { 59 typedef LagrangeDiscreteFunctionSpace< FunctionSpace, GridPart, maxPolOrder, Storage > DiscreteFunctionSpaceType; 60 61 // we take 6 as the maximal polynomial order for the dynamic space 62 static const int maxPolynomialOrder = ( maxPolOrder < 0 ) ? -maxPolOrder : maxPolOrder; 63 static const int minPolynomialOrder = ( maxPolOrder < 0 ) ? 1 : maxPolOrder; 64 65 typedef GridPart GridPartType; 66 typedef GridFunctionSpace< GridPartType, FunctionSpace > FunctionSpaceType; 67 68 typedef IndexSetDofMapper< GridPartType, 69 std::conditional_t< Dune::Capabilities::isCartesian< typename GridPart::GridType >::v, 70 DefaultLocalDofMapping< GridPart >, 71 LagrangeLocalDofMapping< GridPart > > > BlockMapperType; 72 73 typedef Hybrid::IndexRange< int, FunctionSpaceType::dimRange > LocalBlockIndices; 74 75 static const int codimension = 0; 76 77 private: 78 typedef typename GridPartType::template Codim< codimension >::EntityType EntityType; 79 static const int dimLocal = GridPartType::dimension; 80 typedef typename FunctionSpaceType::ScalarFunctionSpaceType ScalarFunctionSpaceType; 81 typedef typename ToNewDimDomainFunctionSpace< ScalarFunctionSpaceType, dimLocal >::Type ShapeFunctionSpaceType; 82 83 public: 84 typedef LagrangeShapeFunctionSet< ShapeFunctionSpaceType, maxPolynomialOrder > LagrangeShapeFunctionSetType; 85 typedef SelectCachingShapeFunctionSet< LagrangeShapeFunctionSetType, Storage > ScalarShapeFunctionSetType; 86 87 typedef ShapeFunctionSetProxy< ScalarShapeFunctionSetType > ScalarShapeFunctionSetProxyType; 88 typedef VectorialShapeFunctionSet< ScalarShapeFunctionSetProxyType, typename FunctionSpaceType::RangeType > ShapeFunctionSetType; 89 90 typedef Dune::Fem::DefaultBasisFunctionSet< EntityType, ShapeFunctionSetType > BasisFunctionSetType; 91 92 template< class DiscreteFunction, class Operation = Dune::Fem::DFCommunicationOperation::Add > 93 struct CommDataHandle 94 { 95 // type of data handle 96 typedef Dune::Fem::DefaultCommunicationHandler< DiscreteFunction, Operation > Type; 97 // type of operation to perform on scatter 98 typedef Operation OperationType; 99 }; 100 }; 101 102 103 104 /** \addtogroup LagrangeDiscreteFunctionSpace 105 * 106 * Provides access to bse function sets for different element types in 107 * one grid and size of function space and maps from local to global dof 108 * number. 109 * 110 * \note This space can only be used with special index sets. If you want 111 * to use the LagrangeDiscreteFunctionSpace with an index set only 112 * supporting the index set interface you will have to use the 113 * IndexSetWrapper class to provide the required functionality. 114 * 115 * \note For adaptive calculations one has to use index sets that are 116 * capable of adaption (i.e. the method adaptive returns true). See also 117 * AdaptiveLeafIndexSet. 118 */ 119 120 // LagrangeDiscreteFunctionSpace 121 // ----------------------------- 122 123 /** \class LagrangeDiscreteFunctionSpace 124 * \ingroup LagrangeDiscreteFunctionSpace 125 * \brief Lagrange discrete function space 126 */ 127 128 template< class FunctionSpace, class GridPart, int maxPolOrder, class Storage > 129 class LagrangeDiscreteFunctionSpace 130 : public DiscreteFunctionSpaceDefault< LagrangeDiscreteFunctionSpaceTraits< FunctionSpace, GridPart, maxPolOrder, Storage > > 131 { 132 133 typedef LagrangeDiscreteFunctionSpace< FunctionSpace, GridPart, maxPolOrder, Storage > ThisType; 134 typedef DiscreteFunctionSpaceDefault< LagrangeDiscreteFunctionSpaceTraits< FunctionSpace, GridPart, maxPolOrder, Storage > > BaseType; 135 136 public: 137 typedef typename BaseType::Traits Traits; 138 static const int maxPolynomialOrder = Traits :: maxPolynomialOrder; 139 static const int minPolynomialOrder = Traits :: minPolynomialOrder; 140 141 static_assert( (maxPolynomialOrder > 0), "LagrangeDiscreteFunctionSpace only defined for polOrder > 0" ); 142 143 typedef typename BaseType::FunctionSpaceType FunctionSpaceType; 144 145 typedef typename BaseType::GridPartType GridPartType; 146 typedef typename BaseType::GridType GridType; 147 typedef typename BaseType::IndexSetType IndexSetType; 148 typedef typename BaseType::EntityType EntityType; 149 typedef typename BaseType::IntersectionType IntersectionType; 150 151 typedef typename BaseType::Traits::ShapeFunctionSetType ShapeFunctionSetType; 152 typedef typename BaseType::BasisFunctionSetType BasisFunctionSetType; 153 154 typedef typename BaseType::BlockMapperType BlockMapperType; 155 156 typedef LagrangePointSet< GridPartType, maxPolynomialOrder > LagrangePointSetType; 157 typedef LagrangeLocalInterpolation< GridPartType, maxPolynomialOrder, BasisFunctionSetType > LocalInterpolationType; 158 typedef LocalInterpolationType InterpolationImplType; 159 160 typedef LocalInterpolationWrapper< ThisType > InterpolationType; 161 162 private: 163 typedef typename Traits::ScalarShapeFunctionSetType ScalarShapeFunctionSetType; 164 typedef BaseSetLocalKeyStorage< ScalarShapeFunctionSetType > ScalarShapeFunctionSetStorageType; 165 166 typedef CompiledLocalKeyContainer< LagrangePointSetType, minPolynomialOrder, maxPolynomialOrder > LagrangePointSetContainerType; 167 typedef typename LagrangePointSetContainerType::LocalKeyStorageType LocalKeyStorageType; 168 169 typedef LagrangeMapperSingletonKey< GridPartType, LocalKeyStorageType > 170 MapperSingletonKeyType; 171 typedef LagrangeMapperSingletonFactory< MapperSingletonKeyType, BlockMapperType > 172 BlockMapperSingletonFactoryType; 173 typedef SingletonList< MapperSingletonKeyType, BlockMapperType, BlockMapperSingletonFactoryType > BlockMapperProviderType; 174 175 // static const InterfaceType defaultInterface = InteriorBorder_InteriorBorder_Interface; 176 static const InterfaceType defaultInterface = GridPart::indexSetInterfaceType; 177 static const CommunicationDirection defaultDirection = ForwardCommunication; 178 179 180 template < int pOrd > 181 struct Initialize 182 { 183 struct ScalarShapeFunctionSetFactory 184 { 185 public: createObjectDune::Fem::LagrangeDiscreteFunctionSpace::Initialize::ScalarShapeFunctionSetFactory186 static ScalarShapeFunctionSetType *createObject ( const GeometryType &type ) 187 { 188 typedef typename Traits :: LagrangeShapeFunctionSetType LagrangeShapeFunctionSetType; 189 return new ScalarShapeFunctionSetType( type, LagrangeShapeFunctionSetType( type, pOrd ) ); 190 } 191 deleteObjectDune::Fem::LagrangeDiscreteFunctionSpace::Initialize::ScalarShapeFunctionSetFactory192 static void deleteObject ( ScalarShapeFunctionSetType *object ) { delete object; } 193 }; 194 applyDune::Fem::LagrangeDiscreteFunctionSpace::Initialize195 static void apply ( ScalarShapeFunctionSetStorageType& scalarShapeFunctionSets, 196 const int polynomialOrder, 197 const GeometryType &type ) 198 { 199 if( polynomialOrder == pOrd ) 200 { 201 typedef SingletonList< const GeometryType, ScalarShapeFunctionSetType, ScalarShapeFunctionSetFactory > SingletonProviderType; 202 scalarShapeFunctionSets.template insert< SingletonProviderType >( type ); 203 } 204 } 205 }; 206 207 208 public: 209 /////////////////////// 210 // Interface methods // 211 /////////////////////// 212 213 using BaseType::order; 214 LagrangeDiscreteFunctionSpace(GridPartType & gridPart,const InterfaceType commInterface,const CommunicationDirection commDirection)215 explicit LagrangeDiscreteFunctionSpace ( GridPartType &gridPart, 216 const InterfaceType commInterface, 217 const CommunicationDirection commDirection ) 218 : LagrangeDiscreteFunctionSpace( gridPart, minPolynomialOrder, commInterface, commDirection ) 219 {} 220 LagrangeDiscreteFunctionSpace(GridPartType & gridPart,const int polOrder=minPolynomialOrder,const InterfaceType commInterface=defaultInterface,const CommunicationDirection commDirection=defaultDirection)221 explicit LagrangeDiscreteFunctionSpace ( GridPartType &gridPart, 222 const int polOrder = minPolynomialOrder, 223 const InterfaceType commInterface = defaultInterface, 224 const CommunicationDirection commDirection = defaultDirection ) 225 : BaseType( gridPart, commInterface, commDirection ), 226 blockMapper_(), 227 lagrangePointSetContainer_( gridPart ), 228 // when min == max polynomial order the dynamic choice is off 229 polynomialOrder_( ( minPolynomialOrder == maxPolynomialOrder ) ? minPolynomialOrder : polOrder ) 230 { 231 const IndexSetType &indexSet = gridPart.indexSet(); 232 233 AllGeomTypes< IndexSetType, GridType > allGeometryTypes( indexSet ); 234 const std::vector< GeometryType > &geometryTypes = allGeometryTypes.geomTypes( 0 ); 235 236 // create shape function sets 237 for( unsigned int i = 0; i < geometryTypes.size(); ++i ) 238 { 239 Fem::ForLoop< Initialize, minPolynomialOrder, maxPolynomialOrder >:: 240 apply( scalarShapeFunctionSets_, polynomialOrder_, geometryTypes[ i ] ); 241 } 242 243 MapperSingletonKeyType key( gridPart, lagrangePointSetContainer_.compiledLocalKeys( polynomialOrder_ ), polynomialOrder_ ); 244 blockMapper_.reset( &BlockMapperProviderType::getObject( key ) ); 245 assert( blockMapper_ ); 246 } 247 248 /** \copydoc Dune::Fem::DiscreteFunctionSpaceInterface::type */ type() const249 DFSpaceIdentifier type () const 250 { 251 return LagrangeSpace_id; 252 } 253 254 /** \copydoc Dune::Fem::DiscreteFunctionSpaceInterface::basisFunctionSet */ basisFunctionSet(const EntityType & entity) const255 const BasisFunctionSetType basisFunctionSet ( const EntityType &entity ) const 256 { 257 return BasisFunctionSetType( entity, shapeFunctionSet( entity ) ); 258 } 259 260 /** \copydoc Dune::Fem::DiscreteFunctionSpaceInterface::continuous */ continuous() const261 bool continuous () const 262 { 263 return true; 264 } 265 266 /** \copydoc Dune::Fem::DiscreteFunctionSpaceInterface::continuous */ continuous(const IntersectionType & intersection) const267 bool continuous ( const IntersectionType &intersection ) const 268 { 269 return intersection.conforming(); 270 } 271 272 /** \copydoc Dune::Fem::DiscreteFunctionSpaceInterface::order */ order() const273 int order () const 274 { 275 return polynomialOrder_; 276 } 277 278 /** \copydoc Dune::Fem::DiscreteFunctionSpaceInterface::blockMapper */ blockMapper() const279 BlockMapperType &blockMapper () const 280 { 281 assert( blockMapper_ ); 282 return *blockMapper_; 283 } 284 285 /////////////////////////// 286 // Non-interface methods // 287 /////////////////////////// 288 289 /** \brief return interpolation object 290 */ interpolation() const291 InterpolationType interpolation () const 292 { 293 return InterpolationType( *this ); 294 } 295 296 /** \brief return local interpolation for given entity 297 * 298 * \param[in] entity grid part entity 299 */ 300 [[deprecated]] interpolation(const EntityType & entity) const301 LocalInterpolationType interpolation ( const EntityType &entity ) const 302 { 303 return LocalInterpolationType( lagrangePointSetContainer_.compiledLocalKey( entity.type(), polynomialOrder_ ), basisFunctionSet( entity ) ); 304 } 305 306 /** \brief return local interpolation for given entity 307 * 308 * \param[in] entity grid part entity 309 */ localInterpolation(const EntityType & entity) const310 LocalInterpolationType localInterpolation ( const EntityType &entity ) const 311 { 312 return LocalInterpolationType( lagrangePointSetContainer_.compiledLocalKey( entity.type(), polynomialOrder_ ), basisFunctionSet( entity ) ); 313 } 314 315 /** \brief return shape function set for given entity 316 * 317 * \param[in] entity entity (of codim 0) for which shape function set 318 * is requested 319 * 320 * \returns ShapeFunctionSetType shape function set 321 */ shapeFunctionSet(const EntityType & entity) const322 ShapeFunctionSetType shapeFunctionSet ( const EntityType &entity ) const 323 { 324 return shapeFunctionSet( entity.type() ); 325 } 326 327 /** \brief return shape unique function set for geometry type 328 * 329 * \param[in] type geometry type (must be a cube) for which 330 * shape function set is requested 331 * 332 * \returns ShapeFunctionSetType shape function set 333 */ shapeFunctionSet(const GeometryType & type) const334 ShapeFunctionSetType shapeFunctionSet ( const GeometryType &type ) const 335 { 336 return ShapeFunctionSetType( &scalarShapeFunctionSets_[ type ] ); 337 } 338 339 LagrangeDiscreteFunctionSpace ( const ThisType & ) = delete; 340 ThisType &operator= ( const ThisType & ) = delete; 341 342 protected: 343 std::unique_ptr< BlockMapperType, typename BlockMapperProviderType::Deleter> blockMapper_; 344 ScalarShapeFunctionSetStorageType scalarShapeFunctionSets_; 345 LagrangePointSetContainerType lagrangePointSetContainer_; 346 const int polynomialOrder_; 347 }; 348 349 } // namespace Fem 350 351 } // namespace Dune 352 353 #endif // #ifndef DUNE_FEM_SPACE_LAGRANGE_SPACE_HH 354