diff --git a/dune/alugrid/3d/entity.hh b/dune/alugrid/3d/entity.hh index 4be3abf797b3dc274be7babe35b506de331a722a..c745a34f6f3810bd2f27e0fcec032cb6e4365e10 100644 --- a/dune/alugrid/3d/entity.hh +++ b/dune/alugrid/3d/entity.hh @@ -46,7 +46,7 @@ class ALU3dGridEntity : public EntityDefaultImplementation <cd,dim,GridImp,ALU3dGridEntity> { // default just returns level - template <class GridType, int cdim> + template <class GridType, int dm, int cdim> struct GetLevel { template <class ItemType> @@ -60,7 +60,7 @@ public EntityDefaultImplementation <cd,dim,GridImp,ALU3dGridEntity> // this the maximum of levels of elements that have this vertex as sub // entity template <class GridType> - struct GetLevel<GridType,dim> + struct GetLevel<GridType,dim,dim> { template <class ItemType> static int getLevel(const GridType & grid, const ItemType & item) diff --git a/dune/alugrid/3d/entity_imp.cc b/dune/alugrid/3d/entity_imp.cc index 1a693ad3c6586c8623e193ed8f53c310766bf627..1aa91c97917f3282edf335af945a71bed1ca7a13 100644 --- a/dune/alugrid/3d/entity_imp.cc +++ b/dune/alugrid/3d/entity_imp.cc @@ -54,7 +54,7 @@ namespace Dune { alu_inline void ALU3dGridEntity<cd,dim,GridImp> :: setElement(const HItemType & item, const GridImp& grid) { - setElement( item, GetLevel<GridImp,cd>::getLevel(grid,item) ); + setElement( item, GetLevel<GridImp,dim,cd>::getLevel(grid,item) ); } template<int cd, int dim, class GridImp> diff --git a/dune/alugrid/3d/iterator.hh b/dune/alugrid/3d/iterator.hh index 056d274a398578117f4382689754c5e5794eb252..394aa98301bea95206ceb736c69eedb27721c848 100644 --- a/dune/alugrid/3d/iterator.hh +++ b/dune/alugrid/3d/iterator.hh @@ -400,7 +400,7 @@ public: typedef typename InternalIteratorType :: val_t val_t; // here the items level will do - template <class GridImp, int codim> + template <class GridImp, int dim, int codim> class GetLevel { public: @@ -413,8 +413,8 @@ public: }; // level is not needed for codim = 0 - template <class GridImp> - class GetLevel<GridImp,0> + template <class GridImp, int dim> + class GetLevel<GridImp,dim,0> { public: template <class ItemType> @@ -424,8 +424,8 @@ public: } }; - template <class GridImp> - class GetLevel<GridImp,3> + template <class GridImp, int dim> + class GetLevel<GridImp, dim, dim> { public: template <class ItemType> @@ -473,7 +473,7 @@ protected: if( item.first ) { it.updateEntityPointer( item.first , - GetLevel<GridImp,codim>::getLevel(grid, *(item.first) , level) ); + GetLevel<GridImp,GridImp::dimension,codim>::getLevel(grid, *(item.first) , level) ); } else it.updateGhostPointer( *item.second ); diff --git a/examples/callback/adaptation.hh b/examples/callback/adaptation.hh index 260bed08567937c6812baaa7753bac80df003d4b..f35e6bf1b407ed8b41f43fcf75c683f38e362bc4 100644 --- a/examples/callback/adaptation.hh +++ b/examples/callback/adaptation.hh @@ -74,6 +74,7 @@ public: balanceStep_( balanceStep ), balanceCounter_( balanceCounter < 0 ? balanceStep-1 : balanceCounter ), adaptTime_( 0.0 ), + restProlTime_( 0.0 ), lbTime_( 0.0 ), commTime_( 0.0 ) {} @@ -89,6 +90,8 @@ public: //! return time spent for the last adapation in sec double adaptationTime() const { return adaptTime_; } + //! return time spent for the last adapation in sec + double restProlTime() const { return restProlTime_; } //! return time spent for the last load balancing in sec double loadBalanceTime() const { return lbTime_; } //! return time spent for the last communication in sec @@ -154,6 +157,7 @@ private: int balanceCounter_; double adaptTime_; + double restProlTime_; double lbTime_; double commTime_; }; @@ -170,6 +174,7 @@ inline void LeafAdaptation< Grid, Vector >::operator() ( Vector &solution ) adaptTime_ = 0.0; lbTime_ = 0.0; commTime_ = 0.0; + restProlTime_ = 0.0; // reset timer adaptTimer_.reset() ; @@ -194,9 +199,15 @@ inline void LeafAdaptation< Grid, Vector >::operator() ( Vector &solution ) hierarchicRestrict( *it, container_ ); } + restProlTime_ = adaptTimer_.elapsed(); + adaptTimer_.reset(); + // adapt grid, returns true if new elements were created const bool refined = grid_.adapt(); + adaptTime_ = adaptTimer_.elapsed(); + adaptTimer_.reset(); + // interpolate all new cells to leaf level if( refined ) { @@ -205,9 +216,15 @@ inline void LeafAdaptation< Grid, Vector >::operator() ( Vector &solution ) for( LevelIterator it = grid_.template lbegin< 0, partition >( 0 ); it != end; ++it ) hierarchicProlong( *it, container_ ); } + + restProlTime_ += adaptTimer_.elapsed(); + #else // CALLBACK_ADAPTATION // callback adaptation, see interface methods above grid_.adapt( *this ); + + adaptTime_ = adaptTimer_.elapsed(); + #endif // CALLBACK_ADAPTATION // reset adaptation information in grid diff --git a/examples/callback/main.cc b/examples/callback/main.cc index 280f2ec14cd00eaba60ec364090f6b22ff8de97b..25d0227972d888fe6782ed0e9525b028f6b0af8f 100644 --- a/examples/callback/main.cc +++ b/examples/callback/main.cc @@ -225,7 +225,8 @@ void method ( int problem, int startLvl, int maxLvl, adaptation.adaptationTime(), // time for adaptation adaptation.loadBalanceTime(), // time for load balance overallTimer.elapsed(), // time step overall time - getMemoryUsage() ); // memory usage + adaptation.restProlTime(), // adaptation restrict prolong + getMemoryUsage() ); // memory usage } diff --git a/examples/communication/main.cc b/examples/communication/main.cc index bc120069dfb805d40d242c8ae81f5c21e67dbcd8..feee11eb168619bd649e9ce178fa373cc80b36e4 100644 --- a/examples/communication/main.cc +++ b/examples/communication/main.cc @@ -239,7 +239,8 @@ void method ( int problem, int startLvl, int maxLvl, adaptation.adaptationTime(), // time for adaptation adaptation.loadBalanceTime(), // time for load balance overallTimer.elapsed(), // time step overall time - getMemoryUsage() ); // memory usage + adaptation.restProlTime(), // adaptation restrict prolong + getMemoryUsage() ); // memory usage } diff --git a/examples/dgf/ffs2d.dgf b/examples/dgf/ffs2d.dgf index c7ae8c29f9afd09fdb39792cea46706f7bc2422e..2a9d62eabf47421922f32ebb6248747f5bde54a3 100644 --- a/examples/dgf/ffs2d.dgf +++ b/examples/dgf/ffs2d.dgf @@ -14,7 +14,8 @@ simplex # boundarydomain -default 3 % default is 3 = reflection -1 -1 -1 0 1 % inflow -2 3 0 4 1 % outflow +default 3 % default is 3 = reflection +1 -1 -1 0 1 % inflow +2 3 0 4 1 % outflow +4 0.6 -1 0.61 0.2 % slip # diff --git a/examples/loadbalance/adaptation.hh b/examples/loadbalance/adaptation.hh index a5d4d702a75a37e9f606a8def2d9f7c9ddd91b90..625acf25726f35a63fc84db536e9056d1471b5d5 100644 --- a/examples/loadbalance/adaptation.hh +++ b/examples/loadbalance/adaptation.hh @@ -73,6 +73,7 @@ public: solution_( 0 ), adaptTimer_(), adaptTime_( 0.0 ), + restProlTime_( 0.0 ), lbTime_( 0.0 ), commTime_( 0.0 ) {} @@ -88,6 +89,8 @@ public: //! return time spent for the last adapation in sec double adaptationTime() const { return adaptTime_; } + //! return time spent for the last adapation in sec + double restProlTime() const { return restProlTime_; } //! return time spent for the last load balancing in sec double loadBalanceTime() const { return lbTime_; } //! return time spent for the last communication in sec @@ -135,6 +138,7 @@ private: Dune :: Timer adaptTimer_ ; double adaptTime_; + double restProlTime_; double lbTime_; double commTime_; }; @@ -151,6 +155,7 @@ inline void LeafAdaptation< Grid, Vector,LoadBalanceHandle >::operator() ( Vecto adaptTime_ = 0.0; lbTime_ = 0.0; commTime_ = 0.0; + restProlTime_ = 0.0; // reset timer adaptTimer_.reset() ; diff --git a/examples/loadbalance/main.cc b/examples/loadbalance/main.cc index 188382f5c30d2da9e4f52e2485e578ad3bb6ae04..b76d978633bbdae26984f9647e597c809d4aad80 100644 --- a/examples/loadbalance/main.cc +++ b/examples/loadbalance/main.cc @@ -257,7 +257,8 @@ void method ( int problem, int startLvl, int maxLvl, adaptation.adaptationTime(), // time for adaptation adaptation.loadBalanceTime(), // time for load balance overallTimer.elapsed(), // time step overall time - getMemoryUsage() ); // memory usage + adaptation.restProlTime(), // adaptation restrict prolong + getMemoryUsage() ); // memory usage } diff --git a/examples/main.cc b/examples/main.cc index b1b30d9a22dc7e0b8a0ba9eb317b1a49915ef2d9..cf00413d75435d461fbad6b627a894568155c1d1 100644 --- a/examples/main.cc +++ b/examples/main.cc @@ -247,8 +247,8 @@ void method ( int problem, int startLvl, int maxLvl, adaptation.adaptationTime(), // time for adaptation adaptation.loadBalanceTime(), // time for load balance overallTimer.elapsed(), // time step overall time - adaptation.restProlTime(), - getMemoryUsage() ); // memory usage + adaptation.restProlTime(), // adaptation restrict prolong + getMemoryUsage() ); // memory usage } diff --git a/examples/problem-euler.hh b/examples/problem-euler.hh index 1dd22514aa8ed73e4851bc33fd2452b66271d94b..ad164452083a7418bced11420f4b25e6ea7f56a5 100644 --- a/examples/problem-euler.hh +++ b/examples/problem-euler.hh @@ -64,10 +64,13 @@ public: int bndType( const DomainType &normal, const DomainType &x, const double time) const { if (normal[0]<-0.1 && x[0]<0.1) - return 1; + return 1; // inflow else if (normal[0]>0.1 && x[0]>1) - return 2; - return 3; + return 2; // outflow + else if (normal[0]>0.1 && x[0] >= 0.59 && x[0] <= 0.61 ) + return 4; // slip + + return 3; // reflection } //! \copydoc ProblemData::endTime @@ -400,6 +403,19 @@ struct EulerModel return *problem_; } + double pressure( const double gamma, const RangeType &u ) const + { + assert( u[0] >= 1e-10 ); + double v = u[1]*u[1] ; + for( int i=2; i<dimDomain+1; ++i ) + v += u[i]*u[i]; + + const double rhoe = u[dimDomain+1]-0.5*(v)/u[0]; + + assert( rhoe>1e-10 ); + return (gamma-1.0)*rhoe; + } + /** \copydoc TransportProblem::numericalFlux */ double numericalFlux ( const DomainType &normal, const double time, @@ -407,18 +423,10 @@ struct EulerModel const RangeType &uLeft, const RangeType &uRight, RangeType &flux ) const { - // calculate unit normal - DomainType unitNormal( normal ); - const double faceVol = normal.two_norm(); - unitNormal *= 1.0 / faceVol; + assert( std::abs( normal.two_norm() - 1.0 ) < 1e-12 ); // apply numerical flux - const double dt = numFlux_.numFlux( uLeft, uRight, unitNormal, flux ); - - // apply face volume - flux *= faceVol; - - return dt * faceVol; + return numFlux_.numFlux( uLeft, uRight, normal, flux ); } /** \brief boundary ids for different types of boundary typical for the @@ -433,8 +441,21 @@ struct EulerModel const RangeType& uLeft, RangeType &flux ) const { + int bndType = problem_->bndType(normal,xGlobal,time); + if( bndType == 4 ) // slip + { + flux = 0; + // slip boundary + const double p = pressure(1.4, uLeft); + for (int i=0;i<dimDomain; ++i) + { + flux[i+1] = normal[i] * p; + } + return 0.; + } + RangeType uRight; - boundaryValue(normal,time,xGlobal,uLeft,uRight); + boundaryValue(bndType, normal,time,xGlobal,uLeft,uRight); return numericalFlux(normal,time,xGlobal,uLeft,uRight,flux); } @@ -454,7 +475,8 @@ struct EulerModel const RangeType& uLeft) const { RangeType uRight; - boundaryValue(normal,time,xGlobal,uLeft,uRight); + int bndType = problem_->bndType(normal,xGlobal,time); + boundaryValue( bndType, normal,time,xGlobal,uLeft,uRight); return indicator( normal,time,xGlobal, uLeft, uRight ); } @@ -462,18 +484,18 @@ private: /** \brief the boundary flux inserts a ghost cell value into the * numerical flux - this function computes these values for the * different boundary types */ - void boundaryValue( const DomainType &normal, + void boundaryValue( const int bndType, + const DomainType &normal, const double time, const DomainType &xGlobal, const RangeType& uLeft, RangeType& uRight) const { - int bndType = problem_->bndType(normal,xGlobal,time); if (bndType == 1) // Dirichlet uRight = problem().boundaryValue(xGlobal,time); else if (bndType == 2) // Neumann uRight = uLeft; - else + else if (bndType == 3) { DomainType unitNormal( normal ); const double faceVol = normal.two_norm(); @@ -488,7 +510,7 @@ private: } } - /** \brief rotate the veclocity field into the face frame of reference */ + /** \brief rotate the velocity field into the face frame of reference */ void rotate( const DomainType &normal, const RangeType &u, RangeType &u_rot) const { diff --git a/examples/testEfficiency/main.cc b/examples/testEfficiency/main.cc index f5c739ffe07bbf8a782b3ae677f6dfa5a55c6610..cea7c290208ce565bc97cc8b73a2ecbea9fb515f 100644 --- a/examples/testEfficiency/main.cc +++ b/examples/testEfficiency/main.cc @@ -245,7 +245,8 @@ void method ( int problem, int startLvl, int maxLvl, adaptation.adaptationTime(), // time for adaptation adaptation.loadBalanceTime(), // time for load balance overallTimer.elapsed(), // time step overall time - getMemoryUsage() ); // memory usage + adaptation.restProlTime(), // adaptation restrict prolong + getMemoryUsage() ); // memory usage }