3 #ifndef DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS_RAVIARTTHOMASSIMPLEX_RAVIARTTHOMASSIMPLEXINTERPOLATION_HH 4 #define DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS_RAVIARTTHOMASSIMPLEX_RAVIARTTHOMASSIMPLEXINTERPOLATION_HH 9 #include <dune/common/exceptions.hh> 11 #include <dune/geometry/topologyfactory.hh> 12 #include <dune/geometry/quadraturerules.hh> 13 #include <dune/geometry/referenceelements.hh> 14 #include <dune/geometry/type.hh> 27 template <
unsigned int dim >
30 template <
unsigned int dim,
class Field >
43 template <
class Setter>
46 setter.setLocalKeys(localKey_);
52 return localKey_[ i ];
57 return localKey_.size();
61 std::vector< LocalKey > localKey_;
69 template <
unsigned int dim >
72 static const unsigned int dimension = dim;
74 typedef unsigned int Key;
83 template <
unsigned int dim >
85 :
public TopologyFactory< RaviartThomasCoefficientsFactoryTraits< dim > >
89 template<
class Topology >
93 if( !supports< Topology >( key ) )
95 typename InterpolationFactory::Object *interpolation = InterpolationFactory::template create< Topology >( key );
97 InterpolationFactory::release( interpolation );
101 template<
class Topology >
119 template<
unsigned int dim,
class Field >
122 static const unsigned int dimension = dim;
127 typedef FieldVector< Field, dimension >
Normal;
136 TestBasisFactory::release( testBasis_ );
137 for( FaceStructure &f : faceStructure_ )
138 TestFaceBasisFactory::release( f.basis_ );
143 GeometryType
type ()
const {
return GeometryType( topologyId(), dimension ); }
145 unsigned int order ()
const {
return order_; }
147 unsigned int faceSize ()
const {
return faceSize_; }
150 TestFaceBasis *
testFaceBasis (
unsigned int f )
const { assert( f < faceSize() );
return faceStructure_[ f ].basis_; }
152 const Normal &
normal (
unsigned int f )
const { assert( f < faceSize() );
return *(faceStructure_[ f ].normal_); }
154 template<
class Topology >
158 topologyId_ = Topology::id;
160 testBasis_ = (order > 0 ? TestBasisFactory::template create< Topology >( order-1 ) :
nullptr);
162 const auto &refElement = ReferenceElements< Field, dimension >::general( type() );
163 faceSize_ = refElement.size( 1 );
164 faceStructure_.reserve( faceSize_ );
165 for(
unsigned int face = 0; face < faceSize_; ++face )
167 TestFaceBasis *faceBasis = Impl::IfTopology< CreateFaceBasis, dimension-1 >::apply( refElement.type( face, 1 ).id(), order );
168 faceStructure_.emplace_back( faceBasis, refElement.integrationOuterNormal( face ) );
170 assert( faceStructure_.size() == faceSize_ );
176 FaceStructure( TestFaceBasis *tfb,
const Normal &n )
177 : basis_( tfb ), normal_( &n )
180 TestFaceBasis *basis_;
181 const Dune::FieldVector< Field, dimension > *normal_;
184 template<
class FaceTopology >
185 struct CreateFaceBasis
187 static TestFaceBasis *apply (
unsigned int order ) {
return TestFaceBasisFactory::template create< FaceTopology >( order ); }
190 std::vector< FaceStructure > faceStructure_;
191 TestBasis *testBasis_ =
nullptr;
192 unsigned int topologyId_, order_, faceSize_;
205 template<
unsigned int dimension,
class F>
220 template<
class Function,
class Fy >
221 void interpolate (
const Function &
function, std::vector< Fy > &coefficients )
const 223 coefficients.resize(
size());
227 template<
class Basis,
class Matrix >
230 matrix.resize(
size(), basis.size() );
243 template <
class Topology>
248 builder_.template build<Topology>(order_);
249 if (builder_.testBasis())
250 size_ += dimension*builder_.testBasis()->size();
251 for (
unsigned int f=0; f<builder_.faceSize(); ++f )
252 if (builder_.testFaceBasis(f))
253 size_ += builder_.testFaceBasis(f)->size();
259 unsigned int row = 0;
260 for (
unsigned int f=0; f<builder_.faceSize(); ++f)
262 if (builder_.faceSize())
263 for (
unsigned int i=0; i<builder_.testFaceBasis(f)->size(); ++i,++row)
266 if (builder_.testBasis())
267 for (
unsigned int i=0; i<builder_.testBasis()->size()*dimension; ++i,++row)
269 assert( row ==
size() );
273 template<
class Func,
class Container,
bool type >
276 const Dune::GeometryType geoType( builder_.topologyId(), dimension );
278 std::vector< Field > testBasisVal;
280 for (
unsigned int i=0; i<
size(); ++i)
281 for (
unsigned int j=0; j<func.size(); ++j)
284 unsigned int row = 0;
287 typedef Dune::QuadratureRule<Field, dimension-1> FaceQuadrature;
288 typedef Dune::QuadratureRules<Field, dimension-1> FaceQuadratureRules;
290 typedef Dune::ReferenceElements< Field, dimension > RefElements;
291 typedef Dune::ReferenceElement< Field, dimension > RefElement;
292 typedef typename RefElement::template Codim< 1 >::Geometry Geometry;
294 const RefElement &refElement = RefElements::general( geoType );
296 for (
unsigned int f=0; f<builder_.faceSize(); ++f)
298 if (!builder_.testFaceBasis(f))
300 testBasisVal.resize(builder_.testFaceBasis(f)->size());
302 const Geometry &geometry = refElement.template geometry< 1 >( f );
303 const Dune::GeometryType subGeoType( geometry.type().id(), dimension-1 );
304 const FaceQuadrature &faceQuad = FaceQuadratureRules::rule( subGeoType, 2*order_+2 );
306 const unsigned int quadratureSize = faceQuad.size();
307 for(
unsigned int qi = 0; qi < quadratureSize; ++qi )
310 builder_.testFaceBasis(f)->template evaluate<0>(faceQuad[qi].position(),testBasisVal);
312 testBasisVal[0] = 1.;
313 fillBnd( row, testBasisVal,
314 func.evaluate( geometry.global( faceQuad[qi].position() ) ),
315 builder_.normal(f), faceQuad[qi].weight(),
319 row += builder_.testFaceBasis(f)->size();
322 if (builder_.testBasis())
324 testBasisVal.resize(builder_.testBasis()->size());
326 typedef Dune::QuadratureRule<Field, dimension> Quadrature;
327 typedef Dune::QuadratureRules<Field, dimension> QuadratureRules;
328 const Quadrature &elemQuad = QuadratureRules::rule( geoType, 2*order_+1 );
330 const unsigned int quadratureSize = elemQuad.size();
331 for(
unsigned int qi = 0; qi < quadratureSize; ++qi )
333 builder_.testBasis()->template evaluate<0>(elemQuad[qi].position(),testBasisVal);
334 fillInterior( row, testBasisVal,
335 func.evaluate(elemQuad[qi].position()),
336 elemQuad[qi].weight(),
340 row += builder_.testBasis()->size()*dimension;
347 template <
class MVal,
class RTVal,
class Matrix>
348 void fillBnd (
unsigned int startRow,
351 const FieldVector<Field,dimension> &normal,
353 Matrix &matrix)
const 355 const unsigned int endRow = startRow+mVal.size();
356 typename RTVal::const_iterator rtiter = rtVal.begin();
357 for (
unsigned int col = 0; col < rtVal.size() ; ++rtiter,++col)
359 Field cFactor = (*rtiter)*normal;
360 typename MVal::const_iterator miter = mVal.begin();
361 for (
unsigned int row = startRow;
362 row!=endRow; ++miter, ++row )
364 matrix.add(row,col, (weight*cFactor)*(*miter) );
366 assert( miter == mVal.end() );
369 template <
class MVal,
class RTVal,
class Matrix>
370 void fillInterior (
unsigned int startRow,
374 Matrix &matrix)
const 376 const unsigned int endRow = startRow+mVal.size()*dimension;
377 typename RTVal::const_iterator rtiter = rtVal.begin();
378 for (
unsigned int col = 0; col < rtVal.size() ; ++rtiter,++col)
380 typename MVal::const_iterator miter = mVal.begin();
381 for (
unsigned int row = startRow;
382 row!=endRow; ++miter,row+=dimension )
384 for (
unsigned int i=0; i<dimension; ++i)
386 matrix.add(row+i,col, (weight*(*miter))*(*rtiter)[i] );
389 assert( miter == mVal.end() );
398 template <
unsigned int dim,
class F >
401 static const unsigned int dimension = dim;
406 template <
unsigned int dim,
class Field >
408 public TopologyFactory< RaviartThomasL2InterpolationFactoryTraits<dim,Field> >
414 template <
class Topology>
417 if ( !supports<Topology>(key) )
419 NonConstObject *interpol =
new NonConstObject();
420 interpol->template build<Topology>(key);
423 template<
class Topology >
432 #endif // #ifndef DUNE_LOCALFUNCTIONS_RAVIARTTHOMAS_RAVIARTTHOMASSIMPLEX_RAVIARTTHOMASSIMPLEXINTERPOLATION_HH
Definition: tensor.hh:165
An L2-based interpolation for Raviart Thomas.
Definition: raviartthomassimplexinterpolation.hh:206
Definition: raviartthomassimplexinterpolation.hh:120
Definition: polynomialbasis.hh:62
const LocalKey & localKey(const unsigned int i) const
Definition: raviartthomassimplexinterpolation.hh:49
static bool supports(const typename Traits::Key &key)
Definition: raviartthomassimplexinterpolation.hh:102
unsigned int Key
Definition: raviartthomassimplexinterpolation.hh:402
LocalCoefficientsContainer(const Setter &setter)
Definition: raviartthomassimplexinterpolation.hh:44
TestBasisFactory::Object TestBasis
Definition: raviartthomassimplexinterpolation.hh:124
RaviartThomasL2InterpolationFactory< dim, F > Factory
Definition: raviartthomassimplexinterpolation.hh:404
~RTL2InterpolationBuilder()
Definition: raviartthomassimplexinterpolation.hh:134
void interpolate(const Function &function, std::vector< Fy > &coefficients) const
Definition: raviartthomassimplexinterpolation.hh:221
TestFaceBasis * testFaceBasis(unsigned int f) const
Definition: raviartthomassimplexinterpolation.hh:150
void interpolate(const Basis &basis, Matrix &matrix) const
Definition: raviartthomassimplexinterpolation.hh:228
const Normal & normal(unsigned int f) const
Definition: raviartthomassimplexinterpolation.hh:152
RaviartThomasL2Interpolation()
Definition: raviartthomassimplexinterpolation.hh:215
void build(unsigned int order)
Definition: raviartthomassimplexinterpolation.hh:244
Definition: brezzidouglasmarini1cube2dlocalbasis.hh:15
OrthonormalBasisFactory< dimension-1, Field > TestFaceBasisFactory
Definition: raviartthomassimplexinterpolation.hh:125
Definition: interpolationhelper.hh:17
unsigned int faceSize() const
Definition: raviartthomassimplexinterpolation.hh:147
void setLocalKeys(std::vector< LocalKey > &keys) const
Definition: raviartthomassimplexinterpolation.hh:256
unsigned int order() const
Definition: raviartthomassimplexinterpolation.hh:235
unsigned int Key
Definition: raviartthomassimplexinterpolation.hh:74
Traits::Object Object
Definition: raviartthomassimplexinterpolation.hh:412
static Traits::Object * createObject(const typename Traits::Key &key)
Definition: raviartthomassimplexinterpolation.hh:415
std::remove_const< Object >::type NonConstObject
Definition: raviartthomassimplexinterpolation.hh:413
FieldVector< Field, dimension > Normal
Definition: raviartthomassimplexinterpolation.hh:127
RTL2InterpolationBuilder< dim, Field > Builder
Definition: raviartthomassimplexinterpolation.hh:411
Definition: raviartthomassimplexinterpolation.hh:31
TestBasis * testBasis() const
Definition: raviartthomassimplexinterpolation.hh:149
unsigned int order() const
Definition: raviartthomassimplexinterpolation.hh:145
const LocalCoefficientsContainer Object
Definition: raviartthomassimplexinterpolation.hh:73
static bool supports(const typename Traits::Key &key)
Definition: raviartthomassimplexinterpolation.hh:424
GeometryType type() const
Definition: raviartthomassimplexinterpolation.hh:143
Definition: interpolationhelper.hh:20
OrthonormalBasisFactory< dimension, Field > TestBasisFactory
Definition: raviartthomassimplexinterpolation.hh:123
RaviartThomasCoefficientsFactoryTraits< dim > Traits
Definition: raviartthomassimplexinterpolation.hh:87
unsigned int size() const
Definition: raviartthomassimplexinterpolation.hh:239
unsigned int size() const
Definition: raviartthomassimplexinterpolation.hh:55
Definition: raviartthomassimplexinterpolation.hh:70
RTL2InterpolationBuilder< dimension, Field > Builder
Definition: raviartthomassimplexinterpolation.hh:214
Definition: raviartthomassimplexinterpolation.hh:28
Definition: interpolationhelper.hh:73
void build(unsigned int order)
Definition: raviartthomassimplexinterpolation.hh:155
F Field
Definition: raviartthomassimplexinterpolation.hh:213
Definition: raviartthomassimplexinterpolation.hh:399
Definition: raviartthomassimplexinterpolation.hh:38
TestFaceBasisFactory::Object TestFaceBasis
Definition: raviartthomassimplexinterpolation.hh:126
static Traits::Object * createObject(const typename Traits::Key &key)
Definition: raviartthomassimplexinterpolation.hh:90
void interpolate(typename Base::template Helper< Func, Container, type > &func) const
Definition: raviartthomassimplexinterpolation.hh:274
Definition: orthonormalbasis.hh:19
RaviartThomasCoefficientsFactory< dim > Factory
Definition: raviartthomassimplexinterpolation.hh:75
Describe position of one degree of freedom.
Definition: localkey.hh:20
RaviartThomasL2InterpolationFactoryTraits< dim, Field > Traits
Definition: raviartthomassimplexinterpolation.hh:410
const RaviartThomasL2Interpolation< dim, F > Object
Definition: raviartthomassimplexinterpolation.hh:403
unsigned int topologyId() const
Definition: raviartthomassimplexinterpolation.hh:141