dune-fem  2.8-git
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 
19 namespace 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  {}
267 
268  bool init ( const EntityType &entity ) { return integrands().init( entity ); }
269  bool init ( const IntersectionType &intersection ) { return integrands().init( intersection ); }
270 
271  bool hasInterior () const { return hasInterior( integrands() ); }
272 
273  template< class Point >
274  RangeValueType interior ( const Point &x, const DomainValueType &u ) const
275  {
276  return interior( integrands(), x, u );
277  }
278 
279  template< class Point >
280  auto linearizedInterior ( const Point &x, const DomainValueType &u ) const
281  {
282  return linearizedInterior( integrands(), x, u );
283  }
284 
285  bool hasBoundary () const { return hasBoundary( integrands() ); }
286 
287  template< class Point >
288  RangeValueType boundary ( const Point &x, const DomainValueType &u ) const
289  {
290  return boundary( integrands(), x, u );
291  }
292 
293  template< class Point >
294  auto linearizedBoundary ( const Point &x, const DomainValueType &u ) const
295  {
296  return linearizedBoundary( integrands(), x, u );
297  }
298 
299  bool hasSkeleton () const { return hasSkeleton( integrands() ); }
300 
301  template< class Point >
302  std::pair< RangeValueType, RangeValueType > skeleton ( const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut ) const
303  {
304  return skeleton( integrands(), xIn, uIn, xOut, uOut );
305  }
306 
307  template< class Point >
308  auto linearizedSkeleton ( const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut ) const
309  {
310  return linearizedSkeleton( integrands(), xIn, uIn, xOut, uOut );
311  }
312 
313  protected:
314  decltype( auto ) integrands () { return std::ref( integrands_ ).get(); }
315  decltype( auto ) integrands () const { return std::ref( integrands_ ).get(); }
316 
317  Integrands integrands_;
318 
319  public:
322  bool hasDirichletBoundary () const
323  {
324  return integrands_.get().hasDirichletBoundary();
325  }
326  bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const
327  {
328  return integrands_.get().isDirichletIntersection(inter,dirichletComponent);
329  }
330  template <class Point>
331  void dirichlet( int bndId, const Point &x, RRangeType &value) const
332  {
333  return integrands_.get().dirichlet(bndId,x,value);
334  }
335  };
336 
337 
338 
339  // VirtualizedIntegrands
340  // ---------------------
341 
342  namespace detail
343  {
344  template <class T>
345  struct GetDimRange;
346  template <class FT,int r>
347  struct GetDimRange<Dune::FieldVector<FT,r>>
348  {
349  typedef Dune::FieldVector<FT,r> type;
350  static const int value = r;
351  };
352  template <class FT,int r,int c>
353  struct GetDimRange<Dune::FieldMatrix<FT,r,c>>
354  {
355  typedef Dune::FieldVector<FT,r> type;
356  static const int value = r;
357  };
358  template <class FT,int r,int c>
359  struct GetDimRange<Dune::Fem::ExplicitFieldVector<Dune::FieldMatrix<FT,c,c>,r>>
360  {
361  typedef Dune::FieldVector<FT,r> type;
362  static const int value = r;
363  };
364  }
365 
366  template< class GridPart, class DomainValue, class RangeValue = DomainValue >
368  {
370 
371  public:
372  typedef GridPart GridPartType;
373  typedef DomainValue DomainValueType;
374  typedef RangeValue RangeValueType;
375 
376  typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
377  typedef typename GridPartType::IntersectionType IntersectionType;
378 
379  using RRangeType = typename detail::GetDimRange<std::tuple_element_t<0,RangeValueType>>::type;
380  typedef std::array<int,RRangeType::dimension> DirichletComponentType;
381  typedef typename EntityType::Geometry::LocalCoordinate DomainType;
382 
383  private:
384  typedef typename EntityType::Geometry::LocalCoordinate LocalCoordinateType;
385 
386  typedef FemPy::CachingPoint< LocalCoordinateType, 0 > InteriorCachingPointType;
387  typedef FemPy::ElementPoint< LocalCoordinateType, 0 > InteriorElementPointType;
388  typedef FemPy::CachingPoint< LocalCoordinateType, 1 > SurfaceCachingPointType;
389  typedef FemPy::ElementPoint< LocalCoordinateType, 1 > SurfaceElementPointType;
390 
391  template< class QP >
392  static Fem::QuadraturePointWrapper< QP > asQP ( const QP &qp )
393  {
394  return static_cast< Fem::QuadraturePointWrapper< QP > >( qp );
395  }
396 
397  template< class R >
398  using Linearization = std::function< R( const DomainValueType & ) >;
399 
400  template< class R >
401  using LinearizationPair = std::pair< Linearization< std::pair< R, R > >, Linearization< std::pair< R, R > > >;
402 
403  struct Interface
404  {
405  virtual ~Interface () {}
406  virtual Interface *clone () const = 0;
407 
408  virtual bool init ( const EntityType &entity ) = 0;
409  virtual bool init ( const IntersectionType &intersection ) = 0;
410 
411  virtual bool hasInterior () const = 0;
412  virtual RangeValueType interior ( const InteriorCachingPointType &x, const DomainValueType &u ) const = 0;
413  virtual RangeValueType interior ( const InteriorElementPointType &x, const DomainValueType &u ) const = 0;
414  virtual Linearization< RangeValueType > linearizedInterior ( const InteriorCachingPointType &x, const DomainValueType &u ) const = 0;
415  virtual Linearization< RangeValueType > linearizedInterior ( const InteriorElementPointType &x, const DomainValueType &u ) const = 0;
416 
417  virtual bool hasBoundary () const = 0;
418  virtual RangeValueType boundary ( const SurfaceCachingPointType &x, const DomainValueType &u ) const = 0;
419  virtual RangeValueType boundary ( const SurfaceElementPointType &x, const DomainValueType &u ) const = 0;
420  virtual Linearization< RangeValueType > linearizedBoundary ( const SurfaceCachingPointType &x, const DomainValueType &u ) const = 0;
421  virtual Linearization< RangeValueType > linearizedBoundary ( const SurfaceElementPointType &x, const DomainValueType &u ) const = 0;
422 
423  virtual bool hasSkeleton () const = 0;
424  virtual std::pair< RangeValueType, RangeValueType > skeleton ( const SurfaceCachingPointType &xIn, const DomainValueType &uIn, const SurfaceCachingPointType &xOut, const DomainValueType &uOut ) const = 0;
425  virtual std::pair< RangeValueType, RangeValueType > skeleton ( const SurfaceElementPointType &xIn, const DomainValueType &uIn, const SurfaceElementPointType &xOut, const DomainValueType &uOut ) const = 0;
426  virtual LinearizationPair< RangeValueType > linearizedSkeleton ( const SurfaceCachingPointType &xIn, const DomainValueType &uIn, const SurfaceCachingPointType &xOut, const DomainValueType &uOut ) const = 0;
427  virtual LinearizationPair< RangeValueType > linearizedSkeleton ( const SurfaceElementPointType &xIn, const DomainValueType &uIn, const SurfaceElementPointType &xOut, const DomainValueType &uOut ) const = 0;
428 
429  virtual bool hasDirichletBoundary () const = 0;
430  virtual bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const = 0;
431  virtual void dirichlet( int bndId, const DomainType &x,RRangeType &value) const = 0;
432  };
433 
434  template< class Impl >
435  struct Implementation final
436  : public Interface
437  {
438  Implementation ( Impl impl ) : impl_( std::move( impl ) ) {}
439  virtual Interface *clone () const override { return new Implementation( *this ); }
440 
441  virtual bool init ( const EntityType &entity ) override { return impl().init( entity ); }
442  virtual bool init ( const IntersectionType &intersection ) override { return impl().init( intersection ); }
443 
444  virtual bool hasInterior () const override { return impl().hasInterior(); }
445  virtual RangeValueType interior ( const InteriorCachingPointType &x, const DomainValueType &u ) const override { return impl().interior( asQP( x ), u ); }
446  virtual RangeValueType interior ( const InteriorElementPointType &x, const DomainValueType &u ) const override { return impl().interior( asQP( x ), u ); }
447  virtual Linearization< RangeValueType > linearizedInterior ( const InteriorCachingPointType &x, const DomainValueType &u ) const override { return impl().linearizedInterior( asQP( x ), u ); }
448  virtual Linearization< RangeValueType > linearizedInterior ( const InteriorElementPointType &x, const DomainValueType &u ) const override { return impl().linearizedInterior( asQP( x ), u ); }
449 
450  virtual bool hasBoundary () const override { return impl().hasBoundary(); }
451  virtual RangeValueType boundary ( const SurfaceCachingPointType &x, const DomainValueType &u ) const override { return impl().boundary( asQP( x ), u ); }
452  virtual RangeValueType boundary ( const SurfaceElementPointType &x, const DomainValueType &u ) const override { return impl().boundary( asQP( x ), u ); }
453  virtual Linearization< RangeValueType > linearizedBoundary ( const SurfaceCachingPointType &x, const DomainValueType &u ) const override { return impl().linearizedBoundary( asQP( x ), u ); }
454  virtual Linearization< RangeValueType > linearizedBoundary ( const SurfaceElementPointType &x, const DomainValueType &u ) const override { return impl().linearizedBoundary( asQP( x ), u ); }
455 
456  virtual bool hasSkeleton () const override { return impl().hasSkeleton(); }
457  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 ); }
458  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 ); }
459  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 ); }
460  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 ); }
461 
462  virtual bool hasDirichletBoundary () const override { return impl().hasDirichletBoundary(); }
463  virtual bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const override { return impl().isDirichletIntersection(inter,dirichletComponent); }
464  virtual void dirichlet( int bndId, const DomainType &x,RRangeType &value) const override { impl().dirichlet(bndId,x,value); }
465 
466  private:
467  const auto &impl () const { return std::cref( impl_ ).get(); }
468  auto &impl () { return std::ref( impl_ ).get(); }
469 
470  Impl impl_;
471  };
472 
473  template< class Integrands >
474  using isVirtualized = std::is_same< std::decay_t< decltype( std::ref( std::declval< Integrands & >() ).get() ) >, This >;
475 
476  public:
477  template< class Integrands, std::enable_if_t< IntegrandsTraits< std::decay_t< Integrands > >::isFull && !isVirtualized< Integrands >::value, int > = 0 >
478  explicit VirtualizedIntegrands ( Integrands integrands )
479  : impl_( new Implementation< Integrands >( std::move( integrands ) ) )
480  {}
481 
482  template< class Integrands, std::enable_if_t< !IntegrandsTraits< Integrands >::isFull && !isVirtualized< Integrands >::value, int > = 0 >
483  explicit VirtualizedIntegrands ( Integrands integrands )
484  : VirtualizedIntegrands( FullIntegrands< std::decay_t< Integrands > >( std::move( integrands ) ) )
485  {}
486 
487  VirtualizedIntegrands ( const This &other ) : impl_( other ? other.impl().clone() : nullptr )
488  {}
489 
490  VirtualizedIntegrands ( This && ) = default;
491 
493  {
494  impl_.reset( other ? other.impl().clone() : nullptr );
495  }
496 
498 
499  explicit operator bool () const { return static_cast< bool >( impl_ ); }
500 
501  bool init ( const EntityType &entity ) { return impl().init( entity ); }
502  bool init ( const IntersectionType &intersection ) { return impl().init( intersection ); }
503 
504  bool hasInterior () const { return impl().hasInterior(); }
505 
506  template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
508  {
509  return impl().interior( InteriorCachingPointType( x ), u );
510  }
511 
512  template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
514  {
515  return impl().interior( InteriorElementPointType( x ), u );
516  }
517 
518  template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
520  {
521  return impl().linearizedInterior( InteriorCachingPointType( x ), u );
522  }
523 
524  template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
526  {
527  return impl().linearizedInterior( InteriorElementPointType( x ), u );
528  }
529 
530  bool hasBoundary () const { return impl().hasBoundary(); }
531 
532  template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
534  {
535  return impl().boundary( SurfaceCachingPointType( 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().boundary( SurfaceElementPointType( x ), u );
542  }
543 
544  template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
546  {
547  return impl().linearizedBoundary( SurfaceCachingPointType( x ), u );
548  }
549 
550  template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
552  {
553  return impl().linearizedBoundary( SurfaceElementPointType( x ), u );
554  }
555 
556  bool hasSkeleton () const { return impl().hasSkeleton(); }
557 
558  template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
559  std::pair< RangeValueType, RangeValueType > skeleton ( const Fem::QuadraturePointWrapper< Quadrature > &xIn, const DomainValueType &uIn, const Fem::QuadraturePointWrapper< Quadrature > &xOut, const DomainValueType &uOut ) const
560  {
561  return impl().skeleton( SurfaceCachingPointType( xIn ), uIn, SurfaceCachingPointType( xOut ), uOut );
562  }
563 
564  template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
565  std::pair< RangeValueType, RangeValueType > skeleton ( const Fem::QuadraturePointWrapper< Quadrature > &xIn, const DomainValueType &uIn, const Fem::QuadraturePointWrapper< Quadrature > &xOut, const DomainValueType &uOut ) const
566  {
567  return impl().skeleton( SurfaceElementPointType( xIn ), uIn, SurfaceElementPointType( xOut ), uOut );
568  }
569 
570  template< class Quadrature, std::enable_if_t< std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
572  {
573  return impl().linearizedSkeleton( SurfaceCachingPointType( xIn ), uIn, SurfaceCachingPointType( xOut ), uOut );
574  }
575 
576  template< class Quadrature, std::enable_if_t< !std::is_convertible< Quadrature, Fem::CachingInterface >::value, int > = 0 >
578  {
579  return impl().linearizedSkeleton( SurfaceElementPointType( xIn ), uIn, SurfaceElementPointType( xOut ), uOut );
580  }
581 
582  bool hasDirichletBoundary () const
583  {
584  return impl().hasDirichletBoundary();
585  }
586  bool isDirichletIntersection( const IntersectionType& inter, DirichletComponentType &dirichletComponent ) const
587  {
588  return impl().isDirichletIntersection(inter,dirichletComponent);
589  }
590  void dirichlet( int bndId, const DomainType &x,RRangeType &value) const
591  {
592  return impl().dirichlet(bndId,x,value);
593  }
594 
595  private:
596  const Interface &impl () const { assert( impl_ ); return *impl_; }
597  Interface &impl () { assert( impl_ ); return *impl_; }
598 
599  std::unique_ptr< Interface > impl_;
600  };
601 
602 
603 
604  // DiffusionModelIntegrands
605  // ------------------------
606 
607  template< class Model >
609  {
610  typedef Model ModelType;
611  typedef typename Model::GridPartType GridPartType;
612 
613  typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
614  typedef typename GridPartType::IntersectionType IntersectionType;
615 
616  typedef typename Model::RangeType RangeType;
617  typedef typename Model::JacobianRangeType JacobianRangeType;
618 
619  typedef std::tuple< RangeType, JacobianRangeType > DomainValueType;
620  typedef std::tuple< RangeType, JacobianRangeType > RangeValueType;
621 
622  explicit DiffusionModelIntegrands ( const Model &model ) : model_( &model ) {}
623 
624  bool init ( const EntityType &entity ) { return model().init( entity ); }
625 
626  bool init ( const IntersectionType &intersection )
627  {
628  return (intersection.boundary() && model().hasNeumanBoundary() && model().init( intersection.inside() ));
629  }
630 
631  template< class Point >
632  RangeValueType interior ( const Point &x, const DomainValueType &u ) const
633  {
634  RangeType source( 0 );
635  model().source( x, std::get< 0 >( u ), std::get< 1 >( u ), source );
636  JacobianRangeType dFlux( 0 );
637  model().diffusiveFlux( x, std::get< 0 >( u ), std::get< 1 >( u ), dFlux );
638  return std::make_tuple( source, dFlux );
639  }
640 
641  template< class Point >
642  auto linearizedInterior ( const Point &x, const DomainValueType &u ) const
643  {
644  return [ this, x, u ] ( const DomainValueType &phi ) {
645  RangeType source( 0 );
646  model().linSource( std::get< 0 >( u ), std::get< 1 >( u ), x, std::get< 0 >( phi ), std::get< 1 >( phi ), source );
647  JacobianRangeType dFlux( 0 );
648  model().linDiffusiveFlux( std::get< 0 >( u ), std::get< 1 >( u ), x, std::get< 0 >( phi ), std::get< 1 >( phi ), dFlux );
649  return std::make_tuple( source, dFlux );
650  };
651  }
652 
653  template< class Point >
654  RangeValueType boundary ( const Point &x, const DomainValueType &u ) const
655  {
656  RangeType alpha( 0 );
657  model().alpha( x, std::get< 0 >( u ), alpha );
658  return std::make_tuple( alpha, 0 );
659  }
660 
661  template< class Point >
662  auto linearizedBoundary ( const Point &x, const DomainValueType &u ) const
663  {
664  return [ this, x, u ] ( const DomainValueType &phi ) {
665  RangeType alpha( 0 );
666  model().linAlpha( std::get< 0 >( u ), x, std::get< 0 >( phi ), alpha );
667  return std::make_tuple( alpha, 0 );
668  };
669  }
670 
671  const Model &model () const { return *model_; }
672 
673  private:
674  const Model *model_;
675  };
676 
677 
678 
679  // DGDiffusionModelIntegrands
680  // --------------------------
681 
682  template< class Model >
684  {
685  typedef Model ModelType;
686  typedef typename Model::GridPartType GridPartType;
687 
688  typedef typename GridPartType::template Codim< 0 >::EntityType EntityType;
689  typedef typename GridPartType::IntersectionType IntersectionType;
690 
691  typedef typename Model::RangeType RangeType;
692  typedef typename Model::JacobianRangeType JacobianRangeType;
693 
694  typedef typename FieldTraits< RangeType >::field_type RangeFieldType;
695 
696  typedef std::tuple< RangeType, JacobianRangeType > DomainValueType;
697  typedef std::tuple< RangeType, JacobianRangeType > RangeValueType;
698 
700  : model_( &model ), penalty_( penalty )
701  {}
702 
703  bool init ( const EntityType &entity )
704  {
705  intersection_ = nullptr;
706  return model().init( entity );
707  }
708 
709  bool init ( const IntersectionType &intersection )
710  {
711  intersection_ = &intersection;
712  if( intersection.boundary() )
713  {
714  const EntityType inside = intersection.inside();
715  beta_ = penalty_ * intersection.geometry().volume() / inside.geometry().volume();
716  return (model().hasNeumanBoundary() && model().init( inside ));
717  }
718  else if( intersection.neighbor() )
719  {
720  const auto volIn = intersection.inside().geometry().volume();
721  const auto volOut = intersection.outside().geometry().volume();
722  beta_ = penalty_ * intersection.geometry().volume() / std::min( volIn, volOut );
723  return true;
724  }
725  }
726 
727  template< class Point >
728  RangeValueType interior ( const Point &x, const DomainValueType &u ) const
729  {
730  RangeType source( 0 );
731  model().source( x, std::get< 0 >( u ), std::get< 1 >( u ), source );
732  JacobianRangeType dFlux( 0 );
733  model().diffusiveFlux( x, std::get< 0 >( u ), std::get< 1 >( u ), dFlux );
734  return RangeValueType( source, dFlux );
735  }
736 
737  template< class Point >
738  auto linearizedInterior ( const Point &x, const DomainValueType &u ) const
739  {
740  return [ this, x, u ] ( const DomainValueType &phi ) {
741  RangeType source( 0 );
742  model().linSource( std::get< 0 >( u ), std::get< 1 >( u ), x, std::get< 0 >( phi ), std::get< 1 >( phi ), source );
743  JacobianRangeType dFlux( 0 );
744  model().linDiffusiveFlux( std::get< 0 >( u ), std::get< 1 >( u ), x, std::get< 0 >( phi ), std::get< 1 >( phi ), dFlux );
745  return RangeValueType( source, dFlux );
746  };
747  }
748 
749  template< class Point >
750  RangeValueType boundary ( const Point &x, const DomainValueType &u ) const
751  {
752  RangeType alpha( 0 );
753  model().alpha( x, std::get< 0 >( u ), alpha );
754  return RangeValueType( alpha, 0 );
755  }
756 
757  template< class Point >
758  auto linearizedBoundary ( const Point &x, const DomainValueType &u ) const
759  {
760  return [ this, x, u ] ( const DomainValueType &phi ) {
761  RangeType alpha( 0 );
762  model().linAlpha( std::get< 0 >( u ), x, std::get< 0 >( phi ), alpha );
763  return RangeValueType( alpha, 0 );
764  };
765  }
766 
767  template< class Point >
768  std::pair< RangeValueType, RangeValueType > skeleton ( const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut ) const
769  {
770  const EntityType inside = intersection().inside();
771  const EntityType outside = intersection().outside();
772 
773  const RangeFieldType half = RangeFieldType( 1 ) / RangeFieldType( 2 );
774  const auto normal = intersection().unitOuterNormal( xIn.localPosition() );
775 
776  DomainValueType uJump( 0, 0 );
777  std::get< 0 >( uJump ) = std::get< 0 >( uOut ) - std::get< 0 >( uIn );
778  for( int i = 0; i < RangeType::dimension; ++i )
779  std::get< 1 >( uJump )[ i ].axpy( std::get< 0 >( uJump )[ i ], normal );
780 
781  model().init( outside );
782  JacobianRangeType dFluxOut( 0 ), dFluxPrimeOut( 0 );
783  model().diffusiveFlux( xOut, std::get< 0 >( uOut ), std::get< 1 >( uJump ), dFluxPrimeOut );
784  model().diffusiveFlux( xOut, std::get< 0 >( uOut ), 0, dFluxOut );
785  dFluxPrimeOut -= dFluxOut;
786  dFluxOut = 0;
787  model().diffusiveFlux( xOut, std::get< 0 >( uOut ), std::get< 1 >( uOut ), dFluxOut );
788 
789  model().init( inside );
790  JacobianRangeType dFluxIn( 0 ), dFluxPrimeIn( 0 );
791  model().diffusiveFlux( xIn, std::get< 0 >( uIn ), std::get< 1 >( uJump ), dFluxPrimeIn );
792  model().diffusiveFlux( xIn, std::get< 0 >( uIn ), 0, dFluxIn );
793  dFluxPrimeIn -= dFluxIn;
794  dFluxIn = 0;
795  model().diffusiveFlux( xIn, std::get< 0 >( uIn ), std::get< 1 >( uIn ), dFluxIn );
796 
797  RangeType int0 = std::get< 0 >( uJump );
798  int0 *= beta_;
799  dFluxIn += dFluxOut;
800  dFluxIn.usmv( -half, normal, int0 );
801 
802  dFluxPrimeIn *= -half;
803  dFluxPrimeOut *= -half;
804 
805  return std::make_pair( RangeValueType( -int0, dFluxPrimeIn ), RangeValueType( int0, dFluxPrimeOut ) );
806  }
807 
808  template< class Point >
809  auto linearizedSkeleton ( const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut ) const
810  {
811  const auto normal = intersection().unitOuterNormal( xIn.localPosition() );
812 
813  DomainValueType uJump( 0, 0 );
814  std::get< 0 >( uJump ) = std::get< 0 >( uOut ) - std::get< 0 >( uIn );
815  for( int i = 0; i < RangeType::dimension; ++i )
816  std::get< 1 >( uJump )[ i ].axpy( std::get< 0 >( uJump )[ i ], normal );
817 
818  auto intIn = [ this, xIn, uIn, xOut, uOut, normal, uJump ] ( const DomainValueType &phiIn ) {
819  const EntityType inside = intersection().inside();
820  const EntityType outside = intersection().outside();
821 
822  const RangeFieldType half = RangeFieldType( 1 ) / RangeFieldType( 2 );
823 
824  DomainValueType phiJump( 0, 0 );
825  std::get< 0 >( phiJump ) -= std::get< 0 >( phiIn );
826  for( int i = 0; i < RangeType::dimension; ++i )
827  std::get< 1 >( phiJump )[ i ].axpy( std::get< 0 >( phiJump )[ i ], normal );
828 
829  model().init( outside );
830  JacobianRangeType dFluxPrimeOut( 0 );
831  model().linDiffusiveFlux( std::get< 0 >( uOut ), std::get< 1 >( uJump ), xOut, 0, std::get< 1 >( phiJump ), dFluxPrimeOut );
832 
833  model().init( inside );
834  JacobianRangeType dFluxIn( 0 ), dFluxPrimeIn( 0 );
835  model().linDiffusiveFlux( std::get< 0 >( uIn ), std::get< 1 >( uJump ), xIn, std::get< 0 >( phiIn ), std::get< 1 >( phiJump ), dFluxPrimeIn );
836  model().linDiffusiveFlux( std::get< 0 >( uIn ), 0, xIn, std::get< 0 >( phiIn ), 0, dFluxIn );
837  dFluxPrimeIn -= dFluxIn;
838  dFluxIn = 0;
839  model().linDiffusiveFlux( std::get< 0 >( uIn ), std::get< 1 >( uIn ), xIn, std::get< 0 >( phiIn ), std::get< 1 >( phiIn ), dFluxIn );
840 
841  RangeType int0 = std::get< 0 >( phiJump );
842  int0 *= beta_;
843  dFluxIn.usmv( -half, normal, int0 );
844 
845  dFluxPrimeIn *= -half;
846  dFluxPrimeOut *= -half;
847 
848  return std::make_pair( RangeValueType( -int0, dFluxPrimeIn ), RangeValueType( int0, dFluxPrimeOut ) );
849  };
850 
851  auto intOut = [ this, xIn, uIn, xOut, uOut, normal, uJump ] ( const DomainValueType &phiOut ) {
852  const EntityType inside = intersection().inside();
853  const EntityType outside = intersection().outside();
854 
855  const RangeFieldType half = RangeFieldType( 1 ) / RangeFieldType( 2 );
856 
857  DomainValueType phiJump( 0, 0 );
858  std::get< 0 >( phiJump ) = std::get< 0 >( phiOut );
859  for( int i = 0; i < RangeType::dimension; ++i )
860  std::get< 1 >( phiJump )[ i ].axpy( std::get< 0 >( phiJump )[ i ], normal );
861 
862  model().init( outside );
863  JacobianRangeType dFluxOut( 0 ), dFluxPrimeOut( 0 );
864  model().linDiffusiveFlux( std::get< 0 >( uOut ), std::get< 1 >( uJump ), xOut, std::get< 0 >( phiOut ), std::get< 1 >( phiJump ), dFluxPrimeOut );
865  model().linDiffusiveFlux( std::get< 0 >( uOut ), 0, xOut, std::get< 0 >( phiOut ), 0, dFluxOut );
866  dFluxPrimeOut -= dFluxOut;
867  dFluxOut = 0;
868  model().linDiffusiveFlux( std::get< 0 >( uOut ), std::get< 1 >( uOut ), xOut, std::get< 0 >( phiOut ), std::get< 1 >( phiOut ), dFluxOut );
869 
870  model().init( inside );
871  JacobianRangeType dFluxPrimeIn( 0 );
872  model().linDiffusiveFlux( std::get< 0 >( uIn ), std::get< 1 >( uJump ), xIn, 0, std::get< 1 >( phiJump ), dFluxPrimeIn );
873 
874  RangeType int0 = std::get< 0 >( phiJump );
875  int0 *= beta_;
876  dFluxOut.usmv( -half, normal, int0 );
877 
878  dFluxPrimeIn *= -half;
879  dFluxPrimeOut *= -half;
880 
881  return std::make_pair( RangeValueType( -int0, dFluxPrimeIn ), RangeValueType( int0, dFluxPrimeOut ) );
882  };
883 
884  return std::make_pair( intIn, intOut );
885  }
886 
887  const Model &model () const { return *model_; }
888 
889  private:
890  const IntersectionType &intersection () const { assert( intersection_ ); return *intersection_; }
891 
892  const Model *model_;
893  RangeFieldType penalty_;
894  const IntersectionType *intersection_ = nullptr;
895  RangeFieldType beta_;
896  };
897 
898  } // namespace Fem
899 
900 } // namespace Dune
901 
902 #endif // #ifndef DUNE_FEM_SCHEMES_INTEGRANDS_HH
Definition: bindguard.hh:11
std::tuple_element< i, Tuple >::type & get(Dune::TypeIndexedTuple< Tuple, Types > &tuple)
Definition: typeindexedtuple.hh:122
static constexpr T min(T a)
Definition: utility.hh:93
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
FullIntegrands(Args &&... args)
Definition: integrands.hh:264
bool init(const IntersectionType &intersection)
Definition: integrands.hh:269
bool hasInterior() const
Definition: integrands.hh:271
bool hasSkeleton() const
Definition: integrands.hh:299
auto linearizedBoundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:294
std::pair< RangeValueType, RangeValueType > skeleton(const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:302
bool hasBoundary() const
Definition: integrands.hh:285
IntegrandsTraits< Integrands >::DomainValueType DomainValueType
Definition: integrands.hh:145
IntegrandsTraits< Integrands >::RRangeType RRangeType
Definition: integrands.hh:320
bool isDirichletIntersection(const IntersectionType &inter, DirichletComponentType &dirichletComponent) const
Definition: integrands.hh:326
bool hasDirichletBoundary() const
Definition: integrands.hh:322
bool init(const EntityType &entity)
Definition: integrands.hh:268
void dirichlet(int bndId, const Point &x, RRangeType &value) const
Definition: integrands.hh:331
IntegrandsTraits< Integrands >::DirichletComponentType DirichletComponentType
Definition: integrands.hh:321
RangeValueType interior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:274
auto linearizedSkeleton(const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:308
IntegrandsTraits< Integrands >::EntityType EntityType
Definition: integrands.hh:149
auto linearizedInterior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:280
IntegrandsTraits< Integrands >::IntersectionType IntersectionType
Definition: integrands.hh:150
RangeValueType boundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:288
IntegrandsTraits< Integrands >::GridPartType GridPartType
Definition: integrands.hh:147
Integrands integrands_
Definition: integrands.hh:317
IntegrandsTraits< Integrands >::RangeValueType RangeValueType
Definition: integrands.hh:146
Definition: integrands.hh:368
RangeValueType boundary(const Fem::QuadraturePointWrapper< Quadrature > &x, const DomainValueType &u) const
Definition: integrands.hh:533
DomainValue DomainValueType
Definition: integrands.hh:373
VirtualizedIntegrands & operator=(const This &other)
Definition: integrands.hh:492
RangeValueType interior(const Fem::QuadraturePointWrapper< Quadrature > &x, const DomainValueType &u) const
Definition: integrands.hh:507
auto linearizedBoundary(const Fem::QuadraturePointWrapper< Quadrature > &x, const DomainValueType &u) const
Definition: integrands.hh:545
GridPart GridPartType
Definition: integrands.hh:372
VirtualizedIntegrands(Integrands integrands)
Definition: integrands.hh:478
VirtualizedIntegrands(This &&)=default
std::array< int, RRangeType::dimension > DirichletComponentType
Definition: integrands.hh:380
typename detail::GetDimRange< std::tuple_element_t< 0, RangeValueType > >::type RRangeType
Definition: integrands.hh:379
VirtualizedIntegrands(const This &other)
Definition: integrands.hh:487
bool hasSkeleton() const
Definition: integrands.hh:556
auto linearizedSkeleton(const Fem::QuadraturePointWrapper< Quadrature > &xIn, const DomainValueType &uIn, const Fem::QuadraturePointWrapper< Quadrature > &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:571
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: integrands.hh:376
EntityType::Geometry::LocalCoordinate DomainType
Definition: integrands.hh:381
bool init(const IntersectionType &intersection)
Definition: integrands.hh:502
bool isDirichletIntersection(const IntersectionType &inter, DirichletComponentType &dirichletComponent) const
Definition: integrands.hh:586
RangeValue RangeValueType
Definition: integrands.hh:374
bool hasInterior() const
Definition: integrands.hh:504
GridPartType::IntersectionType IntersectionType
Definition: integrands.hh:377
bool hasBoundary() const
Definition: integrands.hh:530
auto linearizedInterior(const Fem::QuadraturePointWrapper< Quadrature > &x, const DomainValueType &u) const
Definition: integrands.hh:519
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:559
bool hasDirichletBoundary() const
Definition: integrands.hh:582
bool init(const EntityType &entity)
Definition: integrands.hh:501
void dirichlet(int bndId, const DomainType &x, RRangeType &value) const
Definition: integrands.hh:590
Definition: integrands.hh:609
GridPartType::IntersectionType IntersectionType
Definition: integrands.hh:614
Model::RangeType RangeType
Definition: integrands.hh:616
Model ModelType
Definition: integrands.hh:610
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: integrands.hh:613
std::tuple< RangeType, JacobianRangeType > RangeValueType
Definition: integrands.hh:620
Model::JacobianRangeType JacobianRangeType
Definition: integrands.hh:617
auto linearizedBoundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:662
const Model & model() const
Definition: integrands.hh:671
bool init(const IntersectionType &intersection)
Definition: integrands.hh:626
std::tuple< RangeType, JacobianRangeType > DomainValueType
Definition: integrands.hh:619
bool init(const EntityType &entity)
Definition: integrands.hh:624
DiffusionModelIntegrands(const Model &model)
Definition: integrands.hh:622
RangeValueType boundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:654
RangeValueType interior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:632
auto linearizedInterior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:642
Model::GridPartType GridPartType
Definition: integrands.hh:611
Definition: integrands.hh:684
std::pair< RangeValueType, RangeValueType > skeleton(const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:768
RangeValueType boundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:750
bool init(const IntersectionType &intersection)
Definition: integrands.hh:709
auto linearizedInterior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:738
Model::GridPartType GridPartType
Definition: integrands.hh:686
FieldTraits< RangeType >::field_type RangeFieldType
Definition: integrands.hh:694
RangeValueType interior(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:728
std::tuple< RangeType, JacobianRangeType > RangeValueType
Definition: integrands.hh:697
Model::JacobianRangeType JacobianRangeType
Definition: integrands.hh:692
Model ModelType
Definition: integrands.hh:685
bool init(const EntityType &entity)
Definition: integrands.hh:703
GridPartType::IntersectionType IntersectionType
Definition: integrands.hh:689
DGDiffusionModelIntegrands(const Model &model, RangeFieldType penalty)
Definition: integrands.hh:699
std::tuple< RangeType, JacobianRangeType > DomainValueType
Definition: integrands.hh:696
auto linearizedBoundary(const Point &x, const DomainValueType &u) const
Definition: integrands.hh:758
const Model & model() const
Definition: integrands.hh:887
GridPartType::template Codim< 0 >::EntityType EntityType
Definition: integrands.hh:688
auto linearizedSkeleton(const Point &xIn, const DomainValueType &uIn, const Point &xOut, const DomainValueType &uOut) const
Definition: integrands.hh:809
Model::RangeType RangeType
Definition: integrands.hh:691