Merge branch 'mr_reduce_packinfo_generation' into 'master'

Reduce Packinfo Generation

See merge request !464

python/pystencils_walberla/codegen.py

@@ -146,6 +146,8 @@ def generate_selective_sweep(generation_context, class_name, selection_tree, int
 
 def generate_pack_info_for_field(generation_context, class_name: str, field: Field,
                                  direction_subset: Optional[Tuple[Tuple[int, int, int]]] = None,
+                                 operator=None,
+                                 gl_to_inner=False,
                                  **create_kernel_params):
     """Creates a pack info for a pystencils field assuming a pull-type stencil, packing all cell elements.
 
@@ -155,18 +157,21 @@ def generate_pack_info_for_field(generation_context, class_name: str, field: Fie
         field: pystencils field for which to generate pack info
         direction_subset: optional sequence of directions for which values should be packed
                           otherwise a D3Q27 stencil is assumed
+        operator: optional operator for, e.g., reduction pack infos
+        gl_to_inner: communicates values from ghost layers of sender to interior of receiver
         **create_kernel_params: remaining keyword arguments are passed to `pystencils.create_kernel`
     """
+
     if not direction_subset:
         direction_subset = tuple((i, j, k) for i, j, k in product(*[(-1, 0, 1)] * 3))
 
     all_index_accesses = [field(*ind) for ind in product(*[range(s) for s in field.index_shape])]
-    return generate_pack_info(generation_context, class_name, {direction_subset: all_index_accesses},
-                              **create_kernel_params)
+    return generate_pack_info(generation_context, class_name, {direction_subset: all_index_accesses}, operator=operator,
+                              gl_to_inner=gl_to_inner, **create_kernel_params)
 
 
 def generate_pack_info_from_kernel(generation_context, class_name: str, assignments: Sequence[Assignment],
-                                   kind='pull', **create_kernel_params):
+                                   kind='pull', operator=None, **create_kernel_params):
     """Generates a waLBerla GPU PackInfo from a (pull) kernel.
 
     Args:
@@ -175,6 +180,7 @@ def generate_pack_info_from_kernel(generation_context, class_name: str, assignme
         assignments: list of assignments from the compute kernel - generates PackInfo for "pull" part only
                      i.e. the kernel is expected to only write to the center
         kind: can either be pull or push
+        operator: optional operator for, e.g., reduction pack infos
         **create_kernel_params: remaining keyword arguments are passed to `pystencils.create_kernel`
     """
     assert kind in ('push', 'pull')
@@ -207,12 +213,12 @@ def generate_pack_info_from_kernel(generation_context, class_name: str, assignme
                 spec[(comm_dir,)].add(fa)
     else:
         raise ValueError("Invalid 'kind' parameter")
-    return generate_pack_info(generation_context, class_name, spec, **create_kernel_params)
+    return generate_pack_info(generation_context, class_name, spec, operator=operator, **create_kernel_params)
 
 
 def generate_pack_info(generation_context, class_name: str,
                        directions_to_pack_terms: Dict[Tuple[Tuple], Sequence[Field.Access]],
-                       namespace='pystencils',
+                       namespace='pystencils', operator=None, gl_to_inner=False,
                        **create_kernel_params):
     """Generates a waLBerla GPU PackInfo
 
@@ -222,6 +228,8 @@ def generate_pack_info(generation_context, class_name: str,
         directions_to_pack_terms: maps tuples of directions to read field accesses, specifying which values have to be
                                   packed for which direction
         namespace: inner namespace of the generated class
+        operator: optional operator for, e.g., reduction pack infos
+        gl_to_inner: communicates values from ghost layers of sender to interior of receiver
         **create_kernel_params: remaining keyword arguments are passed to `pystencils.create_kernel`
     """
     items = [(e[0], sorted(e[1], key=lambda x: str(x))) for e in directions_to_pack_terms.items()]
@@ -274,7 +282,10 @@ def generate_pack_info(generation_context, class_name: str,
         pack_ast = create_kernel(pack_assignments, **create_kernel_params, ghost_layers=0)
         pack_ast.function_name = 'pack_{}'.format("_".join(direction_strings))
         pack_ast.assumed_inner_stride_one = create_kernel_params['cpu_vectorize_info']['assume_inner_stride_one']
-        unpack_assignments = [Assignment(term, buffer(i)) for i, term in enumerate(terms)]
+        if operator is None:
+            unpack_assignments = [Assignment(term, buffer(i)) for i, term in enumerate(terms)]
+        else:
+            unpack_assignments = [Assignment(term, operator(term, buffer(i))) for i, term in enumerate(terms)]
         unpack_ast = create_kernel(unpack_assignments, **create_kernel_params, ghost_layers=0)
         unpack_ast.function_name = 'unpack_{}'.format("_".join(direction_strings))
         unpack_ast.assumed_inner_stride_one = create_kernel_params['cpu_vectorize_info']['assume_inner_stride_one']
@@ -296,6 +307,7 @@ def generate_pack_info(generation_context, class_name: str,
         'dtype': dtype,
         'field_name': field_names.pop(),
         'namespace': namespace,
+        'gl_to_inner': gl_to_inner,
     }
     env = Environment(loader=PackageLoader('pystencils_walberla'), undefined=StrictUndefined)
     add_pystencils_filters_to_jinja_env(env)

python/pystencils_walberla/templates/CpuPackInfo.tmpl.cpp

@@ -37,7 +37,11 @@ void {{class_name}}::pack(Direction dir, unsigned char * byte_buffer, IBlock * b
 
     {{fused_kernel|generate_block_data_to_field_extraction(parameters_to_ignore=['buffer'])|indent(4)}}
     CellInterval ci;
+    {% if gl_to_inner -%}
+    {{field_name}}->getGhostRegion(dir, ci, 1, false);
+    {%- else -%}
     {{field_name}}->getSliceBeforeGhostLayer(dir, ci, 1, false);
+    {%- endif %}
 
     switch( dir )
     {
@@ -63,7 +67,11 @@ void {{class_name}}::unpack(Direction dir, unsigned char * byte_buffer, IBlock *
 
     {{fused_kernel|generate_block_data_to_field_extraction(parameters_to_ignore=['buffer'])|indent(4)}}
     CellInterval ci;
+    {% if gl_to_inner -%}
+    {{field_name}}->getSliceBeforeGhostLayer(dir, ci, 1, false);
+    {%- else -%}
     {{field_name}}->getGhostRegion(dir, ci, 1, false);
+    {%- endif %}
     auto communciationDirection = stencil::inverseDir[dir];
 
     switch( communciationDirection )

python/pystencils_walberla/templates/GpuPackInfo.tmpl.cpp

@@ -35,7 +35,11 @@ void {{class_name}}::pack(Direction dir, unsigned char * byte_buffer, IBlock * b
 
     {{fused_kernel|generate_block_data_to_field_extraction(parameters_to_ignore=['buffer'])|indent(4)}}
     CellInterval ci;
+    {% if gl_to_inner -%}
+    {{field_name}}->getGhostRegion(dir, ci, 1, false);
+    {%- else -%}
     {{field_name}}->getSliceBeforeGhostLayer(dir, ci, 1, false);
+    {%- endif %}
 
     switch( dir )
     {
@@ -61,7 +65,11 @@ void {{class_name}}::unpack(Direction dir, unsigned char * byte_buffer, IBlock *
 
     {{fused_kernel|generate_block_data_to_field_extraction(parameters_to_ignore=['buffer'])|indent(4)}}
     CellInterval ci;
+    {% if gl_to_inner -%}
+    {{field_name}}->getSliceBeforeGhostLayer(dir, ci, 1, false);
+    {%- else -%}
     {{field_name}}->getGhostRegion(dir, ci, 1, false);
+    {%- endif %}
     auto communciationDirection = stencil::inverseDir[dir];
 
     switch( communciationDirection )

src/cuda/communication/UniformGPUScheme.impl.h

@@ -118,6 +118,7 @@ UniformGPUScheme<Stencil>::UniformGPUScheme( weak_ptr <StructuredBlockForest> bf
          auto parallelSection = parallelSectionManager_.parallelSection( stream );
          for( auto recvInfo = bufferSystemGPU_.begin(); recvInfo != bufferSystemGPU_.end(); ++recvInfo )
          {
+            recvInfo.buffer().clear();
             for( auto &header : headers_[recvInfo.rank()] )
             {
                auto block = dynamic_cast< Block * >( forest->getBlock( header.blockId ));
@@ -141,6 +142,7 @@ UniformGPUScheme<Stencil>::UniformGPUScheme( weak_ptr <StructuredBlockForest> bf
          {
             auto &gpuBuffer = bufferSystemGPU_.sendBuffer( recvInfo.rank());
 
+            recvInfo.buffer().clear();
             gpuBuffer.clear();
             for( auto &header : headers_[recvInfo.rank()] ) {
                auto block = dynamic_cast< Block * >( forest->getBlock( header.blockId ));

tests/cuda/CMakeLists.txt

@@ -47,4 +47,10 @@ waLBerla_generate_target_from_python(NAME MicroBenchmarkGpuLbmGenerated FILE cod
       OUT_FILES MicroBenchmarkStreamKernel.cu MicroBenchmarkCopyKernel.cu MicroBenchmarkStreamKernel.h MicroBenchmarkCopyKernel.h)
 waLBerla_compile_test( FILES codegen/MicroBenchmarkGpuLbm.cpp DEPENDS MicroBenchmarkGpuLbmGenerated)
 
+waLBerla_generate_target_from_python(NAME CodegenGeneratedGPUFieldPackInfo FILE codegen/GeneratedFieldPackInfoTestGPU.py
+        OUT_FILES ScalarFieldCommunicationGPU.cu ScalarFieldCommunicationGPU.h
+        ScalarFieldPullReductionGPU.cu ScalarFieldPullReductionGPU.h )
+waLBerla_compile_test( FILES codegen/GeneratedFieldPackInfoTestGPU.cpp
+        DEPENDS blockforest core field CodegenGeneratedGPUFieldPackInfo )
+waLBerla_execute_test( NAME GeneratedFieldPackInfoTestGPU )
 endif()
\ No newline at end of file

tests/cuda/codegen/GeneratedFieldPackInfoTestGPU.cpp

+//======================================================================================================================
+//
+//  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 GeneratedFieldPackInfoTestGPU.cpp
+//! \ingroup field
+//! \author Helen Schottenhamml <helen.schottenhamml@fau.de>
+//! \brief Tests if a GPU Field is correctly communicated using generated pack info
+//
+//======================================================================================================================
+
+#include "field/AddToStorage.h"
+#include "field/GhostLayerField.h"
+
+#include "blockforest/Initialization.h"
+
+#include "core/debug/TestSubsystem.h"
+#include "core/Environment.h"
+
+#include "cuda/FieldCopy.h"
+#include "cuda/communication/UniformGPUScheme.h"
+
+#include "stencil/D3Q27.h"
+
+// include generated files
+#include "ScalarFieldCommunicationGPU.h"
+#include "ScalarFieldPullReductionGPU.h"
+
+
+namespace walberla {
+
+using Stencil_T = stencil::D3Q27;
+
+cuda::GPUField<int> * createGPUField( IBlock* const block, StructuredBlockStorage* const storage ) {
+
+   return new cuda::GPUField<int> (
+      storage->getNumberOfXCells( *block ), // number of cells in x direction
+      storage->getNumberOfYCells( *block ), // number of cells in y direction
+      storage->getNumberOfZCells( *block ), // number of cells in z direction
+      1,                                    // fSize
+      1,                                    // number of ghost layers
+      field::fzyx );
+
+}
+
+cuda::GPUField<int> * createSmallGPUField( IBlock * const , StructuredBlockStorage * const ) {
+   return new cuda::GPUField<int> (2, 2, 2, 1, 1, field::fzyx );
+}
+
+
+void testScalarField( std::shared_ptr<blockforest::StructuredBlockForest> & sbf, BlockDataID gpuFieldId ) {
+
+   cuda::communication::UniformGPUScheme< Stencil_T > us{ sbf };
+   us.addPackInfo(std::make_shared< pystencils::ScalarFieldCommunicationGPU >(gpuFieldId));
+
+   for( auto & block : *sbf ) {
+
+      auto & gpuField = *(block.getData< cuda::GPUField< int > >(gpuFieldId));
+
+      field::GhostLayerField< int, 1 > cpuField(gpuField.xSize(), gpuField.ySize(), gpuField.zSize(), 1, 0,
+                                                field::fzyx);
+      cpuField.setWithGhostLayer(0);
+
+      WALBERLA_CHECK_EQUAL(cpuField.xSize(), 2)
+      WALBERLA_CHECK_EQUAL(cpuField.ySize(), 2)
+      WALBERLA_CHECK_EQUAL(cpuField.zSize(), 2)
+
+      // initialize the bottom boundary
+      cpuField(0, 0, 0) = 1;
+      cpuField(0, 1, 0) = 2;
+      cpuField(1, 0, 0) = 3;
+      cpuField(1, 1, 0) = 4;
+
+      cuda::fieldCpy(gpuField, cpuField);
+
+      // communicate
+      us.communicate();
+
+      cuda::fieldCpy(cpuField, gpuField);
+
+      WALBERLA_CHECK_EQUAL(cpuField(0, 0, +2), 1)
+      WALBERLA_CHECK_EQUAL(cpuField(0, 1, +2), 2)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 0, +2), 3)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 1, +2), 4)
+
+   }
+}
+
+void testScalarFieldPullReduction( std::shared_ptr<blockforest::StructuredBlockForest> & sbf, BlockDataID gpuFieldId ) {
+
+   cuda::communication::UniformGPUScheme< Stencil_T > us1{ sbf };
+   us1.addPackInfo(std::make_shared< pystencils::ScalarFieldPullReductionGPU >(gpuFieldId));
+
+   cuda::communication::UniformGPUScheme< Stencil_T > us2{ sbf };
+   us2.addPackInfo(std::make_shared< pystencils::ScalarFieldCommunicationGPU >(gpuFieldId));
+
+   for( auto & block : *sbf ) {
+
+      auto& gpuField = *(block.getData< cuda::GPUField< int > >(gpuFieldId));
+
+      field::GhostLayerField< int, 1 > cpuField(gpuField.xSize(), gpuField.ySize(), gpuField.zSize(), 1, 0,
+                                                field::fzyx);
+      cpuField.setWithGhostLayer(0);
+
+      WALBERLA_CHECK_EQUAL(cpuField.xSize(), 2)
+      WALBERLA_CHECK_EQUAL(cpuField.ySize(), 2)
+      WALBERLA_CHECK_EQUAL(cpuField.zSize(), 2)
+
+      // initialize the bottom ghost layer cells
+      cpuField(0, 0, -1) = 1;
+      cpuField(0, 1, -1) = 2;
+      cpuField(1, 0, -1) = 3;
+      cpuField(1, 1, -1) = 4;
+
+      // initialize the top interior cells
+      cpuField(0, 0, 1) = 1;
+      cpuField(0, 1, 1) = 1;
+      cpuField(1, 0, 1) = 1;
+      cpuField(1, 1, 1) = 1;
+
+      cuda::fieldCpy(gpuField, cpuField);
+
+      // communicate pull += reduction
+      us1.communicate();
+
+      cuda::fieldCpy(cpuField, gpuField);
+
+      // check values in top ghost layer
+      WALBERLA_CHECK_EQUAL(cpuField(0, 0, 2), 0)
+      WALBERLA_CHECK_EQUAL(cpuField(0, 1, 2), 0)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 0, 2), 0)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 1, 2), 0)
+
+      // check values in top interior cells
+      WALBERLA_CHECK_EQUAL(cpuField(0, 0, 1), 2)
+      WALBERLA_CHECK_EQUAL(cpuField(0, 1, 1), 3)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 0, 1), 4)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 1, 1), 5)
+
+      // communicate to sync ghost layers
+      us2.communicate();
+
+      cuda::fieldCpy(cpuField, gpuField);
+
+      // check values in bottom ghost layer
+      WALBERLA_CHECK_EQUAL(cpuField(0, 0, -1), 2)
+      WALBERLA_CHECK_EQUAL(cpuField(0, 1, -1), 3)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 0, -1), 4)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 1, -1), 5)
+
+      // check values in top interior cells
+      WALBERLA_CHECK_EQUAL(cpuField(0, 0, 1), 2)
+      WALBERLA_CHECK_EQUAL(cpuField(0, 1, 1), 3)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 0, 1), 4)
+      WALBERLA_CHECK_EQUAL(cpuField(1, 1, 1), 5)
+
+   }
+}
+
+int main(int argc, char **argv) {
+
+   using blockforest::createUniformBlockGrid;
+
+   debug::enterTestMode();
+   Environment walberlaEnv(argc,argv);
+
+   // Create a BlockForest with 2x2x2 cells per block
+   uint_t processes = uint_c( MPIManager::instance()->numProcesses() );
+   auto blocks = createUniformBlockGrid(processes,1 ,1, //blocks
+                                        2,2,2, //cells
+                                        1, //dx
+                                        true, //one block per process
+                                        true,true,true); //periodicity
+
+   // Create a Field with the same number of cells as the block
+   BlockDataID scalarGPUFieldId = blocks->addStructuredBlockData<cuda::GPUField<int> > ( &createGPUField, "ScalarGPUField" );
+
+   testScalarField( blocks, scalarGPUFieldId );
+
+   // Create a BlockForest with 8x8x8 cells per block
+   blocks = createUniformBlockGrid(processes,1 ,1, //blocks
+                                   8,8,8,          //cells
+                                   1,              //dx
+                                   true,          //one block per process
+                                   true,true,true);//periodicity
+
+   // Create a Field with one quarter as many cells per dimension, i.e. a field with the same size as the one above
+   scalarGPUFieldId = blocks->addStructuredBlockData<cuda::GPUField<int> > ( &createSmallGPUField, "ScalarGPUField" );
+
+   testScalarField( blocks, scalarGPUFieldId );
+
+   testScalarFieldPullReduction( blocks, scalarGPUFieldId );
+
+   return 0;
+
+}
+
+} // namespace walberla
+
+int main( int argc, char* argv[] ) {
+   return walberla::main( argc, argv );
+}
\ No newline at end of file

tests/cuda/codegen/GeneratedFieldPackInfoTestGPU.py

+import operator as op
+import pystencils as ps
+from pystencils_walberla import CodeGeneration, generate_pack_info_for_field
+
+with CodeGeneration() as ctx:
+
+    layout = 'fzyx'
+
+    field = ps.fields("field: int32[3D]", layout=layout)
+
+    # communication
+    generate_pack_info_for_field(ctx, 'ScalarFieldCommunicationGPU', field, target='gpu')
+    generate_pack_info_for_field(ctx, 'ScalarFieldPullReductionGPU', field, target='gpu', operator=op.add,
+                                 gl_to_inner=True)

tests/field/CMakeLists.txt

@@ -77,4 +77,11 @@ waLBerla_generate_target_from_python(NAME CodeGenMultipleFieldSwaps FILE codegen
         OUT_FILES MultipleFieldSwaps.cpp MultipleFieldSwaps.h )
 waLBerla_compile_test( FILES codegen/MultipleFieldSwaps.cpp DEPENDS gui timeloop CodeGenMultipleFieldSwaps)
 waLBerla_execute_test( NAME MultipleFieldSwaps )
+
+waLBerla_generate_target_from_python(NAME CodegenGeneratedCPUFieldPackInfo FILE codegen/GeneratedFieldPackInfoTest.py
+        OUT_FILES ScalarFieldCommunication.cpp ScalarFieldCommunication.h
+                  ScalarFieldPullReduction.cpp ScalarFieldPullReduction.h )
+waLBerla_compile_test( FILES codegen/GeneratedFieldPackInfoTest.cpp
+        DEPENDS blockforest core field CodegenGeneratedCPUFieldPackInfo )
+waLBerla_execute_test( NAME GeneratedFieldPackInfoTest )
 endif()

tests/field/codegen/GeneratedFieldPackInfoTest.cpp

+//======================================================================================================================
+//
+//  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 GeneratedFieldPackInfoTest.cpp
+//! \ingroup field
+//! \author Martin Bauer <martin.bauer@fau.de>
+//! \author Christoph Rettinger <christoph.rettinger@fau.de>
+//! \author Helen Schottenhamml <helen.schottenhamml@fau.de>
+//! \brief Tests if a Field is correctly packed into buffers using generated pack info
+//
+//======================================================================================================================
+
+#include "field/AddToStorage.h"
+#include "field/GhostLayerField.h"
+
+#include "blockforest/Initialization.h"
+
+#include "core/debug/TestSubsystem.h"
+#include "core/Environment.h"
+
+#include "field/communication/PackInfo.h"
+#include "field/communication/UniformPullReductionPackInfo.h"
+
+#include <cstring>
+
+// include generated files
+#include "ScalarFieldCommunication.h"
+#include "ScalarFieldPullReduction.h"
+
+
+namespace walberla {
+
+void testScalarField( IBlock * block, BlockDataID fieldId )
+{
+   GhostLayerField<int,1> & field = *(block->getData<GhostLayerField<int,1> > (fieldId));
+   field.setWithGhostLayer( 0 );
+
+   WALBERLA_CHECK_EQUAL(field.xSize(), 2);
+   WALBERLA_CHECK_EQUAL(field.ySize(), 2);
+   WALBERLA_CHECK_EQUAL(field.zSize(), 2);
+
+   // initialize the bottom boundary
+   field(0,0,0) = 1;
+   field(0,1,0) = 2;
+   field(1,0,0) = 3;
+   field(1,1,0) = 4;
+
+   // -------------- Local Communication Test ----------------------
+
+   // communicate periodic from bottom to top
+   pystencils::ScalarFieldCommunication pi(fieldId);
+   pi.communicateLocal( block, block, stencil::B );
+
+   WALBERLA_CHECK_EQUAL ( field(0,0,+2), 1 );
+   WALBERLA_CHECK_EQUAL ( field(0,1,+2), 2 );
+   WALBERLA_CHECK_EQUAL ( field(1,0,+2), 3 );
+   WALBERLA_CHECK_EQUAL ( field(1,1,+2), 4 );
+
+   // -------------- Buffer Communication Test ---------------------
+
+   // Reset
+   field(0,0,2) = 0;
+   field(0,1,2) = 0;
+   field(1,0,2) = 0;
+   field(1,1,2) = 0;
+
+   mpi::GenericSendBuffer<> sendBuf;
+   pi.packData( block, stencil::B, sendBuf );
+
+   // Manually copy over the send to the receive buffer
+   mpi::GenericRecvBuffer<> recvBuf;
+   recvBuf.resize( sendBuf.size() );
+   memcpy( recvBuf.ptr(), sendBuf.ptr(), sendBuf.size()* sizeof(mpi::GenericSendBuffer<>::ElementType) );
+
+   pi.unpackData( block, stencil::T, recvBuf );
+
+   WALBERLA_CHECK_EQUAL ( field(0,0,+2), 1 );
+   WALBERLA_CHECK_EQUAL ( field(0,1,+2), 2 );
+   WALBERLA_CHECK_EQUAL ( field(1,0,+2), 3 );
+   WALBERLA_CHECK_EQUAL ( field(1,1,+2), 4 );
+}
+
+void testScalarFieldPullReduction( IBlock * block, BlockDataID fieldId )
+{
+   GhostLayerField<int,1> & field = *(block->getData<GhostLayerField<int,1> > (fieldId));
+   field.setWithGhostLayer( 0 );
+
+   WALBERLA_CHECK_EQUAL(field.xSize(), 2);
+   WALBERLA_CHECK_EQUAL(field.ySize(), 2);
+   WALBERLA_CHECK_EQUAL(field.zSize(), 2);
+
+   // initialize the bottom ghost layer cells
+   field(0,0,-1) = 1;
+   field(0,1,-1) = 2;
+   field(1,0,-1) = 3;
+   field(1,1,-1) = 4;
+
+   // initialize the top interior cells
+   field(0,0,1) = 1;
+   field(0,1,1) = 1;
+   field(1,0,1) = 1;
+   field(1,1,1) = 1;
+
+   // communicate periodic from bottom to top with uniform pull scheme
+   pystencils::ScalarFieldPullReduction pi1(fieldId);
+   pi1.communicateLocal( block, block, stencil::B );
+
+   // check values in top ghost layer
+   WALBERLA_CHECK_EQUAL ( field(0,0,2), 0 );
+   WALBERLA_CHECK_EQUAL ( field(0,1,2), 0 );
+   WALBERLA_CHECK_EQUAL ( field(1,0,2), 0 );
+   WALBERLA_CHECK_EQUAL ( field(1,1,2), 0 );
+
+   // check values in top interior cells
+   WALBERLA_CHECK_EQUAL ( field(0,0,1), 2 );
+   WALBERLA_CHECK_EQUAL ( field(0,1,1), 3 );
+   WALBERLA_CHECK_EQUAL ( field(1,0,1), 4 );
+   WALBERLA_CHECK_EQUAL ( field(1,1,1), 5 );
+
+   // communicate periodic from top to bottom with standard form to sync ghost layers
+   pystencils::ScalarFieldCommunication pi2 (fieldId);
+   pi2.communicateLocal( block, block, stencil::T );
+
+   // check values in bottom ghost layer
+   WALBERLA_CHECK_EQUAL ( field(0,0,-1), 2 );
+   WALBERLA_CHECK_EQUAL ( field(0,1,-1), 3 );
+   WALBERLA_CHECK_EQUAL ( field(1,0,-1), 4 );
+   WALBERLA_CHECK_EQUAL ( field(1,1,-1), 5 );
+
+   // check values in top interior cells
+   WALBERLA_CHECK_EQUAL ( field(0,0,1), 2 );
+   WALBERLA_CHECK_EQUAL ( field(0,1,1), 3 );
+   WALBERLA_CHECK_EQUAL ( field(1,0,1), 4 );
+   WALBERLA_CHECK_EQUAL ( field(1,1,1), 5 );
+
+}
+
+int main(int argc, char **argv)
+{
+   using blockforest::createUniformBlockGrid;
+
+   debug::enterTestMode();
+   Environment walberlaEnv(argc,argv);
+
+
+   // Create a BlockForest with 2x2x2 cells per block
+   uint_t processes = uint_c( MPIManager::instance()->numProcesses() );
+   auto blocks = createUniformBlockGrid(processes,1 ,1, //blocks
+                                        2,2,2,          //cells
+                                        1,              //dx
+                                        false,          //one block per process
+                                        true,true,true);//periodicity
+
+   // Create a Field with the same number of cells as the block
+   BlockDataID scalarFieldId = field::addToStorage<GhostLayerField<int,1> > ( blocks, "ScalarField" );
+
+   for( auto blockIt = blocks->begin(); blockIt != blocks->end(); ++blockIt ) // block loop
+      testScalarField( &(*blockIt), scalarFieldId );
+
+   // Create a BlockForest with 8x8x8 cells per block
+   blocks = createUniformBlockGrid(processes,1 ,1, //blocks
+                                   8,8,8,          //cells
+                                   1,              //dx
+                                   false,          //one block per process
+                                   true,true,true);//periodicity
+
+   // Create a Field with one quarter as many cells per dimension, i.e. a field with the same size as the one above
+   auto getSize = []( const shared_ptr< StructuredBlockStorage > &, IBlock * const ) {
+     return Vector3<uint_t>(2,2,2);