Added comments to StencilFieldCoarseningDCA and do test

src/pde/sweeps/Multigrid.impl.h

@@ -153,7 +153,57 @@ void ComputeResidualFixedStencil< Stencil_T >::operator()( IBlock * const block
 }
 
 
+//*************************************************************************************************
+/*! \brief Initialize stencil fields on coarser MG levels dependent to stencils on finest grid
+ *         by multiplying stencils with correct h-dependent scaling factor.
+ *
+ * On each grid level, the stencils are scaled by a scaling factor dependent on order of the pde.
+ * \param stencilFieldId vector of stencil field BlockDataIDs
+ *
+ * \ingroup pde
+ */
+//*************************************************************************************************
 
+template< typename Stencil_T >
+void CoarsenStencilFieldsDCA<Stencil_T >::operator()( const std::vector<BlockDataID> & stencilFieldId ) const
+{
+
+   typedef stencil::Iterator<Stencil_T> stencilIt;
+
+   WALBERLA_ASSERT_EQUAL(numLvl_, stencilFieldId.size(), "The number of stencil fields has to match multigrid levels!");
+
+   const real_t scalingFactor = real_t(1)/real_c( 1 << operatorOrder_ ); // scaling by ( 1/2^operatorOrder ) per MG level
+
+   // Iterate from second finest to coarsest MG level
+   for ( uint_t lvl = 1; lvl < numLvl_; ++lvl )
+   {
+      for( auto block = blocks_->begin( requiredSelectors_, incompatibleSelectors_ ); block != blocks_->end(); ++block )
+      {
+         StencilField_T * const fineStencils   = block-> getData< StencilField_T >( stencilFieldId[lvl-1] ); // template keyword only required for some compilers (e.g. IBM) that otherwise misinterpret "<" as smaller sign
+         StencilField_T * const coarseStencils = block->template getData< StencilField_T >( stencilFieldId[lvl] );
+
+         cell::CellInterval xyz = coarseStencils->xyzSize();
+
+         // Iterate over all cells in CellInterval that contains maximum cell indices in all spatial dimensions
+         for( cell_idx_t z = 0; z <= xyz.zMax(); ++z ) {
+            for( cell_idx_t y = 0; y <= xyz.yMax(); ++y ) {
+               for( cell_idx_t x = 0; x <= xyz.xMax(); ++x )
+               {
+                  for( stencilIt dir = Stencil_T::begin(); dir != Stencil_T::end(); ++dir )
+                     coarseStencils->get(x,y,z, dir.toIdx()) = scalingFactor * fineStencils->get(x,y,z, dir.toIdx());
+               }
+            }
+         }
+
+      }
+   }
+
+}
+
+//*************************************************************************************************
+
+/*
+ // previous version:
 template< typename Stencil_T >
 void CoarsenStencilFieldsDCA<Stencil_T >::operator()( const std::vector<BlockDataID> & stencilFieldId ) const
 {
@@ -175,8 +225,7 @@ void CoarsenStencilFieldsDCA<Stencil_T >::operator()( const std::vector<BlockDat
       }
    }
 }
-
-
+*/
 
 template< >
 void CoarsenStencilFieldsGCA< stencil::D3Q7 >::operator()( const std::vector<BlockDataID> & stencilFieldId ) const