StencilRestrictedPackInfo

- generalized PdfFieldPackInfo to be used with a stencil template
  argument instead of a LatticeModel
- PdfPackInfo now available from Python

src/field/communication/StencilRestrictedPackInfo.h

+//======================================================================================================================
+//
+//  This file is part of waLBerla. waLBerla is free software: you can
+//  redistribute it and/or modify it under the terms of the GNU General Public
+//  License as published by the Free Software Foundation, either version 3 of
+//  the License, or (at your option) any later version.
+//
+//  waLBerla is distributed in the hope that it will be useful, but WITHOUT
+//  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+//  FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+//  for more details.
+//
+//  You should have received a copy of the GNU General Public License along
+//  with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>.
+//
+//! \file PdfFieldPackInfo.h
+//! \ingroup lbm
+//! \author Martin Bauer <martin.bauer@fau.de>
+//! \author Florian Schornbaum <florian.schornbaum@fau.de>
+//! \brief Packs only certain components of a field
+//
+//======================================================================================================================
+
+#pragma once
+
+#include "communication/UniformPackInfo.h"
+#include "core/debug/Debug.h"
+#include "stencil/Directions.h"
+
+
+namespace walberla {
+namespace field {
+namespace communication {
+
+
+/**
+ * \brief PackInfo that packs only the components that point to the neighbor.
+ *
+ * see also documentation for FieldPackInfo
+ *
+ * \warning For fields with nrOfGhostLayers > 1:  only the components pointing towards
+ *          the boundary are communicated, which may not be the desired behavior
+ *          for the 'inner' ghost layers
+ */
+template< typename GhostLayerField_T, typename Stencil_T >
+class StencilRestrictedPackInfo : public walberla::communication::UniformPackInfo
+{
+public:
+   StencilRestrictedPackInfo( const BlockDataID & fieldId ) : fieldId_( fieldId ) {}
+   virtual ~StencilRestrictedPackInfo() {}
+
+   bool constantDataExchange() const { return true; }
+   bool threadsafeReceiving()  const { return true; }
+
+   void unpackData( IBlock * receiver, stencil::Direction dir, mpi::RecvBuffer & buffer );
+
+   void communicateLocal( const IBlock * sender, IBlock * receiver, stencil::Direction dir );
+
+protected:
+
+   void packDataImpl( const IBlock * sender, stencil::Direction dir, mpi::SendBuffer & outBuffer ) const;
+
+   const BlockDataID fieldId_;
+
+   static_assert(GhostLayerField_T::F_SIZE == Stencil_T::Size, "Size of stencil and f size of field have to be equal");
+};
+
+
+template< typename GhostLayerField_T, typename Stencil >
+void StencilRestrictedPackInfo<GhostLayerField_T, Stencil>::unpackData( IBlock * receiver, stencil::Direction dir, mpi::RecvBuffer & buffer )
+{
+   if( Stencil::idx[ stencil::inverseDir[dir] ] >= Stencil::Size )
+      return;
+
+   GhostLayerField_T * pdfField = receiver->getData< GhostLayerField_T >( fieldId_ );
+   WALBERLA_ASSERT_NOT_NULLPTR( pdfField );
+   WALBERLA_ASSERT_EQUAL( pdfField->nrOfGhostLayers(), 1 );
+
+   stencil::Direction packerDirection = stencil::inverseDir[dir];
+
+   for(auto i = pdfField->beginGhostLayerOnlyXYZ(dir); i != pdfField->end(); ++i )
+      for(uint_t f = 0; f < Stencil::d_per_d_length[packerDirection]; ++f)
+         buffer >> i.getF( Stencil::idx[ Stencil::d_per_d[packerDirection][f] ] );
+}
+
+
+
+template< typename GhostLayerField_T, typename Stencil >
+void StencilRestrictedPackInfo<GhostLayerField_T, Stencil>::communicateLocal( const IBlock * sender, IBlock * receiver, stencil::Direction dir )
+{
+   if( Stencil::idx[dir] >= Stencil::Size )
+      return;
+
+   const GhostLayerField_T * sf = sender  ->getData< GhostLayerField_T >( fieldId_ );
+         GhostLayerField_T * rf = receiver->getData< GhostLayerField_T >( fieldId_ );
+
+   WALBERLA_ASSERT_EQUAL( sf->xyzSize(), rf->xyzSize() );
+
+   typename GhostLayerField_T::const_iterator srcIter = sf->beginSliceBeforeGhostLayerXYZ(dir);
+   typename GhostLayerField_T::iterator       dstIter = rf->beginGhostLayerOnlyXYZ(stencil::inverseDir[dir]);
+
+   while( srcIter != sf->end() )
+   {
+      for( uint_t f = 0; f < Stencil::d_per_d_length[dir]; ++f )
+         dstIter.getF( Stencil::idx[ Stencil::d_per_d[dir][f] ] ) = srcIter.getF( Stencil::idx[ Stencil::d_per_d[dir][f] ] );
+
+      ++srcIter;
+      ++dstIter;
+   }
+   WALBERLA_ASSERT( srcIter == sf->end() );
+   WALBERLA_ASSERT( dstIter == rf->end() );
+}
+
+
+
+template< typename GhostLayerField_T, typename Stencil >
+void StencilRestrictedPackInfo<GhostLayerField_T, Stencil>::packDataImpl( const IBlock * sender, stencil::Direction dir, mpi::SendBuffer & outBuffer ) const
+{
+   if( Stencil::idx[dir] >= Stencil::Size )
+      return;
+
+   const GhostLayerField_T * pdfField = sender->getData< GhostLayerField_T >( fieldId_ );
+   WALBERLA_ASSERT_NOT_NULLPTR( pdfField );
+   WALBERLA_ASSERT_EQUAL( pdfField->nrOfGhostLayers(), 1 );
+
+   for( auto i = pdfField->beginSliceBeforeGhostLayerXYZ(dir); i != pdfField->end(); ++i )
+      for(uint_t f = 0; f < Stencil::d_per_d_length[dir]; ++f)
+         outBuffer << i.getF( Stencil::idx[ Stencil::d_per_d[dir][f] ] );
+}
+
+
+
+} // namespace communication
+} // namespace field
+} // namespace walberla

src/field/python/CommunicationExport.impl.h

@@ -26,12 +26,19 @@
 
 
 #include "field/communication/PackInfo.h"
+#include "field/communication/StencilRestrictedPackInfo.h"
 #include "field/communication/UniformMPIDatatypeInfo.h"
 
 #include "python_coupling/helper/MplHelpers.h"
 #include "python_coupling/helper/BoostPythonHelpers.h"
 #include "python_coupling/helper/MplHelpers.h"
 
+#include "stencil/D2Q9.h"
+#include "stencil/D3Q7.h"
+#include "stencil/D3Q15.h"
+#include "stencil/D3Q19.h"
+#include "stencil/D3Q27.h"
+
 
 namespace walberla {
 namespace field {
@@ -39,6 +46,86 @@ namespace field {
 
 namespace internal {
 
+   //===================================================================================================================
+   //
+   //  createStencilRestrictedPackInfo Export
+   //
+   //===================================================================================================================
+
+   template< typename FieldType >
+   typename std::enable_if<FieldType::F_SIZE == 27, boost::python::object>::type
+   createStencilRestrictedPackInfoObject( BlockDataID bdId )
+   {
+      typedef GhostLayerField<typename FieldType::value_type, 27> GlField_T;
+      using field::communication::StencilRestrictedPackInfo;
+      return boost::python::object( make_shared< StencilRestrictedPackInfo<GlField_T, stencil::D3Q27> >( bdId) );
+   }
+
+   template< typename FieldType >
+   typename std::enable_if<FieldType::F_SIZE == 19, boost::python::object>::type
+   createStencilRestrictedPackInfoObject( BlockDataID bdId )
+   {
+      typedef GhostLayerField<typename FieldType::value_type, 19> GlField_T;
+      using field::communication::StencilRestrictedPackInfo;
+      return boost::python::object( make_shared< StencilRestrictedPackInfo<GlField_T, stencil::D3Q19> >( bdId) );
+   }
+
+   template< typename FieldType >
+   typename std::enable_if<FieldType::F_SIZE == 15, boost::python::object>::type
+   createStencilRestrictedPackInfoObject( BlockDataID bdId )
+   {
+      typedef GhostLayerField<typename FieldType::value_type, 15> GlField_T;
+      using field::communication::StencilRestrictedPackInfo;
+      return boost::python::object( make_shared< StencilRestrictedPackInfo<GlField_T, stencil::D3Q15> >( bdId) );
+   }
+
+   template< typename FieldType >
+   typename std::enable_if<FieldType::F_SIZE == 7, boost::python::object>::type
+   createStencilRestrictedPackInfoObject( BlockDataID bdId )
+   {
+      typedef GhostLayerField<typename FieldType::value_type, 7> GlField_T;
+      using field::communication::StencilRestrictedPackInfo;
+      return boost::python::object( make_shared< StencilRestrictedPackInfo<GlField_T, stencil::D3Q7> >( bdId) );
+   }
+
+   template< typename FieldType >
+   typename std::enable_if<FieldType::F_SIZE == 9, boost::python::object>::type
+   createStencilRestrictedPackInfoObject( BlockDataID bdId )
+   {
+      typedef GhostLayerField<typename FieldType::value_type, 9> GlField_T;
+      using field::communication::StencilRestrictedPackInfo;
+      return boost::python::object( make_shared< StencilRestrictedPackInfo<GlField_T, stencil::D2Q9> >( bdId) );
+   }
+
+   template< typename FieldType >
+   typename std::enable_if<!(FieldType::F_SIZE == 9  ||
+                             FieldType::F_SIZE == 7  ||
+                             FieldType::F_SIZE == 15 ||
+                             FieldType::F_SIZE == 19 ||
+                             FieldType::F_SIZE == 27), boost::python::object>::type
+   createStencilRestrictedPackInfoObject( BlockDataID )
+   {
+      PyErr_SetString( PyExc_ValueError, "This works only for fields with fSize in 7, 9, 15, 19 or 27" );
+      throw boost::python::error_already_set();
+   }
+
+   FunctionExporterClass( createStencilRestrictedPackInfoObject, boost::python::object( BlockDataID ) );
+
+   template< typename FieldVector>
+   boost::python::object createStencilRestrictedPackInfo( const shared_ptr<StructuredBlockStorage> & bs,
+                                                          const std::string & blockDataName )
+   {
+      auto bdId = python_coupling::blockDataIDFromString( *bs, blockDataName );
+      if ( bs->begin() == bs->end() ) {
+         // if no blocks are on this field an arbitrary PackInfo can be returned
+         return createStencilRestrictedPackInfoObject< GhostLayerField<real_t,1> > ( bdId );
+      }
+
+      IBlock * firstBlock =  & ( * bs->begin() );
+      python_coupling::Dispatcher<FieldVector, Exporter_createStencilRestrictedPackInfoObject > dispatcher( firstBlock );
+      return dispatcher( bdId )( bdId ) ;
+   }
+
    //===================================================================================================================
    //
    //  createPackInfo Export
@@ -114,6 +201,34 @@ namespace internal {
       return dispatcher( bdId )( bdId, numberOfGhostLayers );
    }
 
+   template< typename T>
+   void exportStencilRestrictedPackInfo()
+   {
+      using field::communication::StencilRestrictedPackInfo;
+      using namespace boost::python;
+
+      {
+         typedef StencilRestrictedPackInfo<GhostLayerField<T, 9>, stencil::D2Q9> Pi;
+         class_< Pi, shared_ptr<Pi>, bases<walberla::communication::UniformPackInfo>,  boost::noncopyable >( "StencilRestrictedPackInfo", no_init );
+      }
+      {
+         typedef StencilRestrictedPackInfo<GhostLayerField<T, 7>, stencil::D3Q7> Pi;
+         class_< Pi, shared_ptr<Pi>, bases<walberla::communication::UniformPackInfo>,  boost::noncopyable >( "StencilRestrictedPackInfo", no_init );
+      }
+      {
+         typedef StencilRestrictedPackInfo<GhostLayerField<T, 15>, stencil::D3Q15> Pi;
+         class_< Pi, shared_ptr<Pi>, bases<walberla::communication::UniformPackInfo>,  boost::noncopyable >( "StencilRestrictedPackInfo", no_init );
+      }
+      {
+         typedef StencilRestrictedPackInfo<GhostLayerField<T, 19>, stencil::D3Q19> Pi;
+         class_< Pi, shared_ptr<Pi>, bases<walberla::communication::UniformPackInfo>,  boost::noncopyable >( "StencilRestrictedPackInfo", no_init );
+      }
+      {
+         typedef StencilRestrictedPackInfo<GhostLayerField<T, 27>, stencil::D3Q27> Pi;
+         class_< Pi, shared_ptr<Pi>, bases<walberla::communication::UniformPackInfo>,  boost::noncopyable >( "StencilRestrictedPackInfo", no_init );
+      }
+
+   }
 
 } // namespace internal
 
@@ -126,8 +241,13 @@ void exportCommunicationClasses()
 {
    using namespace boost::python;
 
+   internal::exportStencilRestrictedPackInfo<float>();
+   internal::exportStencilRestrictedPackInfo<double>();
+
    def( "createMPIDatatypeInfo",&internal::createMPIDatatypeInfo<FieldTypes>, ( arg("blocks"), arg("blockDataName"), arg("numberOfGhostLayers" ) =0 ) );
    def( "createPackInfo",       &internal::createPackInfo<FieldTypes>,        ( arg("blocks"), arg("blockDataName"), arg("numberOfGhostLayers" ) =0 ) );
+   def( "createStencilRestrictedPackInfo", &internal::createStencilRestrictedPackInfo<FieldTypes>,
+        (arg("blocks"), arg("blockDataName") ));
 }
 
 

src/lbm/communication/PdfFieldPackInfo.h

@@ -24,6 +24,7 @@
 #pragma once
 
 #include "lbm/field/PdfField.h"
+#include "field/communication/StencilRestrictedPackInfo.h"
 #include "communication/UniformPackInfo.h"
 #include "core/debug/Debug.h"
 #include "stencil/Directions.h"
@@ -55,96 +56,8 @@ namespace lbm {
  * \ingroup lbm
  */
 template< typename LatticeModel_T >
-class PdfFieldPackInfo : public walberla::communication::UniformPackInfo
-{
-public:
-
-   typedef PdfField<LatticeModel_T>          PdfField_T;
-   typedef typename LatticeModel_T::Stencil  Stencil;
-
-   PdfFieldPackInfo( const BlockDataID & pdfFieldId ) : pdfFieldId_( pdfFieldId ) {}
-   virtual ~PdfFieldPackInfo() {}
-
-   bool constantDataExchange() const { return true; }
-   bool threadsafeReceiving()  const { return true; }
-
-   void unpackData( IBlock * receiver, stencil::Direction dir, mpi::RecvBuffer & buffer );
-
-   void communicateLocal( const IBlock * sender, IBlock * receiver, stencil::Direction dir );
-
-protected:
-
-   void packDataImpl( const IBlock * sender, stencil::Direction dir, mpi::SendBuffer & outBuffer ) const;
-
-
-
-   const BlockDataID pdfFieldId_;
-};
-
-
-
-template< typename LatticeModel_T >
-void PdfFieldPackInfo< LatticeModel_T >::unpackData( IBlock * receiver, stencil::Direction dir, mpi::RecvBuffer & buffer )
-{
-   if( Stencil::idx[ stencil::inverseDir[dir] ] >= Stencil::Size )
-      return;
-
-   PdfField_T * pdfField = receiver->getData< PdfField_T >( pdfFieldId_ );
-   WALBERLA_ASSERT_NOT_NULLPTR( pdfField );
-   WALBERLA_ASSERT_EQUAL( pdfField->nrOfGhostLayers(), 1 );
-
-   stencil::Direction packerDirection = stencil::inverseDir[dir];
-
-   for(auto i = pdfField->beginGhostLayerOnlyXYZ(dir); i != pdfField->end(); ++i )
-      for(uint_t f = 0; f < Stencil::d_per_d_length[packerDirection]; ++f)
-         buffer >> i.getF( Stencil::idx[ Stencil::d_per_d[packerDirection][f] ] );
-}
-
-
-
-template< typename LatticeModel_T >
-void PdfFieldPackInfo< LatticeModel_T >::communicateLocal( const IBlock * sender, IBlock * receiver, stencil::Direction dir )
-{
-   if( Stencil::idx[dir] >= Stencil::Size )
-      return;
-
-   const PdfField_T * sf = sender  ->getData< PdfField_T >( pdfFieldId_ );
-         PdfField_T * rf = receiver->getData< PdfField_T >( pdfFieldId_ );
-
-   WALBERLA_ASSERT_EQUAL( sf->xyzSize(), rf->xyzSize() );
-
-   typename PdfField_T::const_iterator srcIter = sf->beginSliceBeforeGhostLayerXYZ(dir);
-   typename PdfField_T::iterator       dstIter = rf->beginGhostLayerOnlyXYZ(stencil::inverseDir[dir]);
-
-   while( srcIter != sf->end() )
-   {
-      for( uint_t f = 0; f < Stencil::d_per_d_length[dir]; ++f )
-         dstIter.getF( Stencil::idx[ Stencil::d_per_d[dir][f] ] ) = srcIter.getF( Stencil::idx[ Stencil::d_per_d[dir][f] ] );
-
-      ++srcIter;
-      ++dstIter;
-   }
-   WALBERLA_ASSERT( srcIter == sf->end() );
-   WALBERLA_ASSERT( dstIter == rf->end() );
-}
-
-
-
-template< typename LatticeModel_T >
-void PdfFieldPackInfo< LatticeModel_T >::packDataImpl( const IBlock * sender, stencil::Direction dir, mpi::SendBuffer & outBuffer ) const
-{
-   if( Stencil::idx[dir] >= Stencil::Size )
-      return;
-
-   const PdfField_T * pdfField = sender->getData< PdfField_T >( pdfFieldId_ );
-   WALBERLA_ASSERT_NOT_NULLPTR( pdfField );
-   WALBERLA_ASSERT_EQUAL( pdfField->nrOfGhostLayers(), 1 );
-
-   for( auto i = pdfField->beginSliceBeforeGhostLayerXYZ(dir); i != pdfField->end(); ++i )
-      for(uint_t f = 0; f < Stencil::d_per_d_length[dir]; ++f)
-         outBuffer << i.getF( Stencil::idx[ Stencil::d_per_d[dir][f] ] );
-}
-
+using PdfFieldPackInfo = field::communication::StencilRestrictedPackInfo< PdfField<LatticeModel_T>,
+                                                                          typename LatticeModel_T::Stencil>;
 
 
 } // namespace lbm