DUNE-FEM (unstable)

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