1 #ifndef DUNE_FEM_PASS_LOCALDG_HH 2 #define DUNE_FEM_PASS_LOCALDG_HH 38 template <
class DiscreteModelImp,
class PreviousPassImp ,
int passIdImp = -1 >
42 public LocalPass< DiscreteModelImp , PreviousPassImp , passIdImp >
69 typedef typename DiscreteFunctionSpaceType::IteratorType
IteratorType;
72 typedef typename DiscreteFunctionSpaceType::GridType
GridType;
73 typedef typename DiscreteFunctionSpaceType::GridPartType
GridPartType;
74 typedef typename DiscreteFunctionSpaceType::DomainType
DomainType;
75 typedef typename DiscreteFunctionSpaceType::RangeType
RangeType;
82 typedef typename GridType::template Codim<0>::Geometry
Geometry;
91 enum {
dimRange = DiscreteFunctionSpaceType::dimRange };
110 PreviousPassType&
pass,
111 const DiscreteFunctionSpaceType& spc,
112 const int volumeQuadOrd =-1,
113 const int faceQuadOrd=-1,
114 const bool notThreadParallel =
true ) :
132 ( 2 *
spc_.order()) : volumeQuadOrd ),
134 ( 2 *
spc_.order()+1) : faceQuadOrd ),
153 out <<
"LocalDGPass: " 161 const double p = 2 *
spc_.order() + 1;
168 virtual void prepare(
const ArgumentType& arg, DestinationType& dest)
const 170 arg_ =
const_cast<ArgumentType*
>(&arg);
187 for(
int i=0; i<indSize; ++i)
visited_[i] =
false;
194 virtual void finalize(
const ArgumentType& arg, DestinationType& dest)
const 199 spc_.communicate( dest );
225 LocalFunctionType
function =
dest_->localFunction(en);
232 template <
class NeighborChecker>
233 void applyLocal(
const EntityType& en,
const NeighborChecker& nbChecker )
const 246 template <
class NeighborChecker>
248 TemporaryLocalFunctionType& updEn,
249 const NeighborChecker& nbChecker )
const 252 const Geometry & geo = en.geometry();
255 const double vol = geo.volume();
290 IntersectionIteratorType endnit =
gridPart_.iend(en);
291 for (IntersectionIteratorType nit =
gridPart_.ibegin(en); nit != endnit; ++nit)
294 const IntersectionType& intersection = *nit;
298 double wspeedS = 0.0;
299 if (intersection.neighbor())
302 EntityType nb =
make_entity( intersection.outside() );
304 const bool canUpdateNeighbor = nbChecker( en, nb );
310 && !intersection.conforming() )
314 nbvol = applyLocalNeighbor< false >
324 nbvol = applyLocalNeighbor< true >
334 else if( intersection.boundary() )
337 FaceQuadratureType::INSIDE);
343 const size_t faceQuadInner_nop = faceQuadInner.nop();
349 for (
size_t l = 0; l < faceQuadInner_nop; ++l)
355 * faceQuadInner.weight(l);
358 flux *= -faceQuadInner.weight(l);
368 double minvolS =
std::min(vol,nbvol);
380 template <
class LocalFunctionImp>
390 template <
class LocalFunctionImp>
392 LocalFunctionImp& update)
const 395 LocalFunctionType
function =
dest_->localFunction(en);
400 template <
class LocalFunctionImp>
402 const EntityType& en,
403 LocalFunctionImp& update)
const 406 LocalFunctionType
function =
dest_->localFunction(en);
416 template <
class LocalFunctionImp>
418 const Geometry& geo ,
419 const VolumeQuadratureType& volQuad,
420 LocalFunctionImp& updEn )
const 422 const size_t volQuad_nop = volQuad.nop();
428 for (
size_t l = 0; l < volQuad_nop; ++l)
430 JacobianRangeType& flux =
fMatVec_[ l ];
435 const double intel = geo.integrationElement(volQuad.point(l))
443 updEn.axpyQuadrature (volQuad,
fMatVec_ );
449 template <
class LocalFunctionImp>
452 const VolumeQuadratureType& volQuad,
453 LocalFunctionImp& updEn )
const 455 const size_t volQuad_nop = volQuad.nop();
467 for (
size_t l = 0; l < volQuad_nop; ++l)
469 JacobianRangeType& flux =
fMatVec_[ l ];
476 const double intel = geo.integrationElement(volQuad.point(l))
491 template <
bool conforming,
class LocalFunctionImp>
493 const EntityType &en,
494 const EntityType &nb,
495 LocalFunctionImp &updEn,
496 LocalFunctionImp &updNb,
498 const bool canUpdateNeighbor )
const 501 assert( intersection.conforming() == conforming );
505 typedef typename IntersectionQuadratureType :: FaceQuadratureType
QuadratureImp;
511 const QuadratureImp &faceQuadInner = interQuad.inside();
512 const QuadratureImp &faceQuadOuter = interQuad.outside();
518 const Geometry & nbGeo = nb.geometry();
521 const double nbvol = nbGeo.volume();
523 const size_t faceQuadInner_nop = faceQuadInner.
nop();
524 assert( faceQuadInner.
nop() == faceQuadOuter.
nop() );
533 for (
size_t l = 0; l < faceQuadInner_nop; ++l)
543 * faceQuadInner.
weight(l);
546 fluxEn *= -faceQuadInner.
weight(l);
547 fluxNb *= faceQuadOuter.
weight(l);
551 updEn.axpyQuadrature( faceQuadInner,
valEnVec_ );
554 if( canUpdateNeighbor )
559 updNb.axpyQuadrature( faceQuadOuter,
valNbVec_ );
586 const DiscreteFunctionSpaceType&
spc_;
593 mutable TemporaryLocalFunctionType
updEn_;
594 mutable TemporaryLocalFunctionType
updNb_;
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 'true' 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
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