extended pe solver interface

getMaximumPenetration()
getNumberOfContacts()
getNumberOfContactsTreated()

src/pe/cr/DEM.cpp

@@ -54,13 +54,19 @@ DEM::DEM(    const shared_ptr<BodyStorage>&    globalBodyStorage
    , ccdID_(ccdID)
    , fcdID_(fcdID)
    , tt_(tt)
+   , maxPenetration_(0)
+   , numberOfContacts_(0)
+   , numberOfContactsTreated_(0)
 {
 
 }
 
 void DEM::timestep( real_t dt )
 {
-//   auto numContacts = 0;
+   maxPenetration_          = real_c(0.0);
+   numberOfContacts_        = 0;
+   numberOfContactsTreated_ = 0;
+
    for (auto it = blockStorage_->begin(); it != blockStorage_->end(); ++it){
       IBlock & currentBlock = *it;
 
@@ -71,16 +77,18 @@ void DEM::timestep( real_t dt )
       Contacts& cont = fcd->generateContacts( ccd->getPossibleContacts() );
       if (tt_ != NULL) tt_->stop("FCD");
 
-      real_t maxOverlap = real_c(0);
       for (auto cIt = cont.begin(); cIt != cont.end(); ++cIt){
          const real_t overlap( -cIt->getDistance() );
-         if( overlap > maxOverlap )
-            maxOverlap = overlap;
+         if( overlap > maxPenetration_ )
+            maxPenetration_ = overlap;
          if (shouldContactBeTreated( &(*cIt), currentBlock.getAABB() ))
+         {
+            ++numberOfContactsTreated_;
             resolveContact( &(*cIt) );
+         }
       }
 
-//      numContacts += cont.size();
+      numberOfContacts_ += cont.size();
 
       cont.clear();
    }

src/pe/cr/DEM.h

@@ -52,6 +52,10 @@ public:
    void operator()(const real_t dt) { timestep(dt); }
    /// Advances the simulation dt seconds.
    void timestep( const real_t dt );
+
+   virtual inline real_t            getMaximumPenetration()        const WALBERLA_OVERRIDE { return maxPenetration_; }
+   virtual inline size_t            getNumberOfContacts()          const WALBERLA_OVERRIDE { return numberOfContacts_; }
+   virtual inline size_t            getNumberOfContactsTreated()   const WALBERLA_OVERRIDE { return numberOfContactsTreated_; }
 private:
    void resolveContact( ContactID c ) const;
    void move( BodyID id, real_t dt );
@@ -62,6 +66,10 @@ private:
    domain_decomposition::BlockDataID ccdID_;
    domain_decomposition::BlockDataID fcdID_;
    WcTimingTree*                     tt_;
+
+   real_t                            maxPenetration_;
+   size_t                            numberOfContacts_;
+   size_t                            numberOfContactsTreated_;
 };
 
 }  // namespace cr

src/pe/cr/HCSITS.h

@@ -129,8 +129,9 @@ public:
    //**Get functions*******************************************************************************
    /*!\name Get functions */
    //@{
-   inline real_t            getMaximumPenetration() const;
-   inline size_t            getNumberOfContacts()   const;
+   virtual inline real_t            getMaximumPenetration()        const WALBERLA_OVERRIDE;
+   virtual inline size_t            getNumberOfContacts()          const WALBERLA_OVERRIDE;
+   virtual inline size_t            getNumberOfContactsTreated()   const WALBERLA_OVERRIDE;
    inline const std::map<IBlockID::IDType, ContactCache> getContactCache() const { return blockToContactCache_; }
    //@}
    //**********************************************************************************************
@@ -231,6 +232,7 @@ private:
    real_t relaxationParam_;           //!< Parameter specifying underrelaxation of velocity corrections for boundary bodies.
    real_t maximumPenetration_;
    size_t numContacts_;
+   size_t numContactsTreated_;
 
    //**********************************************************************************************
    /*! \cond WALBERLA_INTERNAL */
@@ -294,6 +296,11 @@ inline size_t HardContactSemiImplicitTimesteppingSolvers::getNumberOfContacts()
 {
    return numContacts_;
 }
+
+inline size_t HardContactSemiImplicitTimesteppingSolvers::getNumberOfContactsTreated() const
+{
+   return numContactsTreated_;
+}
 //*************************************************************************************************
 
 

src/pe/cr/HCSITS.impl.h

@@ -102,6 +102,7 @@ inline HardContactSemiImplicitTimesteppingSolvers::HardContactSemiImplicitTimest
    , relaxationParam_  ( real_c(0.9) )
    , maximumPenetration_ ( real_c(0.0) )
    , numContacts_      ( 0 )
+   , numContactsTreated_( 0)
    , requireSync_      ( false )
 {
    // Logging the successful setup of the collision system
@@ -151,6 +152,9 @@ inline void HardContactSemiImplicitTimesteppingSolvers::timestep( const real_t d
 
    const real_t dtinv( real_c(1) / dt );
 
+   numContacts_        = 0;
+   numContactsTreated_ = 0;
+
    if (tt_ != NULL) tt_->start("Simulation Step");
 
    for (auto it = blockStorage_->begin(); it != blockStorage_->end(); ++it)
@@ -194,6 +198,9 @@ inline void HardContactSemiImplicitTimesteppingSolvers::timestep( const real_t d
             contactsMask_[i] = false;
          }
       }
+      numContacts_        += numContacts;
+      numContactsTreated_ += numContactsMasked;
+
       //      WALBERLA_LOG_DEVEL("contact filtering: " << numContactsMasked << "/" << contacts.size());
       if (tt_ != NULL) tt_->stop("Filtering");
       if (tt_ != NULL) tt_->stop("Collision Detection");
@@ -204,7 +211,6 @@ inline void HardContactSemiImplicitTimesteppingSolvers::timestep( const real_t d
 
       {
          maximumPenetration_ = 0;
-         numContacts_ = numContactsMasked;
 
          size_t j = 0;
          for( size_t i = 0; i < numContacts; ++i )

src/pe/cr/ICR.h

@@ -50,11 +50,46 @@ public:
    /// This can be used for example to set a gravitational force.
    inline void        setGlobalLinearAcceleration(const Vec3& acc) { globalLinearAcceleration_ = acc; }
    inline const Vec3& getGlobalLinearAcceleration() { return globalLinearAcceleration_; }
+
+   virtual inline real_t            getMaximumPenetration()        const;
+   virtual inline size_t            getNumberOfContacts()          const;
+   virtual inline size_t            getNumberOfContactsTreated()   const;
 private:
 //   real_t globalLinearDrag_;
    Vec3   globalLinearAcceleration_;
 };
 
+inline real_t ICR::getMaximumPenetration()        const
+{
+   static bool warned = false;
+   if (!warned) {
+      WALBERLA_LOG_WARNING("getMaximumPenetration() is not implemented for this solver!");
+      warned = true;
+   }
+
+   return real_c(0);
+}
+inline size_t ICR::getNumberOfContacts()          const
+{
+   static bool warned = false;
+   if (!warned) {
+      WALBERLA_LOG_WARNING("getMaximumPenetration() is not implemented for this solver!");
+      warned = true;
+   }
+
+   return 0;
+}
+inline size_t ICR::getNumberOfContactsTreated()   const
+{
+   static bool warned = false;
+   if (!warned) {
+      WALBERLA_LOG_WARNING("getMaximumPenetration() is not implemented for this solver!");
+      warned = true;
+   }
+
+   return 0;
+}
+
 }  // namespace cr
 }  // namespace pe
 }  // namespace walberla