diff --git a/dune/fem-dg/operator/limiter/limitpass.hh b/dune/fem-dg/operator/limiter/limitpass.hh
index 3ec362e9aa22367d23db8d22000040f0ebba6ff6..9da1b1df8f4a4e9b2ea483cc2b54e7bd0a375d76 100644
--- a/dune/fem-dg/operator/limiter/limitpass.hh
+++ b/dune/fem-dg/operator/limiter/limitpass.hh
@@ -603,7 +603,6 @@ namespace Dune {
// Types from the base class
typedef typename BaseType::Entity EntityType;
- typedef const EntityType ConstEntityType;
typedef typename BaseType::ArgumentType ArgumentType;
typedef typename PreviousPassType::GlobalArgumentType ArgumentFunctionType;
@@ -659,7 +658,8 @@ namespace Dune {
typedef FieldVector< DomainType , dimRange > DeoModType;
typedef FieldMatrix< DomainFieldType, dimDomain , dimDomain > MatrixType;
- typedef Fem::AllGeomTypes< typename GridPartType :: IndexSetType,
+ typedef typename GridPartType :: IndexSetType IndexSetType;
+ typedef Fem::AllGeomTypes< IndexSetType,
GridType> GeometryInformationType;
typedef Fem::GeometryInformation< GridType, 1 > FaceGeometryInformationType;
@@ -930,6 +930,7 @@ namespace Dune {
dest_(0),
spc_(spc),
gridPart_(spc_.gridPart()),
+ indexSet_( gridPart_.indexSet() ),
localIdSet_( gridPart_.grid().localIdSet()),
lagrangePointSetContainer_(gridPart_),
orderPower_( -((spc_.order()+1.0) * 0.25)),
@@ -956,7 +957,7 @@ namespace Dune {
admissibleFunctions_( getAdmissibleFunctions() ),
usedAdmissibleFunctions_( admissibleFunctions_ )
{
- if( gridPart_.grid().comm().rank() == 0 )
+ if( gridPart_.comm().rank() == 0 )
{
std::cout << "LimitPass: Grid is ";
if( cartesianGrid_ )
@@ -1136,6 +1137,10 @@ namespace Dune {
// calculate maximal indicator (if necessary)
discreteModel_.indicatorMax();
+ // reset visited vector
+ visited_.resize( indexSet_.size( 0 ) );
+ std::fill( visited_.begin(), visited_.end(), false );
+
const int numLevels = gridPart_.grid().maxLevel() + 1;
// check size of matrix cache vec
if( (int) matrixCacheVec_.size() < numLevels )
@@ -1169,13 +1174,53 @@ namespace Dune {
protected:
//! apply local is virtual
- void applyLocal(ConstEntityType& en) const
+ void applyLocal(const EntityType& en) const
{
applyLocalImp(en);
}
+ template <class NeighborChecker>
+ void applyLocalInterior( const EntityType& entity,
+ const NeighborChecker& nbChecker ) const
+ {
+ // check whether on of the intersections is with ghost element
+ // and if so, skip the computation of the limited solution for now
+ const IntersectionIteratorType endnit = gridPart_.iend(entity);
+ for (IntersectionIteratorType nit = gridPart_.ibegin(entity); nit != endnit; ++nit)
+ {
+ const IntersectionType& intersection = *nit;
+ if( intersection.neighbor() )
+ {
+ // get neighbor
+ EntityPointerType outside = intersection.outside();
+ const EntityType & nb = * outside;
+
+ // check whether we have to skip this intersection
+ if( nbChecker.skipIntersection( nb ) )
+ {
+ return ;
+ }
+ }
+ }
+
+ // otherwise apply limiting process
+ applyLocalImp( entity );
+ }
+
+ template <class NeighborChecker>
+ void applyLocalProcessBoundary( const EntityType& entity,
+ const NeighborChecker& nChecker ) const
+ {
+ assert( indexSet_.index( entity ) < int(visited_.size()) );
+ // if entity was already visited, do nothing in this turn
+ if( visited_[ indexSet_.index( entity ) ] ) return ;
+
+ // apply limiter otherwise
+ applyLocalImp( entity );
+ }
+
//! Perform the limitation on all elements.
- void applyLocalImp(ConstEntityType& en) const
+ void applyLocalImp(const EntityType& en) const
{
// timer for shock detection
Timer indiTime;
@@ -1415,6 +1460,10 @@ namespace Dune {
discreteModel_.adaptation( gridPart_.grid() , en, shockIndicator, adaptIndicator );
}
+ // mark entity as finished, even if not limited everything necessary was done
+ assert( indexSet_.index( entity ) < int(visited_.size()) );
+ visited_[ indexSet_.index( en ) ] = true ;
+
// if nothing to limit then just return here
if ( ! limiter ) return ;
@@ -1461,6 +1510,7 @@ namespace Dune {
assert( checkPhysical(en, geo, limitEn) );
stepTime_[1] += indiTime.elapsed();
+
//end limiting process
}
@@ -2549,6 +2599,7 @@ namespace Dune {
const DiscreteFunctionSpaceType& spc_;
GridPartType& gridPart_;
+ const IndexSetType& indexSet_;
const LocalIdSetType& localIdSet_;
LagrangePointSetContainerType lagrangePointSetContainer_;
@@ -2581,6 +2632,9 @@ namespace Dune {
mutable std::vector< RangeType > nbVals_;
mutable std::vector< MatrixCacheType > matrixCacheVec_;
+ // vector for stroing the information which elements have been computed already
+ mutable std::vector< bool > visited_;
+
LocalMassMatrixType localMassMatrix_;
//! true if limiter is used in adaptive scheme
const bool adaptive_;
diff --git a/dune/fem-dg/pass/dgpass.hh b/dune/fem-dg/pass/dgpass.hh
index 2ff8aa79ae6b6074a3292091bbcc364410f5606d..489d3d3cb07d98f91d4539d903e269c8969dcd79 100644
--- a/dune/fem-dg/pass/dgpass.hh
+++ b/dune/fem-dg/pass/dgpass.hh
@@ -142,7 +142,7 @@ namespace Dune {
{
// make sure that either a ghost layer or an overlap layer is there for
// communication of the data, otherwise the scheme will not work
- if( spc_.grid().comm().size() > 1 )
+ if( spc_.gridPart().comm().size() > 1 )
{
if( spc_.grid().ghostSize( 0 ) <= 0 && spc_.grid().overlapSize( 0 ) <= 0 )
{
@@ -396,7 +396,7 @@ namespace Dune {
// only calculate element integral part and interior fluxes
template <class NeighborChecker>
- void interiorIntegral(const EntityType& entity, const NeighborChecker& nbChecker ) const
+ void applyLocalInterior(const EntityType& entity, const NeighborChecker& nbChecker ) const
{
// init local function
initLocalFunction( entity , updEn_ );
@@ -413,8 +413,8 @@ namespace Dune {
// only calculate integral that need data from other processes
template <class NeighborChecker>
- void processBoundaryIntegral(const EntityType& entity,
- const NeighborChecker& nbChecker ) const
+ void applyLocalProcessBoundary(const EntityType& entity,
+ const NeighborChecker& nbChecker ) const
{
// init local function
initLocalFunction( entity , updEn_ );
diff --git a/dune/fem-dg/pass/threadpass.hh b/dune/fem-dg/pass/threadpass.hh
index 35ffcb21d5bb3b9565a372bbd12fe869c5b7f46a..f89737b9245580664aa9e91425e7c4d051d52d1a 100644
--- a/dune/fem-dg/pass/threadpass.hh
+++ b/dune/fem-dg/pass/threadpass.hh
@@ -549,7 +549,7 @@ namespace Dune {
for (Iterator it = iterators_.begin(); it != endit; ++it)
{
assert( iterators_.thread( *it ) == thread );
- myPass.interiorIntegral( *it, nbChecker );
+ myPass.applyLocalInterior( *it, nbChecker );
}
// receive ghost data (only master thread)
@@ -568,7 +568,7 @@ namespace Dune {
for (Iterator it = iterators_.begin(); it != endit; ++it)
{
assert( iterators_.thread( *it ) == thread );
- myPass.processBoundaryIntegral( *it, nbChecker );
+ myPass.applyLocalProcessBoundary( *it, nbChecker );
}
assert( arg_ );