dune-fem  2.4.1-rc
localdg/pass.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_PASS_LOCALDG_HH
2 #define DUNE_FEM_PASS_LOCALDG_HH
3 
12 
13 #include "modelcaller.hh"
14 
15 namespace Dune
16 {
17 
18  namespace Fem
19  {
20 
38  template <class DiscreteModelImp, class PreviousPassImp , int passIdImp = -1 >
41  class LocalDGPass :
42  public LocalPass< DiscreteModelImp , PreviousPassImp , passIdImp >
43  {
45  public:
46 
47  //- Typedefs and enums
50 
52  typedef typename BaseType::PassIds PassIds;
53 
55  typedef DiscreteModelImp DiscreteModelType;
57  typedef PreviousPassImp PreviousPassType;
58 
59  // Types from the base class
60  typedef typename BaseType::Entity EntityType;
62 
63  // Types from the traits
64  typedef typename DiscreteModelType::Traits::DestinationType DestinationType;
65  typedef typename DiscreteModelType::Traits::VolumeQuadratureType VolumeQuadratureType;
66  typedef typename DiscreteModelType::Traits::FaceQuadratureType FaceQuadratureType;
67  typedef typename DiscreteModelType::Traits::DiscreteFunctionSpaceType DiscreteFunctionSpaceType;
69  typedef typename DiscreteFunctionSpaceType::IteratorType IteratorType;
70 
71  // Types extracted from the discrete function space type
72  typedef typename DiscreteFunctionSpaceType::GridType GridType;
73  typedef typename DiscreteFunctionSpaceType::GridPartType GridPartType;
74  typedef typename DiscreteFunctionSpaceType::DomainType DomainType;
75  typedef typename DiscreteFunctionSpaceType::RangeType RangeType;
76  typedef typename DiscreteFunctionSpaceType::JacobianRangeType JacobianRangeType;
77  typedef typename DiscreteFunctionSpaceType::BasisFunctionSetType BasisFunctionSetType;
78 
79  // Types extracted from the underlying grids
80  typedef typename GridPartType::IntersectionIteratorType IntersectionIteratorType;
81  typedef typename IntersectionIteratorType::Intersection IntersectionType;
82  typedef typename GridType::template Codim<0>::Geometry Geometry;
83 
84 
85  // Various other types
86  typedef typename DestinationType::LocalFunctionType LocalFunctionType;
87 
89 
90  // Range of the destination
91  enum { dimRange = DiscreteFunctionSpaceType::dimRange };
92 
93  // type of local id set
94  typedef typename GridPartType::IndexSetType IndexSetType;
96 
99 
100  public:
101  //- Public methods
109  LocalDGPass(DiscreteModelType& problem,
110  PreviousPassType& pass,
111  const DiscreteFunctionSpaceType& spc,
112  const int volumeQuadOrd =-1,
113  const int faceQuadOrd=-1,
114  const bool notThreadParallel = true ) :
115  BaseType(pass, spc),
116  caller_(0),
117  problem_(problem),
118  arg_(0),
119  dest_(0),
120  spc_(spc),
121  gridPart_(spc_.gridPart()),
122  indexSet_(gridPart_.indexSet()),
123  visited_(0),
124  updEn_(spc_),
125  updNb_(spc_),
126  fMatVec_( 20 ),
127  valEnVec_( 20 ),
128  valNbVec_( 20 ),
129  dtMin_(std::numeric_limits<double>::max()),
130  minLimit_(2.0*std::numeric_limits<double>::min()),
131  volumeQuadOrd_( (volumeQuadOrd < 0) ?
132  ( 2 * spc_.order()) : volumeQuadOrd ),
133  faceQuadOrd_( (faceQuadOrd < 0) ?
134  ( 2 * spc_.order()+1) : faceQuadOrd ),
136  notThreadParallel_( notThreadParallel )
137  {
138  fMatVec_.setMemoryFactor( 1.1 );
139  valEnVec_.setMemoryFactor( 1.1 );
140  valNbVec_.setMemoryFactor( 1.1 );
141 
142  assert( volumeQuadOrd_ >= 0 );
143  assert( faceQuadOrd_ >= 0 );
144  }
145 
147  virtual ~LocalDGPass() {
148  }
149 
151  void printTexInfo(std::ostream& out) const {
153  out << "LocalDGPass: "
154  << "\\\\ \n";
155  }
156 
158  double timeStepEstimateImpl() const
159  {
160  // factor for LDG Discretization
161  const double p = 2 * spc_.order() + 1;
162  return dtMin_ / p;
163  }
164 
165  public:
168  virtual void prepare(const ArgumentType& arg, DestinationType& dest) const
169  {
170  arg_ = const_cast<ArgumentType*>(&arg);
171  dest_ = &dest;
172 
173  if( notThreadParallel_ )
174  {
175  // clear destination
176  dest_->clear();
177  }
178 
180  assert( caller_ );
181  caller_->setTime( this->time() );
182 
183  // resize indicator function
184  visited_.resize( indexSet_.size(0) );
185  // set all values to false
186  const int indSize = visited_.size();
187  for(int i=0; i<indSize; ++i) visited_[i] = false;
188 
189  // time initialisation to max value
191  }
192 
194  virtual void finalize(const ArgumentType& arg, DestinationType& dest) const
195  {
196  if( notThreadParallel_ )
197  {
198  // communicate calculated function
199  spc_.communicate( dest );
200  }
201 
202  if( caller_ )
203  delete caller_;
204  caller_ = 0;
205  }
206 
207  size_t numberOfElements() const { return 0; }
208 
210  {
211  bool operator ()(const EntityType& , const EntityType& ) const
212  {
213  return true ;
214  }
215  };
216 
217  void applyLocal( const EntityType& en) const
218  {
219  applyLocal( en, DefaultNBChecker() );
220  }
221 
222  void applyLocalMass( const EntityType& en) const
223  {
224  // get local function and add update
225  LocalFunctionType function = dest_->localFunction(en);
226 
227  // apply local inverse mass matrix
228  localMassMatrix_.applyInverse( caller(), en, function );
229  }
230 
232  template <class NeighborChecker>
233  void applyLocal( const EntityType& en, const NeighborChecker& nbChecker ) const
234  {
235  // init local function
236  initLocalFunction( en , updEn_ );
237 
238  // call real apply local
239  applyLocal(en, updEn_, nbChecker );
240 
241  // add update to real function
243  }
244  protected:
246  template <class NeighborChecker>
247  void applyLocal( const EntityType& en,
248  TemporaryLocalFunctionType& updEn,
249  const NeighborChecker& nbChecker ) const
250  {
251  // only call geometry once, who know what is done in this function
252  const Geometry & geo = en.geometry();
253 
254  // get volume of element
255  const double vol = geo.volume();
256 
257  // only apply volumetric integral if order > 0
258  // otherwise this contribution is zero
259 
260  if( (spc_.order() > 0) || problem_.hasSource() )
261  {
262  // create quadrature
263  VolumeQuadratureType volQuad( en, volumeQuadOrd_ );
264 
265  // set entity and evaluate local functions for quadrature
266  caller().setEntity( en , volQuad );
267 
268  if( problem_.hasSource() )
269  {
270  // evaluate flux and source
271  evalVolumetricPartBoth( en, geo, volQuad, updEn );
272  }
273  else if ( problem_.hasFlux() )
274  {
275  // if only flux, evaluate only flux
276  evalVolumetricPartFlux( en, geo, volQuad, updEn );
277  }
278  }
279  else
280  {
281  // only set entity here without evaluation
282  caller().setEntity( en );
283  }
284 
286  // Surface integral part
288  if ( problem_.hasFlux() )
289  {
290  IntersectionIteratorType endnit = gridPart_.iend(en);
291  for (IntersectionIteratorType nit = gridPart_.ibegin(en); nit != endnit; ++nit)
292  {
293  // get intersection from intersection iterator
294  const IntersectionType& intersection = *nit;
295 
296  //double nbvol;
297  double nbvol = vol;
298  double wspeedS = 0.0;
299  if (intersection.neighbor())
300  {
301  // get neighbor
302  EntityType nb = make_entity( intersection.outside() );
303 
304  const bool canUpdateNeighbor = nbChecker( en, nb );
305 
306  if ( ! visited_[ indexSet_.index( nb ) ] )
307  {
308  // for conforming situations apply Quadrature given
310  && !intersection.conforming() )
311  {
312  // apply neighbor part, return is volume of neighbor which is
313  // needed below
314  nbvol = applyLocalNeighbor< false >
315  (intersection,en,nb,
316  updEn, updNb_ ,
317  wspeedS,
318  canUpdateNeighbor );
319  }
320  else
321  {
322  // apply neighbor part, return is volume of neighbor which is
323  // needed below
324  nbvol = applyLocalNeighbor< true >
325  (intersection,en,nb,
326  updEn, updNb_ ,
327  wspeedS,
328  canUpdateNeighbor );
329  }
330 
331  } // end if do something
332 
333  } // end if neighbor
334  else if( intersection.boundary() )
335  {
336  FaceQuadratureType faceQuadInner(gridPart_, intersection, faceQuadOrd_,
337  FaceQuadratureType::INSIDE);
338 
339  // set neighbor entity to inside entity
340  caller().setBoundary(en, faceQuadInner);
341 
342  // cache number of quadrature points
343  const size_t faceQuadInner_nop = faceQuadInner.nop();
344 
345  if( valEnVec_.size() < faceQuadInner_nop )
346  valEnVec_.resize( faceQuadInner_nop );
347 
348  // loop over quadrature points
349  for (size_t l = 0; l < faceQuadInner_nop; ++l)
350  {
351  RangeType& flux = valEnVec_[ l ];
352 
353  // eval boundary Flux
354  wspeedS += caller().boundaryFlux( intersection, faceQuadInner, l, flux )
355  * faceQuadInner.weight(l);
356 
357  // apply weights
358  flux *= -faceQuadInner.weight(l);
359  }
360 
361  // add factor
362  updEn.axpyQuadrature ( faceQuadInner, valEnVec_ );
363 
364  } // end if boundary
365 
366  if (wspeedS > minLimit_ )
367  {
368  double minvolS = std::min(vol,nbvol);
369  dtMin_ = std::min(dtMin_,minvolS/wspeedS);
370  }
371  } // end intersection loop
372 
373  } // end if problem_.hasFlux()
374 
375  // this entity is finised by now
376  visited_[ indexSet_.index( en ) ] = true ;
377  }
378 
379  // initialize local update function
380  template <class LocalFunctionImp>
381  void initLocalFunction( const EntityType& en, LocalFunctionImp& update) const
382  {
383  // init local function
384  update.init( en );
385  // clear dof values
386  update.clear();
387  }
388 
390  template <class LocalFunctionImp>
391  void updateFunction( const EntityType& en,
392  LocalFunctionImp& update) const
393  {
394  // get local function and add update
395  LocalFunctionType function = dest_->localFunction(en);
396  function += update;
397  }
398 
400  template <class LocalFunctionImp>
402  const EntityType& en,
403  LocalFunctionImp& update) const
404  {
405  // get local function and add update
406  LocalFunctionType function = dest_->localFunction(en);
407  function += update;
408 
409  // apply local inverse mass matrix
410  localMassMatrix_.applyInverse( caller(), en, function );
411  }
412 
414  // Volumetric integral part only flux
416  template <class LocalFunctionImp>
417  void evalVolumetricPartFlux( const EntityType& en,
418  const Geometry& geo ,
419  const VolumeQuadratureType& volQuad,
420  LocalFunctionImp& updEn ) const
421  {
422  const size_t volQuad_nop = volQuad.nop();
423  if( fMatVec_.size() < volQuad_nop )
424  {
425  fMatVec_.resize( volQuad_nop );
426  }
427 
428  for (size_t l = 0; l < volQuad_nop; ++l)
429  {
430  JacobianRangeType& flux = fMatVec_[ l ];
431 
432  // evaluate analytical flux and source
433  caller().analyticalFlux(en, volQuad, l, flux );
434 
435  const double intel = geo.integrationElement(volQuad.point(l))
436  * volQuad.weight(l);
437 
438  // apply integration weights
439  flux *= intel;
440  }
441 
442  // add values to local function
443  updEn.axpyQuadrature (volQuad, fMatVec_ );
444  }
445 
447  // Volumetric integral part only flux
449  template <class LocalFunctionImp>
450  void evalVolumetricPartBoth( const EntityType& en,
451  const Geometry& geo,
452  const VolumeQuadratureType& volQuad,
453  LocalFunctionImp& updEn ) const
454  {
455  const size_t volQuad_nop = volQuad.nop();
456 
457  if( fMatVec_.size() < volQuad_nop )
458  {
459  fMatVec_.resize( volQuad_nop );
460  }
461 
462  if( valEnVec_.size() < volQuad_nop )
463  {
464  valEnVec_.resize( volQuad_nop );
465  }
466 
467  for (size_t l = 0; l < volQuad_nop; ++l)
468  {
469  JacobianRangeType& flux = fMatVec_[ l ];
470  RangeType& source = valEnVec_[ l ];
471 
472  // evaluate analytical flux and source
473  const double dtEst =
474  caller().analyticalFluxAndSource(en, volQuad, l, flux, source );
475 
476  const double intel = geo.integrationElement(volQuad.point(l))
477  * volQuad.weight(l);
478 
479  // apply integration weights
480  source *= intel;
481  flux *= intel;
482 
483  if( dtEst > minLimit_ )
484  dtMin_ = std::min(dtMin_, dtEst);
485  }
486 
487  // add values to local function
488  updEn.axpyQuadrature (volQuad, valEnVec_, fMatVec_ );
489  }
490 
491  template <bool conforming, class LocalFunctionImp>
492  double applyLocalNeighbor ( const IntersectionType &intersection,
493  const EntityType &en,
494  const EntityType &nb,
495  LocalFunctionImp &updEn,
496  LocalFunctionImp &updNb,
497  double &wspeedS,
498  const bool canUpdateNeighbor ) const
499  {
500  // make sure we got the right conforming statement
501  assert( intersection.conforming() == conforming );
502 
503  // use IntersectionQuadrature to create appropriate face quadratures
504  typedef IntersectionQuadrature< FaceQuadratureType, conforming > IntersectionQuadratureType;
505  typedef typename IntersectionQuadratureType :: FaceQuadratureType QuadratureImp;
506 
507  // create intersection quadrature
508  IntersectionQuadratureType interQuad( gridPart_, intersection, faceQuadOrd_ );
509 
510  // get appropriate references
511  const QuadratureImp &faceQuadInner = interQuad.inside();
512  const QuadratureImp &faceQuadOuter = interQuad.outside();
513 
514  // make Entity known in caller
515  caller().setNeighbor(nb, faceQuadInner, faceQuadOuter);
516 
517  // get goemetry of neighbor
518  const Geometry & nbGeo = nb.geometry();
519 
520  // get neighbors volume
521  const double nbvol = nbGeo.volume();
522 
523  const size_t faceQuadInner_nop = faceQuadInner.nop();
524  assert( faceQuadInner.nop() == faceQuadOuter.nop() );
525 
526  // check valNbVec_ here, valEnVec_ might have been resized
527  if( valNbVec_.size() < faceQuadInner_nop )
528  {
529  valEnVec_.resize( faceQuadInner_nop );
530  valNbVec_.resize( faceQuadInner_nop );
531  }
532 
533  for (size_t l = 0; l < faceQuadInner_nop; ++l)
534  {
535  RangeType& fluxEn = valEnVec_[ l ];
536  RangeType& fluxNb = valNbVec_[ l ];
537 
538  wspeedS += caller().numericalFlux(intersection,
539  faceQuadInner,
540  faceQuadOuter,
541  l,
542  fluxEn, fluxNb )
543  * faceQuadInner.weight(l);
544 
545  // apply weights
546  fluxEn *= -faceQuadInner.weight(l);
547  fluxNb *= faceQuadOuter.weight(l);
548  }
549 
550  // add values to local functions
551  updEn.axpyQuadrature( faceQuadInner, valEnVec_ );
552 
553  // update on neighbor
554  if( canUpdateNeighbor )
555  {
556  // init local function
557  initLocalFunction( nb, updNb );
558  // add fluxes
559  updNb.axpyQuadrature( faceQuadOuter, valNbVec_ );
560  // add update to real function
561  updateFunction(nb, updNb );
562  }
563 
564  return nbvol;
565  }
566 
567  protected:
568  DiscreteModelCallerType &caller () const
569  {
570  assert( caller_ );
571  return *caller_;
572  }
573 
574  private:
575  LocalDGPass();
576  LocalDGPass(const LocalDGPass&);
577  LocalDGPass& operator=(const LocalDGPass&);
578 
579  protected:
580  mutable DiscreteModelCallerType *caller_;
581  DiscreteModelType& problem_;
582 
583  mutable ArgumentType* arg_;
584  mutable DestinationType* dest_;
585 
586  const DiscreteFunctionSpaceType& spc_;
587  const GridPartType & gridPart_;
588  const IndexSetType& indexSet_;
589 
590  // indicator for grid walk
592 
593  mutable TemporaryLocalFunctionType updEn_;
594  mutable TemporaryLocalFunctionType updNb_;
595 
600 
601  mutable double dtMin_;
602  const double minLimit_;
603 
605  LocalMassMatrixType localMassMatrix_;
606  const bool notThreadParallel_;
607  };
609 
610  } // namespace Fem
611 
612 } // namespace Dune
613 
614 #endif // #ifndef DUNE_FEM_PASS_LOCALDG_HH
TemporaryLocalFunctionType updNb_
Definition: localdg/pass.hh:594
GridType::template Codim< 0 >::Geometry Geometry
Definition: localdg/pass.hh:82
MutableArray< RangeType > valNbVec_
Definition: localdg/pass.hh:599
void applyLocal(const EntityType &en, const NeighborChecker &nbChecker) const
local integration
Definition: localdg/pass.hh:233
DiscreteFunctionSpaceType::IteratorType IteratorType
Iterator over the space.
Definition: localdg/pass.hh:69
DiscreteModelImp DiscreteModelType
Repetition of template arguments.
Definition: localdg/pass.hh:55
double dtMin_
Definition: localdg/pass.hh:601
GridPartType::IntersectionIteratorType IntersectionIteratorType
Definition: localdg/pass.hh:80
size_t size() const
return number of enties of array
Definition: arrays.hh:213
static constexpr T min(T a)
Definition: utility.hh:81
const GridPartType & gridPart_
Definition: localdg/pass.hh:587
void applyLocalMass(const EntityType &en) const
Definition: localdg/pass.hh:222
double applyLocalNeighbor(const IntersectionType &intersection, const EntityType &en, const EntityType &nb, LocalFunctionImp &updEn, LocalFunctionImp &updNb, double &wspeedS, const bool canUpdateNeighbor) const
Definition: localdg/pass.hh:492
TemporaryLocalFunctionType updEn_
Definition: localdg/pass.hh:593
const DiscreteFunctionSpaceType & spc_
Definition: localdg/pass.hh:586
DiscreteModelType & problem_
Definition: localdg/pass.hh:581
void analyticalFlux(const EntityType &entity, const VolumeQuadratureType &quadrature, const int qp, JacobianRangeType &flux)
Definition: modelcaller.hh:145
void updateFunction(const EntityType &en, LocalFunctionImp &update) const
add update to destination
Definition: localdg/pass.hh:391
void resize(size_t nsize)
Definition: arrays.hh:491
static constexpr T max(T a)
Definition: utility.hh:65
GridPartType::IndexSetType IndexSetType
Definition: localdg/pass.hh:94
size_t nop() const
obtain the number of integration points
Definition: quadratureimp.hh:108
const double minLimit_
Definition: localdg/pass.hh:602
void setTime(double time)
Definition: modelcaller.hh:92
DiscreteModelCallerType * caller_
Definition: localdg/pass.hh:580
const FieldType & weight(size_t i) const
obtain weight of i-th integration point
Definition: quadratureimp.hh:242
LocalMassMatrixType localMassMatrix_
Definition: localdg/pass.hh:605
DGDiscreteModelCaller< DiscreteModelType, ArgumentType, PassIds > DiscreteModelCallerType
Definition: localdg/pass.hh:88
const bool notThreadParallel_
Definition: localdg/pass.hh:606
Definition: localdg/pass.hh:91
DiscreteFunctionSpaceType::JacobianRangeType JacobianRangeType
Definition: localdg/pass.hh:76
specialize with &#39;true&#39; if implementation guarantees conforming level grids. (default=false) ...
Definition: gridpart/common/capabilities.hh:79
IntersectionQuadrature is a helper class for creating the appropriate face quadratures for integratin...
Definition: intersectionquadrature.hh:25
void printTexInfo(std::ostream &out) const
printTex info of operator
Definition: common/pass.hh:210
const int faceQuadOrd_
Definition: localdg/pass.hh:604
MutableArray< bool > visited_
Definition: localdg/pass.hh:591
void updateFunctionAndApplyMass(const EntityType &en, LocalFunctionImp &update) const
add update to destination
Definition: localdg/pass.hh:401
static double max(const Double &v, const double p)
Definition: double.hh:387
EntityType Entity
Definition: local.hh:57
void applyInverse(MassCaller &caller, const EntityType &entity, LocalFunction &lf) const
Definition: localmassmatrix.hh:231
BaseType::ArgumentType ArgumentType
Definition: localdg/pass.hh:61
double boundaryFlux(const IntersectionType &intersection, const FaceQuadratureType &quadrature, const int qp, RangeType &gLeft)
Definition: modelcaller.hh:190
double numericalFlux(const IntersectionType &intersection, const QuadratureType &inside, const QuadratureType &outside, const int qp, RangeType &gLeft, RangeType &gRight)
Definition: modelcaller.hh:180
void evalVolumetricPartBoth(const EntityType &en, const Geometry &geo, const VolumeQuadratureType &volQuad, LocalFunctionImp &updEn) const
Definition: localdg/pass.hh:450
bool operator()(const EntityType &, const EntityType &) const
Definition: localdg/pass.hh:211
MutableArray< RangeType > valEnVec_
Definition: localdg/pass.hh:598
Dune::EntityPointer< Grid, Implementation >::Entity make_entity(const Dune::EntityPointer< Grid, Implementation > &entityPointer)
Definition: compatibility.hh:23
void setBoundary(const EntityType &entity, const QuadratureType &quadrature)
Definition: modelcaller.hh:138
virtual ~LocalDGPass()
Destructor.
Definition: localdg/pass.hh:147
void evalVolumetricPartFlux(const EntityType &en, const Geometry &geo, const VolumeQuadratureType &volQuad, LocalFunctionImp &updEn) const
Definition: localdg/pass.hh:417
void pass(const GlobalArgumentType &arg) const
Definition: common/pass.hh:267
Definition: localdg/pass.hh:209
MutableArray< JacobianRangeType > fMatVec_
Some helper variables.
Definition: localdg/pass.hh:597
DiscreteFunctionSpaceType::DomainType DomainType
Definition: localdg/pass.hh:74
Definition: coordinate.hh:4
virtual void prepare(const ArgumentType &arg, DestinationType &dest) const
Definition: localdg/pass.hh:168
void applyLocal(const EntityType &en, TemporaryLocalFunctionType &updEn, const NeighborChecker &nbChecker) const
local integration
Definition: localdg/pass.hh:247
DiscreteFunctionSpaceType::GridPartType GridPartType
Definition: localdg/pass.hh:73
TemporaryLocalFunction< DiscreteFunctionSpaceType > TemporaryLocalFunctionType
Definition: localdg/pass.hh:95
LocalMassMatrix< DiscreteFunctionSpaceType, VolumeQuadratureType > LocalMassMatrixType
type of local mass matrix
Definition: localdg/pass.hh:98
DiscreteModelType::Traits::DiscreteFunctionSpaceType DiscreteFunctionSpaceType
Definition: localdg/pass.hh:67
double time() const
return current time of calculation
Definition: common/pass.hh:254
STL namespace.
BaseType::Entity EntityType
Definition: localdg/pass.hh:60
void setMemoryFactor(const double memFactor)
set memory factor
Definition: arrays.hh:465
const int volumeQuadOrd_
Definition: localdg/pass.hh:604
BaseType::TotalArgumentType ArgumentType
The type of the argument (and destination) type of the overall operator.
Definition: local.hh:45
ArgumentType * arg_
Definition: localdg/pass.hh:583
LocalDGPass(DiscreteModelType &problem, PreviousPassType &pass, const DiscreteFunctionSpaceType &spc, const int volumeQuadOrd=-1, const int faceQuadOrd=-1, const bool notThreadParallel=true)
Definition: localdg/pass.hh:109
model caller for local DG pass
Definition: modelcaller.hh:30
DiscreteModelType::Traits::DestinationType DestinationType
Definition: localdg/pass.hh:64
void applyLocal(const EntityType &en) const
Definition: localdg/pass.hh:217
DiscreteModelType::Traits::VolumeQuadratureType VolumeQuadratureType
Definition: localdg/pass.hh:65
Generic implementation of a Dune quadrature.
Definition: quadratureimp.hh:178
virtual void finalize(const ArgumentType &arg, DestinationType &dest) const
Some timestep size management.
Definition: localdg/pass.hh:194
IntersectionIteratorType::Intersection IntersectionType
Definition: localdg/pass.hh:81
LocalPass< DiscreteModelImp, PreviousPassImp, passIdImp > BaseType
Base class.
Definition: localdg/pass.hh:49
DestinationType * dest_
Definition: localdg/pass.hh:584
DiscreteFunctionSpaceType::BasisFunctionSetType BasisFunctionSetType
Definition: localdg/pass.hh:77
DiscreteFunctionSpaceType::RangeType RangeType
Definition: localdg/pass.hh:75
PreviousPassImp PreviousPassType
Repetition of template arguments.
Definition: localdg/pass.hh:57
double timeStepEstimateImpl() const
Estimate for the timestep size.
Definition: localdg/pass.hh:158
void printTexInfo(std::ostream &out) const
print tex info
Definition: localdg/pass.hh:151
DestinationType::LocalFunctionType LocalFunctionType
Definition: localdg/pass.hh:86
DiscreteModelCallerType & caller() const
Definition: localdg/pass.hh:568
DiscreteModelType::Traits::FaceQuadratureType FaceQuadratureType
Definition: localdg/pass.hh:66
void initLocalFunction(const EntityType &en, LocalFunctionImp &update) const
Definition: localdg/pass.hh:381
const IndexSetType & indexSet_
Definition: localdg/pass.hh:588
void axpyQuadrature(const QuadratureType &quad, const VectorType &values)
evaluate all basisfunctions for all quadrature points, multiply with the given factor and add the res...
Definition: localfunction.hh:361
double analyticalFluxAndSource(const EntityType &entity, const VolumeQuadratureType &quadrature, const int qp, JacobianRangeType &flux, RangeType &source)
Definition: modelcaller.hh:166
static double min(const Double &v, const double p)
Definition: double.hh:375
void setEntity(const EntityType &entity)
Definition: modelcaller.hh:97
Specialisation of Pass which provides a grid walk-through, but leaves open what needs to be done on e...
Definition: local.hh:32
size_t numberOfElements() const
return number of elements visited during operator computation
Definition: localdg/pass.hh:207
DiscreteFunctionSpaceType::GridType GridType
Definition: localdg/pass.hh:72
void setNeighbor(const EntityType &entity)
Definition: modelcaller.hh:104
BaseType::PassIds PassIds
pass ids up to here (tuple of integral constants)
Definition: localdg/pass.hh:52
Dune::PushBackTuple< typename PreviousPassType::PassIds, std::integral_constant< int, passIdImp > >::type PassIds
pass ids up to here (tuple of integral constants)
Definition: common/pass.hh:157
Definition: localdg/pass.hh:41