dune-fem 2.8.0
Loading...
Searching...
No Matches
integrands.hh
Go to the documentation of this file.
1#ifndef DUNE_FEM_SCHEMES_INTEGRANDS_HH
2#define DUNE_FEM_SCHEMES_INTEGRANDS_HH
3
4#include <cassert>
5
6#include <algorithm>
7#include <functional>
8#include <tuple>
9#include <type_traits>
10#include <utility>
11
12#include <dune/common/ftraits.hh>
13
14#include <dune/fempy/quadrature/cachingpoint.hh>
15#include <dune/fempy/quadrature/elementpoint.hh>
18
19namespace Dune
20{
21
22 namespace Fem
23 {
24
25 // IntegrandsTraits
26 // ----------------
27
28 namespace Impl
29 {
30
31 namespace IntegrandsTraits
32 {
33
34 template< class Integrands >
35 using DomainValueType = typename std::decay_t< decltype( std::ref( std::declval< const Integrands & >() ).get() ) >::DomainValueType;
36
37 template< class Integrands >
38 using RangeValueType = typename std::decay_t< decltype( std::ref( std::declval< const Integrands & >() ).get() ) >::RangeValueType;
39
40 template< class Integrands >
41 using GridPartType = typename std::decay_t< decltype( std::ref( std::declval< const Integrands & >() ).get() ) >::GridPartType;
42
43
44 template< class Integrands >
45 using EntityType = typename GridPartType< Integrands >::template Codim< 0 >::EntityType;
46
47 template< class Integrands >
48 using IntersectionType = typename GridPartType< Integrands >::IntersectionType;
49
50
51 template< class Integrands >
52 using InteriorQuadratureType = CachingQuadrature< GridPartType< Integrands >, 0 >;
53
54 template< class Integrands >
55 using SurfaceQuadratureType = CachingQuadrature< GridPartType< Integrands >, 1 >;
56
57
58 template< class Integrands >
59 using InteriorQuadraturePointType = QuadraturePointWrapper< InteriorQuadratureType< Integrands > >;
60
61 template< class Integrands >
62 using SurfaceQuadraturePointType = QuadraturePointWrapper< SurfaceQuadratureType< Integrands > >;
63
64
65 template< class Integrands >
66 std::true_type interior ( const Integrands &, decltype( std::declval< const Integrands & >().interior( std::declval< const InteriorQuadraturePointType< Integrands > & >(), std::declval< const DomainValueType< Integrands > & >() ) ) * = nullptr );
67
68 std::false_type interior ( ... );
69
70 template< class Integrands, std::enable_if_t< std::is_same< decltype( std::declval< const Integrands & >().hasInterior() ), bool >::value, int > = 0 >
71 std::true_type hasInterior ( const Integrands & );
72
73 std::false_type hasInterior ( ... );
74
75
76 template< class Integrands >
77 std::true_type boundary ( const Integrands &, decltype( std::declval< const Integrands & >().boundary( std::declval< const SurfaceQuadraturePointType< Integrands > & >(), std::declval< const DomainValueType< Integrands > & >() ) ) * = nullptr );
78
79 std::false_type boundary ( ... );
80
81 template< class Integrands, std::enable_if_t< std::is_same< decltype( std::declval< const Integrands & >().hasBoundary() ), bool >::value, int > = 0 >
82 std::true_type hasBoundary ( const Integrands & );
83
84 std::false_type hasBoundary ( ... );
85
86
87 template< class Integrands >
88 std::true_type skeleton ( const Integrands &, decltype( std::declval< const Integrands & >().skeleton( std::declval< const SurfaceQuadraturePointType< Integrands > & >(), std::declval< const DomainValueType< Integrands > & >(), std::declval< const SurfaceQuadraturePointType< Integrands > & >(), std::declval< const DomainValueType< Integrands > & >() ) ) * = nullptr );
89
90 std::false_type skeleton ( ... );
91
92 template< class Integrands, std::enable_if_t< std::is_same< decltype( std::declval< const Integrands & >().hasSkeleton() ), bool >::value, int > = 0 >
93 std::true_type hasSkeleton ( const Integrands & );
94
95 std::false_type hasSkeleton ( ... );
96
97 template< class Integrands >
98 using Get = std::decay_t< decltype( std::ref( std::declval< const Integrands & >() ).get() ) >;
99
100 template< class Integrands >
101 using RRangeType = typename Get<Integrands>::RRangeType;
102 template< class Integrands >
103 using DirichletComponentType = typename Get<Integrands>::DirichletComponentType;
104 } // namespace IntegrandsTraits
105
106 } // namespace Impl
107
108 template< class Integrands >
110 {
111 typedef Impl::IntegrandsTraits::DomainValueType< Integrands > DomainValueType;
112 typedef Impl::IntegrandsTraits::RangeValueType< Integrands > RangeValueType;
113 typedef Impl::IntegrandsTraits::GridPartType< Integrands > GridPartType;
114
115 typedef Impl::IntegrandsTraits::EntityType< Integrands > EntityType;
116 typedef Impl::IntegrandsTraits::IntersectionType< Integrands > IntersectionType;
117
118 static const bool interior = decltype( Impl::IntegrandsTraits::interior( std::declval< const Integrands & >() ) )::value;
119 static const bool hasInterior = decltype( Impl::IntegrandsTraits::hasInterior( std::declval< const Integrands & >() ) )::value;
120
121 static const bool boundary = decltype( Impl::IntegrandsTraits::boundary( std::declval< const Integrands & >() ) )::value;
122 static const bool hasBoundary = decltype( Impl::IntegrandsTraits::hasBoundary( std::declval< const Integrands & >() ) )::value;
123
124 static const bool skeleton = decltype( Impl::IntegrandsTraits::skeleton( std::declval< const Integrands & >() ) )::value;
125 static const bool hasSkeleton = decltype( Impl::IntegrandsTraits::hasSkeleton( std::declval< const Integrands & >() ) )::value;
126
127 static_assert( (!hasInterior || interior), "Existence of method 'hasInterior' implies existence of method interior." );
128 static_assert( (!hasBoundary || boundary), "Existence of method 'hasBoundary' implies existence of method boundary." );
129 static_assert( (!hasSkeleton || skeleton), "Existence of method 'hasSkeleton' implies existence of method skeleton." );
130
131 static const bool isFull = hasInterior && hasBoundary && hasSkeleton;
132
133 typedef Impl::IntegrandsTraits::RRangeType< Integrands > RRangeType;
134 typedef Impl::IntegrandsTraits::DirichletComponentType< Integrands > DirichletComponentType;
135 };
136
137
138
139 // FullIntegrands
140 // --------------
141
142 template< class Integrands >
144 {
148
151
152 private:
153 template< class T, std::enable_if_t< IntegrandsTraits< T >::hasInterior, int > = 0 >
154 static bool hasInterior ( const T &integrands )
155 {
156 return integrands.hasInterior();
157 }
158
159 template< class T, std::enable_if_t< !IntegrandsTraits< T >::hasInterior, int > = 0 >
160 static bool hasInterior ( const T &integrands )
161 {
163 }
164
165 template< class T, class Point, std::enable_if_t< IntegrandsTraits< T >::interior, int > = 0 >
166 static RangeValueType interior ( const T &integrands, const Point &x, const DomainValueType &u )
167 {
168 return integrands.interior( x, u );
169 }
170
171 template< class T, class Point, std::enable_if_t< !IntegrandsTraits< T >::interior, int > = 0 >
172 static RangeValueType interior ( const T &integrands, const Point &x, const DomainValueType &u )
173 {
174 return RangeValueType();
175 }
176
177 template< class T, class Point, std::enable_if_t< IntegrandsTraits< T >::interior, int > = 0 >
178 static auto linearizedInterior ( const T &integrands, const Point &x, const DomainValueType &u )
179 {
180 return integrands.linearizedInterior( x, u );
181 }
182
183 template< class T, class Point, std::enable_if_t< !IntegrandsTraits< T >::interior, int > = 0 >
184 static auto linearizedInterior ( const T &integrands, const Point &x, const DomainValueType &u )
185 {
186 return [] ( const DomainValueType & ) { return RangeValueType(); };
187 }
188
189 template< class T, std::enable_if_t< IntegrandsTraits< T >::hasBoundary, int > = 0 >
190 static bool hasBoundary ( const T &integrands )
191 {
192 return integrands.hasBoundary();
193 }
194
195 template< class T, std::enable_if_t< !IntegrandsTraits< T >::hasBoundary, int > = 0 >
196 static bool hasBoundary ( const T &integrands )
197 {
199 }
200
201 template< class T, class Point, std::enable_if_t< IntegrandsTraits< T >::boundary, int > = 0 >
202 static RangeValueType boundary ( const T &integrands, const Point &x, const DomainValueType &u )
203 {
204 return integrands.boundary( x, u );
205 }
206
207 template< class T, class Point, std::enable_if_t< !IntegrandsTraits< T >::boundary, int > = 0 >
208 static RangeValueType boundary ( const T &integrands, const Point &x, const DomainValueType &u )
209 {
210 return RangeValueType();
211 }
212
213 template< class T, class Point, std::enable_if_t< IntegrandsTraits< T >::boundary, int > = 0 >
214 static auto linearizedBoundary ( const T &integrands, const Point &x, const DomainValueType &u )
215 {
216 return integrands.linearizedBoundary( x, u );
217 }
218
219 template< class T, class Point, std::enable_if_t< !IntegrandsTraits< T >::boundary, int > = 0 >
220 static auto linearizedBoundary ( const T &integrands, const Point &x, const DomainValueType &u )
221 {
222 return [] ( const DomainValueType & ) { return RangeValueType(); };
223 }
224
225 template< class T, std::enable_if_t< IntegrandsTraits< T >::hasSkeleton, int > = 0 >
226 static bool hasSkeleton ( const T &integrands )
227 {
228 return integrands.hasSkeleton();
229 }
230
231 template< class T, std::enable_if_t< !IntegrandsTraits< T >::hasSkeleton, int > = 0 >
232 static bool hasSkeleton ( const T &integrands )
233 {
235 }
236
237 template< class T, class Point, std::enable_if_t< IntegrandsTraits< T >::skeleton, int > = 0 >
238 static std::pair< RangeValueType, RangeValueType > skeleton ( const T &integrands, const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut )
239 {
240 return integrands.skeleton( xIn, uIn, xOut, uOut );
241 }
242
243 template< class T, class Point, std::enable_if_t< !IntegrandsTraits< T >::skeleton, int > = 0 >
244 static std::pair< RangeValueType, RangeValueType > skeleton ( const T &integrands, const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut )
245 {
246 return std::make_pair( RangeValueType(), RangeValueType() );
247 }
248
249 template< class T, class Point, std::enable_if_t< IntegrandsTraits< T >::skeleton, int > = 0 >
250 static auto linearizedSkeleton ( const T &integrands, const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut )
251 {
252 return integrands.linearizedSkeleton( xIn, uIn, xOut, uOut );
253 }
254
255 template< class T, class Point, std::enable_if_t< !IntegrandsTraits< T >::skeleton, int > = 0 >
256 static auto linearizedSkeleton ( const T &integrands, const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut )
257 {
258 auto zero = [] ( const DomainValueType & ) { return std::make_pair( RangeValueType(), RangeValueType() ); };
259 return std::make_pair( zero, zero );
260 }
261
262 public:
263 template< class... Args >
264 explicit FullIntegrands ( Args &&... args )
265 : integrands_( std::forward< Args >( args )... ),
266 rInt_( std::ref( integrands_ ).get() )
267 {
268 }
269
270 bool init ( const EntityType &entity ) { return integrands().init( entity ); }
271 bool init ( const IntersectionType &intersection ) { return integrands().init( intersection ); }
272 void unbind ( ) { integrands().unbind( ); }
273
274 bool hasInterior () const { return hasInterior( integrands() ); }
275
276 template< class Point >
277 RangeValueType interior ( const Point &x, const DomainValueType &u ) const
278 {
279 return interior( integrands(), x, u );
280 }
281
282 template< class Point >
283 auto linearizedInterior ( const Point &x, const DomainValueType &u ) const
284 {
285 return linearizedInterior( integrands(), x, u );
286 }
287
288 bool hasBoundary () const { return hasBoundary( integrands() ); }
289
290 template< class Point >
291 RangeValueType boundary ( const Point &x, const DomainValueType &u ) const
292 {
293 return boundary( integrands(), x, u );
294 }
295
296 template< class Point >
297 auto linearizedBoundary ( const Point &x, const DomainValueType &u ) const
298 {
299 return linearizedBoundary( integrands(), x, u );
300 }
301
302 bool hasSkeleton () const { return hasSkeleton( integrands() ); }
303
304 template< class Point >
305 std::pair< RangeValueType, RangeValueType > skeleton ( const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut ) const
306 {
307 return skeleton( integrands(), xIn, uIn, xOut, uOut );
308 }
309
310 template< class Point >
311 auto linearizedSkeleton ( const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut ) const
312 {
313 return linearizedSkeleton( integrands(), xIn, uIn, xOut, uOut );
314 }
315
316 protected:
317 typedef typename Integrands::type RealIntegrands;
318 //decltype( auto ) integrands () { return std::ref( integrands_ ).get(); }
319 //decltype( auto ) integrands () const { return std::ref( integrands_ ).get(); }
321 const RealIntegrands& integrands () const { return rInt_; }
322
323
324
325 Integrands integrands_;
327
328 public:
332 {
333 return integrands().hasDirichletBoundary();
334 }
335 bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const
336 {
337 return integrands().isDirichletIntersection(inter,dirichletComponent);
338 }
339 template <class Point>
340 void dirichlet( int bndId, const Point &x, RRangeType &value) const
341 {
342 return integrands().dirichlet(bndId,x,value);
343 }
344 };
345
346
347
348 // VirtualizedIntegrands
349 // ---------------------
350
351 namespace detail
352 {
353 template <class T>
354 struct GetDimRange;
355 template <class FT,int r>
356 struct GetDimRange<Dune::FieldVector<FT,r>>
357 {
358 typedef Dune::FieldVector<FT,r> type;
359 static const int value = r;
360 };
361 template <class FT,int r,int c>
362 struct GetDimRange<Dune::FieldMatrix<FT,r,c>>
363 {
364 typedef Dune::FieldVector<FT,r> type;
365 static const int value = r;
366 };
367 template <class FT,int r,int c>
368 struct GetDimRange<Dune::Fem::ExplicitFieldVector<Dune::FieldMatrix<FT,c,c>,r>>
369 {
370 typedef Dune::FieldVector<FT,r> type;
371 static const int value = r;
372 };
373 }
374
375 template< class GridPart, class DomainValue, class RangeValue = DomainValue >
377 {
379
380 public:
381 typedef GridPart GridPartType;
382 typedef DomainValue DomainValueType;
383 typedef RangeValue RangeValueType;
384
385 typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
386 typedef typename GridPartType::IntersectionType IntersectionType;
387
388 using RRangeType = typename detail::GetDimRange<std::tuple_element_t<0,RangeValueType>>::type;
389 typedef std::array<int,RRangeType::dimension> DirichletComponentType;
390 typedef typename EntityType::Geometry::LocalCoordinate DomainType;
391
392 private:
393 typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
394
395 typedef FemPy::CachingPoint< LocalCoordinateType, 0 > InteriorCachingPointType;
396 typedef FemPy::ElementPoint< LocalCoordinateType, 0 > InteriorElementPointType;
397 typedef FemPy::CachingPoint< LocalCoordinateType, 1 > SurfaceCachingPointType;
398 typedef FemPy::ElementPoint< LocalCoordinateType, 1 > SurfaceElementPointType;
399
400 template< class QP >
401 static Fem::QuadraturePointWrapper< QP > asQP ( const QP &qp )
402 {
403 return static_cast< Fem::QuadraturePointWrapper< QP > >( qp );
404 }
405
406 template< class R >
407 using Linearization = std::function< R( const DomainValueType & ) >;
408
409 template< class R >
410 using LinearizationPair = std::pair< Linearization< std::pair< R, R > >, Linearization< std::pair< R, R > > >;
411
412 struct Interface
413 {
414 virtual ~Interface () {}
415 virtual Interface *clone () const = 0;
416
417 virtual bool init ( const EntityType &entity ) = 0;
418 virtual bool init ( const IntersectionType &intersection ) = 0;
419 virtual void unbind ( ) = 0;
420
421 virtual bool hasInterior () const = 0;
422 virtual RangeValueType interior ( const InteriorCachingPointType &x, const DomainValueType &u ) const = 0;
423 virtual RangeValueType interior ( const InteriorElementPointType &x, const DomainValueType &u ) const = 0;
424 virtual Linearization< RangeValueType > linearizedInterior ( const InteriorCachingPointType &x, const DomainValueType &u ) const = 0;
425 virtual Linearization< RangeValueType > linearizedInterior ( const InteriorElementPointType &x, const DomainValueType &u ) const = 0;
426
427 virtual bool hasBoundary () const = 0;
428 virtual RangeValueType boundary ( const SurfaceCachingPointType &x, const DomainValueType &u ) const = 0;
429 virtual RangeValueType boundary ( const SurfaceElementPointType &x, const DomainValueType &u ) const = 0;
430 virtual Linearization< RangeValueType > linearizedBoundary ( const SurfaceCachingPointType &x, const DomainValueType &u ) const = 0;
431 virtual Linearization< RangeValueType > linearizedBoundary ( const SurfaceElementPointType &x, const DomainValueType &u ) const = 0;
432
433 virtual bool hasSkeleton () const = 0;
434 virtual std::pair< RangeValueType, RangeValueType > skeleton ( const SurfaceCachingPointType &xIn, const DomainValueType &uIn, const SurfaceCachingPointType &xOut, const DomainValueType &uOut ) const = 0;
435 virtual std::pair< RangeValueType, RangeValueType > skeleton ( const SurfaceElementPointType &xIn, const DomainValueType &uIn, const SurfaceElementPointType &xOut, const DomainValueType &uOut ) const = 0;
436 virtual LinearizationPair< RangeValueType > linearizedSkeleton ( const SurfaceCachingPointType &xIn, const DomainValueType &uIn, const SurfaceCachingPointType &xOut, const DomainValueType &uOut ) const = 0;
437 virtual LinearizationPair< RangeValueType > linearizedSkeleton ( const SurfaceElementPointType &xIn, const DomainValueType &uIn, const SurfaceElementPointType &xOut, const DomainValueType &uOut ) const = 0;
438
439 virtual bool hasDirichletBoundary () const = 0;
440 virtual bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const = 0;
441 virtual void dirichlet( int bndId, const DomainType &x,RRangeType &value) const = 0;
442 };
443
444 template< class Impl >
445 struct Implementation final
446 : public Interface
447 {
448 Implementation ( Impl impl ) : impl_( std::move( impl ) )
449 {
450 }
451 virtual Interface *clone () const override { return new Implementation( *this ); }
452
453 virtual bool init ( const EntityType &entity ) override { return impl().init( entity ); }
454 virtual bool init ( const IntersectionType &intersection ) override { return impl().init( intersection ); }
455 virtual void unbind ( ) override { impl().unbind( ); }
456
457 virtual bool hasInterior () const override { return impl().hasInterior(); }
458 virtual RangeValueType interior ( const InteriorCachingPointType &x, const DomainValueType &u ) const override { return impl().interior( asQP( x ), u ); }
459 virtual RangeValueType interior ( const InteriorElementPointType &x, const DomainValueType &u ) const override { return impl().interior( asQP( x ), u ); }
460 virtual Linearization< RangeValueType > linearizedInterior ( const InteriorCachingPointType &x, const DomainValueType &u ) const override { return impl().linearizedInterior( asQP( x ), u ); }
461 virtual Linearization< RangeValueType > linearizedInterior ( const InteriorElementPointType &x, const DomainValueType &u ) const override { return impl().linearizedInterior( asQP( x ), u ); }
462
463 virtual bool hasBoundary () const override { return impl().hasBoundary(); }
464 virtual RangeValueType boundary ( const SurfaceCachingPointType &x, const DomainValueType &u ) const override { return impl().boundary( asQP( x ), u ); }
465 virtual RangeValueType boundary ( const SurfaceElementPointType &x, const DomainValueType &u ) const override { return impl().boundary( asQP( x ), u ); }
466 virtual Linearization< RangeValueType > linearizedBoundary ( const SurfaceCachingPointType &x, const DomainValueType &u ) const override { return impl().linearizedBoundary( asQP( x ), u ); }
467 virtual Linearization< RangeValueType > linearizedBoundary ( const SurfaceElementPointType &x, const DomainValueType &u ) const override { return impl().linearizedBoundary( asQP( x ), u ); }
468
469 virtual bool hasSkeleton () const override { return impl().hasSkeleton(); }
470 virtual std::pair< RangeValueType, RangeValueType > skeleton ( const SurfaceCachingPointType &xIn, const DomainValueType &uIn, const SurfaceCachingPointType &xOut, const DomainValueType &uOut ) const override { return impl().skeleton( asQP( xIn ), uIn, asQP( xOut ), uOut ); }
471 virtual std::pair< RangeValueType, RangeValueType > skeleton ( const SurfaceElementPointType &xIn, const DomainValueType &uIn, const SurfaceElementPointType &xOut, const DomainValueType &uOut ) const override { return impl().skeleton( asQP( xIn ), uIn, asQP( xOut ), uOut ); }
472 virtual LinearizationPair< RangeValueType > linearizedSkeleton ( const SurfaceCachingPointType &xIn, const DomainValueType &uIn, const SurfaceCachingPointType &xOut, const DomainValueType &uOut ) const override { return impl().linearizedSkeleton( asQP( xIn ), uIn, asQP( xOut ), uOut ); }
473 virtual LinearizationPair< RangeValueType > linearizedSkeleton ( const SurfaceElementPointType &xIn, const DomainValueType &uIn, const SurfaceElementPointType &xOut, const DomainValueType &uOut ) const override { return impl().linearizedSkeleton( asQP( xIn ), uIn, asQP( xOut ), uOut ); }
474
475 virtual bool hasDirichletBoundary () const override { return impl().hasDirichletBoundary(); }
476 virtual bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const override { return impl().isDirichletIntersection(inter,dirichletComponent); }
477 virtual void dirichlet( int bndId, const DomainType &x,RRangeType &value) const override { impl().dirichlet(bndId,x,value); }
478
479 private:
480 const auto &impl () const { return std::cref( impl_ ).get(); }
481 auto &impl () { return std::ref( impl_ ).get(); }
482
483 Impl impl_;
484 };
485
486 template< class Integrands >
487 using isVirtualized = std::is_same< std::decay_t< decltype( std::ref( std::declval< Integrands & >() ).get() ) >, This >;
488
489 public:
490 template< class Integrands, std::enable_if_t< IntegrandsTraits< std::decay_t< Integrands > >::isFull && !isVirtualized< Integrands >::value, int > = 0 >
491 explicit VirtualizedIntegrands ( Integrands integrands )
492 : impl_( new Implementation< Integrands >( std::move( integrands ) ) )
493 {}
494
495 template< class Integrands, std::enable_if_t< !IntegrandsTraits< Integrands >::isFull && !isVirtualized< Integrands >::value, int > = 0 >
496 explicit VirtualizedIntegrands ( Integrands integrands )
497 : VirtualizedIntegrands( FullIntegrands< std::decay_t< Integrands > >( std::move( integrands ) ) )
498 {}
499
500 VirtualizedIntegrands ( const This &other ) : impl_( other ? other.impl().clone() : nullptr )
501 {}
502
503 VirtualizedIntegrands ( This && ) = default;
504
506 {
507 impl_.reset( other ? other.impl().clone() : nullptr );
508 }
509
511
512 explicit operator bool () const { return static_cast< bool >( impl_ ); }
513
514 bool init ( const EntityType &entity ) { return impl().init( entity ); }
515 bool init ( const IntersectionType &intersection ) { return impl().init( intersection ); }
516 void unbind ( ) { impl().unbind( ); }
517
518 bool hasInterior () const { return impl().hasInterior(); }
519
520 template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
522 {
523 return impl().interior( InteriorCachingPointType( x ), u );
524 }
525
526 template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
528 {
529 return impl().interior( InteriorElementPointType( x ), u );
530 }
531
532 template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
534 {
535 return impl().linearizedInterior( InteriorCachingPointType( x ), u );
536 }
537
538 template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
540 {
541 return impl().linearizedInterior( InteriorElementPointType( x ), u );
542 }
543
544 bool hasBoundary () const { return impl().hasBoundary(); }
545
546 template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
548 {
549 return impl().boundary( SurfaceCachingPointType( x ), u );
550 }
551
552 template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
554 {
555 return impl().boundary( SurfaceElementPointType( x ), u );
556 }
557
558 template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
560 {
561 return impl().linearizedBoundary( SurfaceCachingPointType( x ), u );
562 }
563
564 template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
566 {
567 return impl().linearizedBoundary( SurfaceElementPointType( x ), u );
568 }
569
570 bool hasSkeleton () const { return impl().hasSkeleton(); }
571
572 template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
573 std::pair< RangeValueType, RangeValueType > skeleton ( const Fem::QuadraturePointWrapper< Quadrature > &xIn, const DomainValueType &uIn, const Fem::QuadraturePointWrapper< Quadrature > &xOut, const DomainValueType &uOut ) const
574 {
575 return impl().skeleton( SurfaceCachingPointType( xIn ), uIn, SurfaceCachingPointType( xOut ), uOut );
576 }
577
578 template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
579 std::pair< RangeValueType, RangeValueType > skeleton ( const Fem::QuadraturePointWrapper< Quadrature > &xIn, const DomainValueType &uIn, const Fem::QuadraturePointWrapper< Quadrature > &xOut, const DomainValueType &uOut ) const
580 {
581 return impl().skeleton( SurfaceElementPointType( xIn ), uIn, SurfaceElementPointType( xOut ), uOut );
582 }
583
584 template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
586 {
587 return impl().linearizedSkeleton( SurfaceCachingPointType( xIn ), uIn, SurfaceCachingPointType( xOut ), uOut );
588 }
589
590 template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
592 {
593 return impl().linearizedSkeleton( SurfaceElementPointType( xIn ), uIn, SurfaceElementPointType( xOut ), uOut );
594 }
595
597 {
598 return impl().hasDirichletBoundary();
599 }
600 bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const
601 {
602 return impl().isDirichletIntersection(inter,dirichletComponent);
603 }
604 void dirichlet( int bndId, const DomainType &x,RRangeType &value) const
605 {
606 return impl().dirichlet(bndId,x,value);
607 }
608
609 private:
610 const Interface &impl () const { assert( impl_ ); return *impl_; }
611 Interface &impl () { assert( impl_ ); return *impl_; }
612
613 std::unique_ptr< Interface > impl_;
614 };
615
616
617
618 // ConservationLawModelIntegrands
619 // ------------------------------
620
621 template< class Model >
623 {
624 typedef Model ModelType;
625 typedef typename Model::GridPartType GridPartType;
626
627 typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
628 typedef typename GridPartType::IntersectionType IntersectionType;
629
630 typedef typename Model::RangeType RangeType;
631 typedef typename Model::JacobianRangeType JacobianRangeType;
632
633 typedef std::tuple< RangeType, JacobianRangeType > DomainValueType;
634 typedef std::tuple< RangeType, JacobianRangeType > RangeValueType;
635
636 explicit ConservationLawModelIntegrands ( const Model &model ) : model_( &model ) {}
637
638 bool init ( const EntityType &entity ) { return model().init( entity ); }
639
640 bool init ( const IntersectionType &intersection )
641 {
642 return (intersection.boundary() && model().hasNeumanBoundary() && model().init( intersection.inside() ));
643 }
644
645 void unbind ( ) { model().unbind( ); }
646
647 template< class Point >
648 RangeValueType interior ( const Point &x, const DomainValueType &u ) const
649 {
650 RangeType source( 0 );
651 model().source( x, std::get< 0 >( u ), std::get< 1 >( u ), source );
652 JacobianRangeType dFlux( 0 );
653 model().flux( x, std::get< 0 >( u ), std::get< 1 >( u ), dFlux );
654 return std::make_tuple( source, dFlux );
655 }
656
657 template< class Point >
658 auto linearizedInterior ( const Point &x, const DomainValueType &u ) const
659 {
660 return [ this, x, u ] ( const DomainValueType &phi ) {
661 RangeType source( 0 );
662 model().linSource( std::get< 0 >( u ), std::get< 1 >( u ), x, std::get< 0 >( phi ), std::get< 1 >( phi ), source );
663 JacobianRangeType dFlux( 0 );
664 model().linFlux( std::get< 0 >( u ), std::get< 1 >( u ), x, std::get< 0 >( phi ), std::get< 1 >( phi ), dFlux );
665 return std::make_tuple( source, dFlux );
666 };
667 }
668
669 template< class Point >
670 RangeValueType boundary ( const Point &x, const DomainValueType &u ) const
671 {
672 RangeType alpha( 0 );
673 model().alpha( x, std::get< 0 >( u ), alpha );
674 return std::make_tuple( alpha, 0 );
675 }
676
677 template< class Point >
678 auto linearizedBoundary ( const Point &x, const DomainValueType &u ) const
679 {
680 return [ this, x, u ] ( const DomainValueType &phi ) {
681 RangeType alpha( 0 );
682 model().linAlpha( std::get< 0 >( u ), x, std::get< 0 >( phi ), alpha );
683 return std::make_tuple( alpha, 0 );
684 };
685 }
686
687 const Model &model () const { return *model_; }
688
689 private:
690 const Model *model_ = nullptr;
691 };
692
693
694
695 // DGConservationLawModelIntegrands
696 // --------------------------------
697
698 template< class Model >
700 {
701 typedef Model ModelType;
702 typedef typename Model::GridPartType GridPartType;
703
704 typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
705 typedef typename GridPartType::IntersectionType IntersectionType;
706
707 typedef typename Model::RangeType RangeType;
708 typedef typename Model::JacobianRangeType JacobianRangeType;
709
710 typedef typename FieldTraits< RangeType >::field_type RangeFieldType;
711
712 typedef std::tuple< RangeType, JacobianRangeType > DomainValueType;
713 typedef std::tuple< RangeType, JacobianRangeType > RangeValueType;
714
716 : model_( &model ), penalty_( penalty )
717 {}
718
719 bool init ( const EntityType &entity )
720 {
721 intersection_ = nullptr;
722 return model().init( entity );
723 }
724
725 bool init ( const IntersectionType &intersection )
726 {
727 intersection_ = &intersection;
728 if( intersection.boundary() )
729 {
730 const EntityType inside = intersection.inside();
731 beta_ = penalty_ * intersection.geometry().volume() / inside.geometry().volume();
732 return (model().hasNeumanBoundary() && model().init( inside ));
733 }
734 else if( intersection.neighbor() )
735 {
736 const auto volIn = intersection.inside().geometry().volume();
737 const auto volOut = intersection.outside().geometry().volume();
738 beta_ = penalty_ * intersection.geometry().volume() / std::min( volIn, volOut );
739 return true;
740 }
741 }
742
743 void unbind ( )
744 {
745 model().unbind( );
746 }
747
748 template< class Point >
749 RangeValueType interior ( const Point &x, const DomainValueType &u ) const
750 {
751 RangeType source( 0 );
752 model().source( x, std::get< 0 >( u ), std::get< 1 >( u ), source );
753 JacobianRangeType dFlux( 0 );
754 model().flux( x, std::get< 0 >( u ), std::get< 1 >( u ), dFlux );
755 return RangeValueType( source, dFlux );
756 }
757
758 template< class Point >
759 auto linearizedInterior ( const Point &x, const DomainValueType &u ) const
760 {
761 return [ this, x, u ] ( const DomainValueType &phi ) {
762 RangeType source( 0 );
763 model().linSource( std::get< 0 >( u ), std::get< 1 >( u ), x, std::get< 0 >( phi ), std::get< 1 >( phi ), source );
764 JacobianRangeType dFlux( 0 );
765 model().linFlux( std::get< 0 >( u ), std::get< 1 >( u ), x, std::get< 0 >( phi ), std::get< 1 >( phi ), dFlux );
766 return RangeValueType( source, dFlux );
767 };
768 }
769
770 template< class Point >
771 RangeValueType boundary ( const Point &x, const DomainValueType &u ) const
772 {
773 RangeType alpha( 0 );
774 model().alpha( x, std::get< 0 >( u ), alpha );
775 return RangeValueType( alpha, 0 );
776 }
777
778 template< class Point >
779 auto linearizedBoundary ( const Point &x, const DomainValueType &u ) const
780 {
781 return [ this, x, u ] ( const DomainValueType &phi ) {
782 RangeType alpha( 0 );
783 model().linAlpha( std::get< 0 >( u ), x, std::get< 0 >( phi ), alpha );
784 return RangeValueType( alpha, 0 );
785 };
786 }
787
788 template< class Point >
789 std::pair< RangeValueType, RangeValueType > skeleton ( const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut ) const
790 {
791 const EntityType inside = intersection().inside();
792 const EntityType outside = intersection().outside();
793
794 const RangeFieldType half = RangeFieldType( 1 ) / RangeFieldType( 2 );
795 const auto normal = intersection().unitOuterNormal( xIn.localPosition() );
796
797 DomainValueType uJump( 0, 0 );
798 std::get< 0 >( uJump ) = std::get< 0 >( uOut ) - std::get< 0 >( uIn );
799 for( int i = 0; i < RangeType::dimension; ++i )
800 std::get< 1 >( uJump )[ i ].axpy( std::get< 0 >( uJump )[ i ], normal );
801
802 model().init( outside );
803 JacobianRangeType dFluxOut( 0 ), dFluxPrimeOut( 0 );
804 model().flux( xOut, std::get< 0 >( uOut ), std::get< 1 >( uJump ), dFluxPrimeOut );
805 model().flux( xOut, std::get< 0 >( uOut ), 0, dFluxOut );
806 dFluxPrimeOut -= dFluxOut;
807 dFluxOut = 0;
808 model().flux( xOut, std::get< 0 >( uOut ), std::get< 1 >( uOut ), dFluxOut );
809
810 model().init( inside );
811 JacobianRangeType dFluxIn( 0 ), dFluxPrimeIn( 0 );
812 model().flux( xIn, std::get< 0 >( uIn ), std::get< 1 >( uJump ), dFluxPrimeIn );
813 model().flux( xIn, std::get< 0 >( uIn ), 0, dFluxIn );
814 dFluxPrimeIn -= dFluxIn;
815 dFluxIn = 0;
816 model().flux( xIn, std::get< 0 >( uIn ), std::get< 1 >( uIn ), dFluxIn );
817
818 RangeType int0 = std::get< 0 >( uJump );
819 int0 *= beta_;
820 dFluxIn += dFluxOut;
821 dFluxIn.usmv( -half, normal, int0 );
822
823 dFluxPrimeIn *= -half;
824 dFluxPrimeOut *= -half;
825
826 return std::make_pair( RangeValueType( -int0, dFluxPrimeIn ), RangeValueType( int0, dFluxPrimeOut ) );
827 }
828
829 template< class Point >
830 auto linearizedSkeleton ( const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut ) const
831 {
832 const auto normal = intersection().unitOuterNormal( xIn.localPosition() );
833
834 DomainValueType uJump( 0, 0 );
835 std::get< 0 >( uJump ) = std::get< 0 >( uOut ) - std::get< 0 >( uIn );
836 for( int i = 0; i < RangeType::dimension; ++i )
837 std::get< 1 >( uJump )[ i ].axpy( std::get< 0 >( uJump )[ i ], normal );
838
839 auto intIn = [ this, xIn, uIn, xOut, uOut, normal, uJump ] ( const DomainValueType &phiIn ) {
840 const EntityType inside = intersection().inside();
841 const EntityType outside = intersection().outside();
842
843 const RangeFieldType half = RangeFieldType( 1 ) / RangeFieldType( 2 );
844
845 DomainValueType phiJump( 0, 0 );
846 std::get< 0 >( phiJump ) -= std::get< 0 >( phiIn );
847 for( int i = 0; i < RangeType::dimension; ++i )
848 std::get< 1 >( phiJump )[ i ].axpy( std::get< 0 >( phiJump )[ i ], normal );
849
850 model().init( outside );
851 JacobianRangeType dFluxPrimeOut( 0 );
852 model().linFlux( std::get< 0 >( uOut ), std::get< 1 >( uJump ), xOut, 0, std::get< 1 >( phiJump ), dFluxPrimeOut );
853
854 model().init( inside );
855 JacobianRangeType dFluxIn( 0 ), dFluxPrimeIn( 0 );
856 model().linFlux( std::get< 0 >( uIn ), std::get< 1 >( uJump ), xIn, std::get< 0 >( phiIn ), std::get< 1 >( phiJump ), dFluxPrimeIn );
857 model().linFlux( std::get< 0 >( uIn ), 0, xIn, std::get< 0 >( phiIn ), 0, dFluxIn );
858 dFluxPrimeIn -= dFluxIn;
859 dFluxIn = 0;
860 model().linFlux( std::get< 0 >( uIn ), std::get< 1 >( uIn ), xIn, std::get< 0 >( phiIn ), std::get< 1 >( phiIn ), dFluxIn );
861
862 RangeType int0 = std::get< 0 >( phiJump );
863 int0 *= beta_;
864 dFluxIn.usmv( -half, normal, int0 );
865
866 dFluxPrimeIn *= -half;
867 dFluxPrimeOut *= -half;
868
869 return std::make_pair( RangeValueType( -int0, dFluxPrimeIn ), RangeValueType( int0, dFluxPrimeOut ) );
870 };
871
872 auto intOut = [ this, xIn, uIn, xOut, uOut, normal, uJump ] ( const DomainValueType &phiOut ) {
873 const EntityType inside = intersection().inside();
874 const EntityType outside = intersection().outside();
875
876 const RangeFieldType half = RangeFieldType( 1 ) / RangeFieldType( 2 );
877
878 DomainValueType phiJump( 0, 0 );
879 std::get< 0 >( phiJump ) = std::get< 0 >( phiOut );
880 for( int i = 0; i < RangeType::dimension; ++i )
881 std::get< 1 >( phiJump )[ i ].axpy( std::get< 0 >( phiJump )[ i ], normal );
882
883 model().init( outside );
884 JacobianRangeType dFluxOut( 0 ), dFluxPrimeOut( 0 );
885 model().linFlux( std::get< 0 >( uOut ), std::get< 1 >( uJump ), xOut, std::get< 0 >( phiOut ), std::get< 1 >( phiJump ), dFluxPrimeOut );
886 model().linFlux( std::get< 0 >( uOut ), 0, xOut, std::get< 0 >( phiOut ), 0, dFluxOut );
887 dFluxPrimeOut -= dFluxOut;
888 dFluxOut = 0;
889 model().linFlux( std::get< 0 >( uOut ), std::get< 1 >( uOut ), xOut, std::get< 0 >( phiOut ), std::get< 1 >( phiOut ), dFluxOut );
890
891 model().init( inside );
892 JacobianRangeType dFluxPrimeIn( 0 );
893 model().linFlux( std::get< 0 >( uIn ), std::get< 1 >( uJump ), xIn, 0, std::get< 1 >( phiJump ), dFluxPrimeIn );
894
895 RangeType int0 = std::get< 0 >( phiJump );
896 int0 *= beta_;
897 dFluxOut.usmv( -half, normal, int0 );
898
899 dFluxPrimeIn *= -half;
900 dFluxPrimeOut *= -half;
901
902 return std::make_pair( RangeValueType( -int0, dFluxPrimeIn ), RangeValueType( int0, dFluxPrimeOut ) );
903 };
904
905 return std::make_pair( intIn, intOut );
906 }
907
908 const Model &model () const { return *model_; }
909
910 private:
911 const IntersectionType &intersection () const { assert( intersection_ ); return *intersection_; }
912
913 const Model *model_ = nullptr;
914 RangeFieldType penalty_;
915 const IntersectionType *intersection_ = nullptr;
916 RangeFieldType beta_;
917 };
918
919 } // namespace Fem
920
921} // namespace Dune
922
923#endif // #ifndef DUNE_FEM_SCHEMES_INTEGRANDS_HH
STL namespace.
double min(const Dune::Fem::Double &v, const double p)
Definition: double.hh:953
Definition: bindguard.hh:11
std::tuple_element< i, Tuple >::type & get(Dune::TypeIndexedTuple< Tuple, Types > &tuple)
Definition: typeindexedtuple.hh:122
wrapper for a (Quadrature,int) pair
Definition: quadrature.hh:43
Definition: integrands.hh:110
static const bool boundary
Definition: integrands.hh:121
static const bool hasBoundary
Definition: integrands.hh:122
Impl::IntegrandsTraits::DomainValueType< Integrands > DomainValueType
Definition: integrands.hh:111
static const bool hasSkeleton
Definition: integrands.hh:125
Impl::IntegrandsTraits::EntityType< Integrands > EntityType
Definition: integrands.hh:115
Impl::IntegrandsTraits::RRangeType< Integrands > RRangeType
Definition: integrands.hh:133
static const bool interior
Definition: integrands.hh:118
Impl::IntegrandsTraits::GridPartType< Integrands > GridPartType
Definition: integrands.hh:113
Impl::IntegrandsTraits::IntersectionType< Integrands > IntersectionType
Definition: integrands.hh:116
static const bool hasInterior
Definition: integrands.hh:119
static const bool isFull
Definition: integrands.hh:131
Impl::IntegrandsTraits::RangeValueType< Integrands > RangeValueType
Definition: integrands.hh:112
static const bool skeleton
Definition: integrands.hh:124
Impl::IntegrandsTraits::DirichletComponentType< Integrands > DirichletComponentType
Definition: integrands.hh:134
Definition: integrands.hh:144
Integrands::type RealIntegrands
Definition: integrands.hh:317
FullIntegrands(Args &&... args)
Definition: integrands.hh:264
IntegrandsTraits< Integrands >::GridPartType GridPartType
Definition: integrands.hh:147
RealIntegrands rInt_
Definition: integrands.hh:326
void unbind()
Definition: integrands.hh:272
bool init(const IntersectionType &intersection)
Definition: integrands.hh:271
const RealIntegrands & integrands() const
Definition: integrands.hh:321
bool hasInterior() const
Definition: integrands.hh:274
bool hasSkeleton() const
Definition: integrands.hh:302
auto linearizedBoundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:297
std::pair< RangeValueType, RangeValueType > skeleton(const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:305
bool hasBoundary() const
Definition: integrands.hh:288
RealIntegrands & integrands()
Definition: integrands.hh:320
bool isDirichletIntersection(const IntersectionType &inter, DirichletComponentType &dirichletComponent) const
Definition: integrands.hh:335
IntegrandsTraits< Integrands >::IntersectionType IntersectionType
Definition: integrands.hh:150
IntegrandsTraits< Integrands >::EntityType EntityType
Definition: integrands.hh:149
bool hasDirichletBoundary() const
Definition: integrands.hh:331
bool init(const EntityType &entity)
Definition: integrands.hh:270
IntegrandsTraits< Integrands >::DirichletComponentType DirichletComponentType
Definition: integrands.hh:330
void dirichlet(int bndId, const Point &x, RRangeType &value) const
Definition: integrands.hh:340
IntegrandsTraits< Integrands >::RRangeType RRangeType
Definition: integrands.hh:329
RangeValueType interior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:277
auto linearizedSkeleton(const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:311
auto linearizedInterior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:283
RangeValueType boundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:291
IntegrandsTraits< Integrands >::DomainValueType DomainValueType
Definition: integrands.hh:145
IntegrandsTraits< Integrands >::RangeValueType RangeValueType
Definition: integrands.hh:146
Integrands integrands_
Definition: integrands.hh:325
Definition: integrands.hh:377
RangeValueType boundary(const Fem::QuadraturePointWrapper< Quadrature > &x, const DomainValueType &u) const
Definition: integrands.hh:547
DomainValue DomainValueType
Definition: integrands.hh:382
RangeValueType interior(const Fem::QuadraturePointWrapper< Quadrature > &x, const DomainValueType &u) const
Definition: integrands.hh:521
VirtualizedIntegrands & operator=(const This &other)
Definition: integrands.hh:505
auto linearizedBoundary(const Fem::QuadraturePointWrapper< Quadrature > &x, const DomainValueType &u) const
Definition: integrands.hh:559
GridPart GridPartType
Definition: integrands.hh:381
VirtualizedIntegrands(Integrands integrands)
Definition: integrands.hh:491
VirtualizedIntegrands(This &&)=default
std::array< int, RRangeType::dimension > DirichletComponentType
Definition: integrands.hh:389
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: integrands.hh:385
void unbind()
Definition: integrands.hh:516
typename detail::GetDimRange< std::tuple_element_t< 0, RangeValueType > >::type RRangeType
Definition: integrands.hh:388
VirtualizedIntegrands(const This &other)
Definition: integrands.hh:500
bool hasSkeleton() const
Definition: integrands.hh:570
auto linearizedSkeleton(const Fem::QuadraturePointWrapper< Quadrature > &xIn, const DomainValueType &uIn, const Fem::QuadraturePointWrapper< Quadrature > &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:585
EntityType::Geometry::LocalCoordinate DomainType
Definition: integrands.hh:390
bool init(const IntersectionType &intersection)
Definition: integrands.hh:515
bool isDirichletIntersection(const IntersectionType &inter, DirichletComponentType &dirichletComponent) const
Definition: integrands.hh:600
RangeValue RangeValueType
Definition: integrands.hh:383
bool hasInterior() const
Definition: integrands.hh:518
GridPartType::IntersectionType IntersectionType
Definition: integrands.hh:386
bool hasBoundary() const
Definition: integrands.hh:544
auto linearizedInterior(const Fem::QuadraturePointWrapper< Quadrature > &x, const DomainValueType &u) const
Definition: integrands.hh:533
std::pair< RangeValueType, RangeValueType > skeleton(const Fem::QuadraturePointWrapper< Quadrature > &xIn, const DomainValueType &uIn, const Fem::QuadraturePointWrapper< Quadrature > &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:573
bool hasDirichletBoundary() const
Definition: integrands.hh:596
bool init(const EntityType &entity)
Definition: integrands.hh:514
void dirichlet(int bndId, const DomainType &x, RRangeType &value) const
Definition: integrands.hh:604
Definition: integrands.hh:623
Model::RangeType RangeType
Definition: integrands.hh:630
std::tuple< RangeType, JacobianRangeType > DomainValueType
Definition: integrands.hh:633
auto linearizedInterior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:658
bool init(const EntityType &entity)
Definition: integrands.hh:638
void unbind()
Definition: integrands.hh:645
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: integrands.hh:627
RangeValueType boundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:670
bool init(const IntersectionType &intersection)
Definition: integrands.hh:640
Model ModelType
Definition: integrands.hh:624
GridPartType::IntersectionType IntersectionType
Definition: integrands.hh:628
Model::GridPartType GridPartType
Definition: integrands.hh:625
RangeValueType interior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:648
const Model & model() const
Definition: integrands.hh:687
ConservationLawModelIntegrands(const Model &model)
Definition: integrands.hh:636
std::tuple< RangeType, JacobianRangeType > RangeValueType
Definition: integrands.hh:634
Model::JacobianRangeType JacobianRangeType
Definition: integrands.hh:631
auto linearizedBoundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:678
Definition: integrands.hh:700
auto linearizedBoundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:779
Model::JacobianRangeType JacobianRangeType
Definition: integrands.hh:708
bool init(const IntersectionType &intersection)
Definition: integrands.hh:725
std::tuple< RangeType, JacobianRangeType > RangeValueType
Definition: integrands.hh:713
const Model & model() const
Definition: integrands.hh:908
GridPartType::IntersectionType IntersectionType
Definition: integrands.hh:705
FieldTraits< RangeType >::field_type RangeFieldType
Definition: integrands.hh:710
Model::RangeType RangeType
Definition: integrands.hh:707
RangeValueType boundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:771
RangeValueType interior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:749
Model::GridPartType GridPartType
Definition: integrands.hh:702
bool init(const EntityType &entity)
Definition: integrands.hh:719
auto linearizedSkeleton(const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:830
std::pair< RangeValueType, RangeValueType > skeleton(const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:789
auto linearizedInterior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:759
std::tuple< RangeType, JacobianRangeType > DomainValueType
Definition: integrands.hh:712
DGConservationLawModelIntegrands(const Model &model, RangeFieldType penalty)
Definition: integrands.hh:715
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: integrands.hh:704
Model ModelType
Definition: integrands.hh:701
void unbind()
Definition: integrands.hh:743