1 #ifndef DUNE_FEM_LOCALDG_DISCRETEMODEL_HH 2 #define DUNE_FEM_LOCALDG_DISCRETEMODEL_HH 6 #include <dune/common/bartonnackmanifcheck.hh> 32 template <
class DGDiscreteModelTraits>
36 typedef DGDiscreteModelTraits
Traits;
45 typedef typename FunctionSpaceType::DomainType
DomainType;
47 typedef typename FunctionSpaceType::RangeType
RangeType;
54 typedef typename DiscreteFunctionSpaceType::GridPartType
GridPartType;
56 typedef typename DiscreteFunctionSpaceType::EntityType
EntityType;
64 enum {
dimRange = DiscreteFunctionSpaceType::dimRange };
69 template<
template <
class,
int,
int >
class MatrixType,
70 class ctype,
int dimD,
int dimR
74 typedef MatrixType< RangeFieldType, dimRange, dimRange >
Type;
119 template <
class ArgumentTuple,
class FaceDomainType>
122 const FaceDomainType &x,
123 const ArgumentTuple &uLeft,
124 const ArgumentTuple &uRight,
128 CHECK_INTERFACE_IMPLEMENTATION(
asImp().
numericalFlux( intersection, time, x, uLeft, uRight, gLeft, gRight) );
129 return asImp().numericalFlux( intersection, time, x, uLeft, uRight, gLeft, gRight );
147 template <
class ArgumentTuple,
class FaceDomainType>
150 const FaceDomainType &x,
151 const ArgumentTuple &uLeft,
154 CHECK_INTERFACE_IMPLEMENTATION(
asImp().
boundaryFlux( intersection, time, x, uLeft, gLeft ) );
155 return asImp().boundaryFlux( intersection, time, x, uLeft, gLeft );
169 template <
class ArgumentTuple>
172 const LocalCoordinateType& x,
173 const ArgumentTuple& u,
174 JacobianRangeType& f )
176 CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
193 template <
class ArgumentTuple,
class JacobianTuple>
196 const LocalCoordinateType& x,
197 const ArgumentTuple& u,
198 const JacobianTuple& jac,
201 CHECK_INTERFACE_IMPLEMENTATION(
203 return asImp().source(en, time, x, u, jac, s);
216 template <
class ArgumentTuple>
217 void mass (
const EntityType& en,
219 const LocalCoordinateType& x,
220 const ArgumentTuple& u,
223 CHECK_AND_CALL_INTERFACE_IMPLEMENTATION(
250 DGDiscreteModelType&
asImp() {
return static_cast<DGDiscreteModelType&
>(*this); }
251 const DGDiscreteModelType&
asImp()
const {
return static_cast<const DGDiscreteModelType&
>(*this); }
265 template<
class DGDiscreteModelTraits
319 inline bool hasFlux ()
const {
return false; }
331 inline bool hasMass ()
const {
return false; }
336 template <
class ArgumentTuple,
class FaceDomainType>
339 const FaceDomainType& x,
340 const ArgumentTuple& uLeft,
341 const ArgumentTuple& uRight,
354 template <
class ArgumentTuple,
class FaceDomainType>
357 const FaceDomainType &x,
358 const ArgumentTuple &uLeft,
369 template<
class ArgumentTuple >
372 const LocalCoordinateType &x,
373 const ArgumentTuple &u,
374 JacobianRangeType &f )
383 template <
class ArgumentTuple,
class JacobianTuple>
386 const LocalCoordinateType& x,
387 const ArgumentTuple& u,
388 const JacobianTuple& jac,
399 template <
class ArgumentTuple>
400 void mass (
const EntityType& en,
402 const LocalCoordinateType& x,
403 const ArgumentTuple& u,
406 enum { rows = MassFactorType :: rows };
407 enum { cols = MassFactorType :: cols };
409 for(
int i=0; i<rows; ++i)
414 for(
int j=i+1; j<cols; ++j)
415 m[i][j] = m[j][i] = 0;
435 template <
class Adaptation>
437 const double weight = 1.0 ) {}
440 template <
class Adaptation,
class ThreadFilter>
443 const double weight = 1.0 ) {}
456 template <
class DGDiscreteModelTraits
468 :
public DGDiscreteModelDefault< DGDiscreteModelTraits, N1, N2, N3, N4, N5, N6, N7, N8, N9 >,
477 : enVol_(-1.0) , nbVol_(-1.0) , en_(0) , nb_(0)
481 : enVol_(other.enVol_) , nbVol_(other.nbVol_),
482 en_(other.en_) , nb_(other.nb_)
492 enVol_ = en.geometry().volume();
502 nbVol_ = nb.geometry().volume();
528 assert(enVol_ > 0.0);
535 assert( nbVol_ > 0.0 );
543 const EntityType* en_;
544 const EntityType* nb_;
551 #endif // #ifndef DUNE_FEM_LOCALDG_DISCRETEMODEL_HH bool hasSource() const
Definition: discretemodel.hh:325
double source(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, const JacobianTuple &jac, RangeType &s)
Empty implementation that fails if problem claims to have a source term.
Definition: discretemodel.hh:384
Dune::Fem::Selector< N1, N2, N3, N4, N5, N6, N7, N8, N9 >::Type Selector
Definition: discretemodel.hh:305
void setEntity(const EntityType &en)
empty implementation
Definition: discretemodel.hh:420
double enVolume() const
return volume of entity
Definition: discretemodel.hh:526
BaseType::EntityType EntityType
Definition: discretemodel.hh:474
void setNeighbor(EntityType &nb)
Passes the active neigbor entity to the model. This can be used, to set local functions required as d...
Definition: discretemodel.hh:244
Definition: discretemodel.hh:467
Traits::DGDiscreteModelType DGDiscreteModelType
Implementation type for Barton-Nackman trick.
Definition: discretemodel.hh:38
void analyticalFlux(const EntityType &entity, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, JacobianRangeType &f)
Empty implementation that fails if problem claims to have a flux contribution.
Definition: discretemodel.hh:370
void removeAdaptation()
remove pointer to adaptation handle
Definition: discretemodel.hh:446
void setEntity(const EntityType &en)
method setting pointer of inside entity and getting volume
Definition: discretemodel.hh:489
const DGDiscreteModelType & asImp() const
Definition: discretemodel.hh:251
ExtractMatrix< JacobianRangeType >::Type MassFactorType
Type of mass factor.
Definition: discretemodel.hh:78
void analyticalFlux(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, JacobianRangeType &f)
Computes the analytical flux of the problem. Analytical flux of the problem.
Definition: discretemodel.hh:170
Mass factor type.
Definition: discretemodel.hh:67
BaseType::JacobianRangeType JacobianRangeType
Definition: discretemodel.hh:312
DiscreteFunctionSpaceType::EntityType EntityType
Element (codim 0 entity) of the grid.
Definition: discretemodel.hh:56
void setNeighbor(const EntityType &nb)
method seting pointer of outside entity and getting volume
Definition: discretemodel.hh:499
BaseType::EntityType EntityType
Definition: discretemodel.hh:308
double numericalFlux(const IntersectionType &it, const double time, const FaceDomainType &x, const ArgumentTuple &uLeft, const ArgumentTuple &uRight, RangeType &gLeft, RangeType &gRight)
Empty implementation that fails if problem claims to have a flux contribution.
Definition: discretemodel.hh:337
Interface class for problem definition in the LDG context.
Definition: discretemodel.hh:33
FunctionSpaceType::RangeType RangeType
Vector type of the function space's range.
Definition: discretemodel.hh:47
DGDiscreteModelType & asImp()
Definition: discretemodel.hh:250
Default implementation of the DGDiscreteModelInterface where methods for the fluxes and the source te...
Definition: discretemodel.hh:276
BaseType::IntersectionType IntersectionType
Definition: discretemodel.hh:309
bool hasMass() const
Returns true if problem has a mass matrix factor disctinct from the identity.
Definition: discretemodel.hh:94
bool hasSource() const
Returns true if problem has a source term. If you program this method to return true, make sure to implement method source well.
Definition: discretemodel.hh:91
BaseType::LocalCoordinateType LocalCoordinateType
Definition: discretemodel.hh:310
void setEntity(EntityType &en)
Passes the active entity to the model. This can be used, to set local functions required as data func...
Definition: discretemodel.hh:233
Traits::DiscreteFunctionSpaceType DiscreteFunctionSpaceType
Discrete function space.
Definition: discretemodel.hh:41
DiscreteFunctionSpaceType::GridPartType GridPartType
Type of GridPart.
Definition: discretemodel.hh:54
double numericalFlux(const IntersectionType &intersection, const double time, const FaceDomainType &x, const ArgumentTuple &uLeft, const ArgumentTuple &uRight, RangeType &gLeft, RangeType &gRight)
Computes the numerical flux at a cell interface.
Definition: discretemodel.hh:120
FunctionSpaceType::DomainType DomainType
Coordinate type (world coordinates)
Definition: discretemodel.hh:45
Definition: coordinate.hh:4
FunctionSpaceType::JacobianRangeType JacobianRangeType
Matrix type of the function space's jacobian.
Definition: discretemodel.hh:51
DGDiscreteModelTraits Traits
Traits class defined by the user.
Definition: discretemodel.hh:36
bool hasMass() const
Definition: discretemodel.hh:331
void setAdaptation(Adaptation &, const double weight=1.0)
default method for setting adaptation handle to discrete model
Definition: discretemodel.hh:436
Dune::Fem::ElementTuple< N1, N2, N3, N4, N5, N6, N7, N8, N9 >::Type Type
tuple consisting of Dune::integral_constant< int, N_i >
Definition: selector.hh:71
DGDiscreteModelDefaultWithInsideOutside(const DGDiscreteModelDefaultWithInsideOutside &other)
Definition: discretemodel.hh:480
double boundaryFlux(const IntersectionType &intersection, const double time, const FaceDomainType &x, const ArgumentTuple &uLeft, RangeType &gLeft)
Computes the flux at the boundary Special kind of numerical flux. The intersection iterator provides ...
Definition: discretemodel.hh:148
DiscreteFunctionSpaceType::FunctionSpaceType FunctionSpaceType
Function space type.
Definition: discretemodel.hh:43
Definition: threadfilter.hh:20
double boundaryFlux(const IntersectionType &intersection, const double time, const FaceDomainType &x, const ArgumentTuple &uLeft, RangeType &gLeft)
Empty implementation that fails if problem claims to have a flux contribution.
Definition: discretemodel.hh:355
double source(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, const JacobianTuple &jac, RangeType &s)
Implements the source term of the problem. Source term contribution. The source term can contain part...
Definition: discretemodel.hh:194
BaseType::MassFactorType MassFactorType
Definition: discretemodel.hh:313
DGDiscreteModelDefaultWithInsideOutside()
Definition: discretemodel.hh:476
EntityType::Geometry::LocalCoordinate LocalCoordinateType
Type of local coordinate.
Definition: discretemodel.hh:61
const EntityType & inside() const
method returning reference to inside entity
Definition: discretemodel.hh:509
bool hasFlux() const
Definition: discretemodel.hh:319
void mass(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, MassFactorType &m)
Implements the mass factor term of the problem. Mass term contribution.
Definition: discretemodel.hh:217
Definition: discretemodel.hh:64
void mass(const EntityType &en, const double time, const LocalCoordinateType &x, const ArgumentTuple &u, MassFactorType &m)
empty implementation for mass factor default implementation sets this factor to 1.0
Definition: discretemodel.hh:400
bool hasFlux() const
Returns true if problem has a flux contribution. If you program this method to return true...
Definition: discretemodel.hh:85
BaseType::RangeType RangeType
Definition: discretemodel.hh:311
const EntityType & outside() const
method returning reference to outside entity
Definition: discretemodel.hh:519
DiscreteFunctionSpaceType::IntersectionType IntersectionType
Intersection type.
Definition: discretemodel.hh:58
double nbVolume() const
return volume of neighbor
Definition: discretemodel.hh:533
void setAdaptation(Adaptation &, const ThreadFilter &, const double weight=1.0)
default method for setting adaptation handle and thead filter to discrete model
Definition: discretemodel.hh:441
FunctionSpaceType::RangeFieldType RangeFieldType
Vector type of the function space's range field type.
Definition: discretemodel.hh:49
void setNeighbor(const EntityType &nb)
empty implementation
Definition: discretemodel.hh:424
Definition: discretemodel.hh:431
MatrixType< RangeFieldType, dimRange, dimRange > Type
Definition: discretemodel.hh:74