Merge branch 'mesa' into 'master'

Updates for MESA-PD

See merge request walberla/walberla!228

apps/benchmarks/GranularGas/CMakeLists.txt

@@ -2,13 +2,25 @@ waLBerla_link_files_to_builddir( *.cfg )
 waLBerla_link_files_to_builddir( *.py )
 
 waLBerla_add_executable ( NAME PE_GranularGas
-   FILES PE_GranularGas.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp
+                          FILES PE_GranularGas.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp
                           DEPENDS blockforest core pe postprocessing sqlite )
 
+waLBerla_add_executable ( NAME PE_LoadBalancing
+                          FILES PE_LoadBalancing.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp
+                          DEPENDS blockforest core pe postprocessing sqlite )
+
+waLBerla_add_executable ( NAME MESA_PD_LoadBalancing
+                          FILES MESA_PD_LoadBalancing.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp sortParticleStorage.cpp CreateParticles.cpp
+                          DEPENDS blockforest core pe mesa_pd postprocessing sqlite vtk )
+
 waLBerla_add_executable ( NAME MESA_PD_GranularGas
-   FILES MESA_PD_GranularGas.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp sortParticleStorage.cpp CreateParticles.cpp
-                          DEPENDS blockforest core pe mesa_pd postprocessing vtk sqlite )
+                          FILES MESA_PD_GranularGas.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp sortParticleStorage.cpp CreateParticles.cpp
+                          DEPENDS blockforest core pe mesa_pd postprocessing sqlite vtk )
 
 waLBerla_add_executable ( NAME MESA_PD_KernelBenchmark
-   FILES MESA_PD_KernelBenchmark.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp sortParticleStorage.cpp CreateParticles.cpp
-                          DEPENDS blockforest core pe mesa_pd postprocessing vtk sqlite )
+                          FILES MESA_PD_KernelBenchmark.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp sortParticleStorage.cpp CreateParticles.cpp
+                          DEPENDS blockforest core pe mesa_pd postprocessing sqlite vtk )
+
+waLBerla_add_executable ( NAME MESA_PD_KernelLoadBalancing
+                          FILES MESA_PD_KernelLoadBalancing.cpp SQLProperties.cpp Parameters.cpp NodeTimings.cpp sortParticleStorage.cpp CreateParticles.cpp
+                          DEPENDS blockforest core pe mesa_pd postprocessing sqlite vtk )

apps/benchmarks/GranularGas/GenerateModule.py

@@ -49,6 +49,9 @@ if __name__ == '__main__':
    ps.addProperty("flags",            "walberla::mesa_pd::data::particle_flags::FlagT", defValue="", syncMode="COPY")
    ps.addProperty("nextParticle",     "int",                     defValue="-1",        syncMode="NEVER")
 
+   ps.addInclude("blockforest/BlockForest.h")
+   ps.addProperty("currentBlock",     "blockforest::Block*",     defValue="nullptr",   syncMode="NEVER")
+
    kernels = []
    kernels.append( kernel.DoubleCast(shapes) )
    kernels.append( kernel.ExplicitEuler() )
@@ -74,6 +77,7 @@ if __name__ == '__main__':
    comm.append(mpi.ClearNextNeighborSync())
    comm.append(mpi.ReduceContactHistory())
    comm.append(mpi.ReduceProperty())
+   comm.append(mpi.SyncGhostOwners(ps))
    comm.append(mpi.SyncNextNeighbors(ps))
 
 

apps/benchmarks/GranularGas/GranularGas.cfg

@@ -2,11 +2,21 @@ GranularGas
 {
    simulationCorner < 0, 0, 0 >;
    simulationDomain < 40, 40, 40 >;
-   blocks < 2,2,2 >;
-   isPeriodic < 1, 1, 1 >;
-   initialRefinementLevel 0;
-   sorting none;
+   blocks < 1,1,1 >;
+   isPeriodic < 0, 0, 0 >;
+   initialRefinementLevel 1;
+   sorting linear;
 
+   LBAlgorithm Morton;
+   baseWeight 1;
+
+   recalculateBlockLevelsInRefresh 1;
+   reevaluateMinTargetLevelsAfterForcedRefinement 1;
+   allowRefreshChangingDepth 1;
+   regridMin 2000;
+   regridMax 100;
+
+   normal  <1,1,1>;
    radius  0.6;
    spacing 1.0;
    vMax    0.0;

apps/benchmarks/GranularGas/MESA_PD_GranularGas.cpp

@@ -36,6 +36,7 @@
 #include <mesa_pd/data/ParticleStorage.h>
 #include <mesa_pd/data/ShapeStorage.h>
 #include <mesa_pd/domain/BlockForestDomain.h>
+#include <mesa_pd/kernel/AssocToBlock.h>
 #include <mesa_pd/kernel/DoubleCast.h>
 #include <mesa_pd/kernel/ExplicitEulerWithShape.h>
 #include <mesa_pd/kernel/InsertParticleIntoLinkedCells.h>
@@ -44,6 +45,7 @@
 #include <mesa_pd/mpi/ContactFilter.h>
 #include <mesa_pd/mpi/ReduceProperty.h>
 #include <mesa_pd/mpi/SyncNextNeighbors.h>
+#include <mesa_pd/mpi/SyncNextNeighborsBlockForest.h>
 
 #include <mesa_pd/mpi/notifications/ForceTorqueNotification.h>
 
@@ -101,7 +103,7 @@ int main( int argc, char ** argv )
       WALBERLA_LOG_INFO_ON_ROOT( "No BlockForest created ... exiting!");
       return EXIT_SUCCESS;
    }
-   domain::BlockForestDomain domain(forest);
+   auto domain = std::make_shared<domain::BlockForestDomain> (forest);
 
    auto simulationDomain = forest->getDomain();
    auto localDomain = forest->begin()->getAABB();
@@ -123,7 +125,7 @@ int main( int argc, char ** argv )
    for (auto& iBlk : *forest)
    {
       for (auto pt : grid_generator::SCGrid(iBlk.getAABB(),
-                                            Vector3<real_t>(params.spacing) * real_c(0.5),
+                                            Vector3<real_t>(params.spacing) * real_c(0.5) + params.shift,
                                             params.spacing))
       {
          WALBERLA_CHECK(iBlk.getAABB().contains(pt));
@@ -158,9 +160,9 @@ int main( int argc, char ** argv )
    WALBERLA_LOG_INFO_ON_ROOT("*** SETUP - END ***");
 
    WALBERLA_LOG_INFO_ON_ROOT("*** VTK ***");
-   auto vtkDomainOutput = walberla::vtk::createVTKOutput_DomainDecomposition( forest, "domain_decomposition", 1, "vtk_out", "simulation_step" );
+   auto vtkDomainOutput = walberla::vtk::createVTKOutput_DomainDecomposition( forest, "domain_decomposition", 1, params.vtk_out, "simulation_step" );
    auto vtkOutput       = make_shared<mesa_pd::vtk::ParticleVtkOutput>(ps) ;
-   auto vtkWriter       = walberla::vtk::createVTKOutput_PointData(vtkOutput, "Bodies", 1, "vtk", "simulation_step", false, false);
+   auto vtkWriter       = walberla::vtk::createVTKOutput_PointData(vtkOutput, "Bodies", 1, params.vtk_out, "simulation_step", false, false);
    vtkOutput->addOutput<SelectRank>("rank");
    vtkOutput->addOutput<data::SelectParticleOwner>("owner");
    //   vtkDomainOutput->write();
@@ -175,13 +177,15 @@ int main( int argc, char ** argv )
    dem.setDampingT (0, 0, real_t(0));
    dem.setFriction (0, 0, real_t(0));
    collision_detection::AnalyticContactDetection              acd;
+   kernel::AssocToBlock                  assoc(forest);
    kernel::DoubleCast                    double_cast;
    mpi::ContactFilter                    contact_filter;
    mpi::ReduceProperty                   RP;
-   mpi::SyncNextNeighbors                SNN;
+   mpi::SyncNextNeighborsBlockForest     SNN;
 
    // initial sync
-   SNN(*ps, domain);
+   ps->forEachParticle(false, kernel::SelectLocal(), accessor, assoc, accessor);
+   SNN(*ps, forest, domain);
    sortParticleStorage(*ps, params.sorting, lc.domain_, uint_c(lc.numCellsPerDim_[0]));
 //   vtkWriter->write();
 
@@ -211,6 +215,11 @@ int main( int argc, char ** argv )
          //         vtkWriter->write();
          //      }
 
+         tp["AssocToBlock"].start();
+         ps->forEachParticle(false, kernel::SelectLocal(), accessor, assoc, accessor);
+         if (params.bBarrier) WALBERLA_MPI_BARRIER();
+         tp["AssocToBlock"].end();
+
          tp["GenerateLinkedCells"].start();
          lc.clear();
          ps->forEachParticle(true, kernel::SelectAll(), accessor, ipilc, accessor, lc);
@@ -230,7 +239,7 @@ int main( int argc, char ** argv )
             if (double_cast(idx1, idx2, ac, acd, ac ))
             {
                ++contactsDetected;
-               if (contact_filter(acd.getIdx1(), acd.getIdx2(), ac, acd.getContactPoint(), domain))
+               if (contact_filter(acd.getIdx1(), acd.getIdx2(), ac, acd.getContactPoint(), *domain))
                {
                   ++contactsTreated;
                   dem(acd.getIdx1(), acd.getIdx2(), ac, acd.getContactPoint(), acd.getContactNormal(), acd.getPenetrationDepth());
@@ -253,7 +262,7 @@ int main( int argc, char ** argv )
          tp["Euler"].end();
 
          tp["SNN"].start();
-         SNN(*ps, domain);
+         SNN(*ps, forest, domain);
          if (params.bBarrier) WALBERLA_MPI_BARRIER();
          tp["SNN"].end();
       }
@@ -329,9 +338,9 @@ int main( int argc, char ** argv )
       walberla::mpi::reduceInplace(linkedCellsVolume, walberla::mpi::SUM);
       size_t numLinkedCells = lc.cells_.size();
       walberla::mpi::reduceInplace(numLinkedCells, walberla::mpi::SUM);
-      size_t local_aabbs         = domain.getNumLocalAABBs();
-      size_t neighbor_subdomains = domain.getNumNeighborSubdomains();
-      size_t neighbor_processes  = domain.getNumNeighborProcesses();
+      size_t local_aabbs         = domain->getNumLocalAABBs();
+      size_t neighbor_subdomains = domain->getNumNeighborSubdomains();
+      size_t neighbor_processes  = domain->getNumNeighborProcesses();
       walberla::mpi::reduceInplace(local_aabbs, walberla::mpi::SUM);
       walberla::mpi::reduceInplace(neighbor_subdomains, walberla::mpi::SUM);
       walberla::mpi::reduceInplace(neighbor_processes, walberla::mpi::SUM);
@@ -347,6 +356,11 @@ int main( int argc, char ** argv )
          stringProperties["tag"]                  = "mesa_pd";
          integerProperties["mpi_num_processes"]   = mpiManager->numProcesses();
          integerProperties["omp_max_threads"]     = omp_get_max_threads();
+         realProperties["PUpS"]                   = double_c(PUpS);
+         realProperties["timer_min"]              = timer_reduced->min();
+         realProperties["timer_max"]              = timer_reduced->max();
+         realProperties["timer_average"]          = timer_reduced->average();
+         realProperties["timer_total"]            = timer_reduced->total();
          integerProperties["outerIteration"]      = int64_c(outerIteration);
          integerProperties["num_particles"]       = numParticles;
          integerProperties["num_ghost_particles"] = numGhostParticles;

apps/benchmarks/GranularGas/MESA_PD_KernelBenchmark.cpp

@@ -128,7 +128,7 @@ int main( int argc, char ** argv )
    for (auto& iBlk : *forest)
    {
       for (auto pt : grid_generator::SCGrid(iBlk.getAABB(),
-                                            Vector3<real_t>(params.spacing) * real_c(0.5),
+                                            Vector3<real_t>(params.spacing) * real_c(0.5) + params.shift,
                                             params.spacing))
       {
          WALBERLA_CHECK(iBlk.getAABB().contains(pt));
@@ -163,9 +163,9 @@ int main( int argc, char ** argv )
    WALBERLA_LOG_INFO_ON_ROOT("*** SETUP - END ***");
 
    WALBERLA_LOG_INFO_ON_ROOT("*** VTK ***");
-   auto vtkDomainOutput = walberla::vtk::createVTKOutput_DomainDecomposition( forest, "domain_decomposition", 1, "vtk_out", "simulation_step" );
+   auto vtkDomainOutput = walberla::vtk::createVTKOutput_DomainDecomposition( forest, "domain_decomposition", 1, params.vtk_out, "simulation_step" );
    auto vtkOutput       = make_shared<mesa_pd::vtk::ParticleVtkOutput>(ps) ;
-   auto vtkWriter       = walberla::vtk::createVTKOutput_PointData(vtkOutput, "Bodies", 1, "vtk", "simulation_step", false, false);
+   auto vtkWriter       = walberla::vtk::createVTKOutput_PointData(vtkOutput, "Bodies", 1, params.vtk_out, "simulation_step", false, false);
    vtkOutput->addOutput<SelectRank>("rank");
    vtkOutput->addOutput<data::SelectParticleOwner>("owner");
    vtkOutput->addOutput<SelectIdx>("idx");

apps/benchmarks/GranularGas/PE_GranularGas.cpp

@@ -90,7 +90,7 @@ int main( int argc, char ** argv )
    auto cfg = env.config();
    if (cfg == nullptr) WALBERLA_ABORT("No config specified!");
    const Config::BlockHandle mainConf  = cfg->getBlock( "GranularGas" );
-   Parameters params;
+   mesa_pd::Parameters params;
    loadFromConfig(params, mainConf);
 
    WALBERLA_LOG_INFO_ON_ROOT("*** GLOBALBODYSTORAGE ***");
@@ -174,9 +174,9 @@ int main( int argc, char ** argv )
    }
 
    WALBERLA_LOG_INFO_ON_ROOT("*** VTK ***");
-   auto vtkDomainOutput = vtk::createVTKOutput_DomainDecomposition( forest, "domain_decomposition", 1, "vtk_out", "simulation_step" );
+   auto vtkDomainOutput = vtk::createVTKOutput_DomainDecomposition( forest, "domain_decomposition", 1, params.vtk_out, "simulation_step" );
    auto vtkSphereHelper = make_shared<SphereVtkOutput>(storageID, *forest) ;
-   auto vtkSphereOutput = vtk::createVTKOutput_PointData(vtkSphereHelper, "Bodies", 1, "vtk_out", "simulation_step", false, false);
+   auto vtkSphereOutput = vtk::createVTKOutput_PointData(vtkSphereHelper, "Bodies", 1, params.vtk_out, "simulation_step", false, false);
 
    WALBERLA_LOG_INFO_ON_ROOT("*** SETUP - START ***");
    //const real_t   static_cof  ( real_c(0.1) / 2 );   // Coefficient of static friction. Note: pe doubles the input coefficient of friction for material-material contacts.
@@ -190,7 +190,7 @@ int main( int argc, char ** argv )
    for (auto& currentBlock : *forest)
    {
       for (auto it = grid_generator::SCIterator(currentBlock.getAABB().getIntersection(generationDomain),
-                                                Vector3<real_t>(params.spacing) * real_c(0.5),
+                                                Vector3<real_t>(params.spacing) * real_c(0.5) + params.shift,
                                                 params.spacing);
            it != grid_generator::SCIterator();
            ++it)

apps/benchmarks/GranularGas/Parameters.cpp

@@ -29,15 +29,22 @@
 #include <core/logging/Logging.h>
 
 namespace walberla {
+namespace mesa_pd {
 
 void loadFromConfig(Parameters& params, const Config::BlockHandle& cfg)
 {
    params.sorting = cfg.getParameter<std::string>("sorting", "none" );
    WALBERLA_LOG_INFO_ON_ROOT("sorting: " << params.sorting);
    
+   params.normal = cfg.getParameter<Vec3>("normal", Vec3(real_t(1.0), real_t(1.0), real_t(1.0)) );
+   WALBERLA_LOG_INFO_ON_ROOT("normal: " << params.normal);
+   
    params.spacing = cfg.getParameter<real_t>("spacing", real_t(1.0) );
    WALBERLA_LOG_INFO_ON_ROOT("spacing: " << params.spacing);
    
+   params.shift = cfg.getParameter<Vec3>("shift", Vec3(real_t(0.1), real_t(0.1), real_t(0.1)) );
+   WALBERLA_LOG_INFO_ON_ROOT("shift: " << params.shift);
+   
    params.radius = cfg.getParameter<real_t>("radius", real_t(0.5) );
    WALBERLA_LOG_INFO_ON_ROOT("radius: " << params.radius);
    
@@ -65,12 +72,60 @@ void loadFromConfig(Parameters& params, const Config::BlockHandle& cfg)
    params.visSpacing = cfg.getParameter<int64_t>("visSpacing", 1000 );
    WALBERLA_LOG_INFO_ON_ROOT("visSpacing: " << params.visSpacing);
    
-   params.path = cfg.getParameter<std::string>("path", "vtk_out" );
-   WALBERLA_LOG_INFO_ON_ROOT("path: " << params.path);
+   params.vtk_out = cfg.getParameter<std::string>("vtk_out", "vtk_out" );
+   WALBERLA_LOG_INFO_ON_ROOT("vtk_out: " << params.vtk_out);
    
    params.sqlFile = cfg.getParameter<std::string>("sqlFile", "benchmark.sqlite" );
    WALBERLA_LOG_INFO_ON_ROOT("sqlFile: " << params.sqlFile);
    
+   params.recalculateBlockLevelsInRefresh = cfg.getParameter<bool>("recalculateBlockLevelsInRefresh", false );
+   WALBERLA_LOG_INFO_ON_ROOT("recalculateBlockLevelsInRefresh: " << params.recalculateBlockLevelsInRefresh);
+   
+   params.alwaysRebalanceInRefresh = cfg.getParameter<bool>("alwaysRebalanceInRefresh", true );
+   WALBERLA_LOG_INFO_ON_ROOT("alwaysRebalanceInRefresh: " << params.alwaysRebalanceInRefresh);
+   
+   params.reevaluateMinTargetLevelsAfterForcedRefinement = cfg.getParameter<bool>("reevaluateMinTargetLevelsAfterForcedRefinement", false );
+   WALBERLA_LOG_INFO_ON_ROOT("reevaluateMinTargetLevelsAfterForcedRefinement: " << params.reevaluateMinTargetLevelsAfterForcedRefinement);
+   
+   params.allowRefreshChangingDepth = cfg.getParameter<bool>("allowRefreshChangingDepth", false );
+   WALBERLA_LOG_INFO_ON_ROOT("allowRefreshChangingDepth: " << params.allowRefreshChangingDepth);
+   
+   params.allowMultipleRefreshCycles = cfg.getParameter<bool>("allowMultipleRefreshCycles", false );
+   WALBERLA_LOG_INFO_ON_ROOT("allowMultipleRefreshCycles: " << params.allowMultipleRefreshCycles);
+   
+   params.checkForEarlyOutInRefresh = cfg.getParameter<bool>("checkForEarlyOutInRefresh", true );
+   WALBERLA_LOG_INFO_ON_ROOT("checkForEarlyOutInRefresh: " << params.checkForEarlyOutInRefresh);
+   
+   params.checkForLateOutInRefresh = cfg.getParameter<bool>("checkForLateOutInRefresh", true );
+   WALBERLA_LOG_INFO_ON_ROOT("checkForLateOutInRefresh: " << params.checkForLateOutInRefresh);
+   
+   params.regridMin = cfg.getParameter<uint_t>("regridMin", uint_c(100) );
+   WALBERLA_LOG_INFO_ON_ROOT("regridMin: " << params.regridMin);
+   
+   params.regridMax = cfg.getParameter<uint_t>("regridMax", uint_c(1000) );
+   WALBERLA_LOG_INFO_ON_ROOT("regridMax: " << params.regridMax);
+   
+   params.maxBlocksPerProcess = cfg.getParameter<int>("maxBlocksPerProcess", int_c(1000) );
+   WALBERLA_LOG_INFO_ON_ROOT("maxBlocksPerProcess: " << params.maxBlocksPerProcess);
+   
+   params.baseWeight = cfg.getParameter<real_t>("baseWeight", real_t(10.0) );
+   WALBERLA_LOG_INFO_ON_ROOT("baseWeight: " << params.baseWeight);
+   
+   params.metisipc2redist = cfg.getParameter<real_t>("metisipc2redist", real_t(1000.0) );
+   WALBERLA_LOG_INFO_ON_ROOT("metisipc2redist: " << params.metisipc2redist);
+   
+   params.LBAlgorithm = cfg.getParameter<std::string>("LBAlgorithm", "Hilbert" );
+   WALBERLA_LOG_INFO_ON_ROOT("LBAlgorithm: " << params.LBAlgorithm);
+   
+   params.metisAlgorithm = cfg.getParameter<std::string>("metisAlgorithm", "PART_GEOM_KWAY" );
+   WALBERLA_LOG_INFO_ON_ROOT("metisAlgorithm: " << params.metisAlgorithm);
+   
+   params.metisWeightsToUse = cfg.getParameter<std::string>("metisWeightsToUse", "BOTH_WEIGHTS" );
+   WALBERLA_LOG_INFO_ON_ROOT("metisWeightsToUse: " << params.metisWeightsToUse);
+   
+   params.metisEdgeSource = cfg.getParameter<std::string>("metisEdgeSource", "EDGES_FROM_EDGE_WEIGHTS" );
+   WALBERLA_LOG_INFO_ON_ROOT("metisEdgeSource: " << params.metisEdgeSource);
+   
 }
 
 void saveToSQL(const Parameters& params,
@@ -80,8 +135,10 @@ void saveToSQL(const Parameters& params,
 {
    stringProperties["sorting"] = params.sorting;
    
+   
    realProperties["spacing"] = double_c(params.spacing);
    
+   
    realProperties["radius"] = double_c(params.radius);
    
    
@@ -97,10 +154,33 @@ void saveToSQL(const Parameters& params,
    
    integerProperties["visSpacing"] = params.visSpacing;
    
-   stringProperties["path"] = params.path;
+   stringProperties["vtk_out"] = params.vtk_out;
    
    stringProperties["sqlFile"] = params.sqlFile;
    
+   
+   
+   
+   
+   
+   
+   
+   
+   
+   
+   realProperties["baseWeight"] = double_c(params.baseWeight);
+   
+   realProperties["metisipc2redist"] = double_c(params.metisipc2redist);
+   
+   stringProperties["LBAlgorithm"] = params.LBAlgorithm;
+   
+   stringProperties["metisAlgorithm"] = params.metisAlgorithm;
+   
+   stringProperties["metisWeightsToUse"] = params.metisWeightsToUse;
+   
+   stringProperties["metisEdgeSource"] = params.metisEdgeSource;
+   
 }
 
+} //namespace mesa_pd
 } //namespace walberla
\ No newline at end of file

apps/benchmarks/GranularGas/Parameters.h

@@ -27,16 +27,19 @@
 #pragma once
 
 #include <core/config/Config.h>
-#include <core/DataTypes.h>
+#include <mesa_pd/data/DataTypes.h>
 
 #include <string>
 
 namespace walberla {
+namespace mesa_pd {
 
 struct Parameters
 {
    std::string sorting = "none";
+   Vec3 normal = Vec3(real_t(1.0), real_t(1.0), real_t(1.0));
    real_t spacing = real_t(1.0);
+   Vec3 shift = Vec3(real_t(0.1), real_t(0.1), real_t(0.1));
    real_t radius = real_t(0.5);
    bool bBarrier = false;
    bool storeNodeTimings = false;
@@ -46,8 +49,24 @@ struct Parameters
    int64_t simulationSteps = 10;
    real_t dt = real_t(0.01);
    int64_t visSpacing = 1000;
-   std::string path = "vtk_out";
+   std::string vtk_out = "vtk_out";
    std::string sqlFile = "benchmark.sqlite";
+   bool recalculateBlockLevelsInRefresh = false;
+   bool alwaysRebalanceInRefresh = true;
+   bool reevaluateMinTargetLevelsAfterForcedRefinement = false;
+   bool allowRefreshChangingDepth = false;
+   bool allowMultipleRefreshCycles = false;
+   bool checkForEarlyOutInRefresh = true;
+   bool checkForLateOutInRefresh = true;
+   uint_t regridMin = uint_c(100);
+   uint_t regridMax = uint_c(1000);
+   int maxBlocksPerProcess = int_c(1000);
+   real_t baseWeight = real_t(10.0);
+   real_t metisipc2redist = real_t(1000.0);
+   std::string LBAlgorithm = "Hilbert";
+   std::string metisAlgorithm = "PART_GEOM_KWAY";
+   std::string metisWeightsToUse = "BOTH_WEIGHTS";
+   std::string metisEdgeSource = "EDGES_FROM_EDGE_WEIGHTS";
 };
 
 void loadFromConfig(Parameters& params,
@@ -58,4 +77,5 @@ void saveToSQL(const Parameters& params,
                std::map< std::string, double >&            realProperties,
                std::map< std::string, std::string >&       stringProperties );
 
+} //namespace mesa_pd
 } //namespace walberla
\ No newline at end of file

apps/benchmarks/GranularGas/Parameters.templ.cpp

@@ -29,6 +29,7 @@
 #include <core/logging/Logging.h>
 
 namespace walberla {
+namespace mesa_pd {
 
 void loadFromConfig(Parameters& params, const Config::BlockHandle& cfg)
 {
@@ -56,4 +57,5 @@ void saveToSQL(const Parameters& params,
    {% endfor %}
 }
 
+} //namespace mesa_pd
 } //namespace walberla

apps/benchmarks/GranularGas/Parameters.templ.h

@@ -27,11 +27,12 @@
 #pragma once
 
 #include <core/config/Config.h>
-#include <core/DataTypes.h>
+#include <mesa_pd/data/DataTypes.h>
 
 #include <string>
 
 namespace walberla {
+namespace mesa_pd {
 
 struct Parameters
 {
@@ -48,4 +49,5 @@ void saveToSQL(const Parameters& params,
                std::map< std::string, double >&            realProperties,
                std::map< std::string, std::string >&       stringProperties );
 
+} //namespace mesa_pd
 } //namespace walberla

apps/benchmarks/GranularGas/SQLProperties.cpp

@@ -108,6 +108,7 @@ void addSlurmPropertiesToSQL( std::map< std::string, int64_t > &        /*intege
    stringProperties["SLURM_NTASKS_PER_CORE"]    = envToString(std::getenv( "SLURM_NTASKS_PER_CORE" ));
    stringProperties["SLURM_NTASKS_PER_NODE"]    = envToString(std::getenv( "SLURM_NTASKS_PER_NODE" ));
    stringProperties["SLURM_NTASKS_PER_SOCKET"]  = envToString(std::getenv( "SLURM_NTASKS_PER_SOCKET" ));
+   stringProperties["SLURM_CPU_BIND_TYPE"]      = envToString(std::getenv( "SLURM_CPU_BIND_TYPE" ));
 }
 
 } //namespace mesa_pd

apps/benchmarks/GranularGas/generateConfig.py

@@ -5,7 +5,9 @@ from ConfigGenerator import Config
 
 cfg = Config()
 cfg.addParameter("sorting",                "std::string", '"none"')
+cfg.addParameter("normal",                 "Vec3",        "Vec3(real_t(1.0), real_t(1.0), real_t(1.0))")
 cfg.addParameter("spacing",                "real_t",      "real_t(1.0)")
+cfg.addParameter("shift",                  "Vec3",        "Vec3(real_t(0.1), real_t(0.1), real_t(0.1))")
 cfg.addParameter("radius",                 "real_t",      "real_t(0.5)")
 cfg.addParameter("bBarrier",               "bool",        "false")
 cfg.addParameter("storeNodeTimings",       "bool",        "false")
@@ -15,7 +17,26 @@ cfg.addParameter("initialRefinementLevel", "int64_t",     "0")
 cfg.addParameter("simulationSteps",        "int64_t",     "10")
 cfg.addParameter("dt",                     "real_t",      "real_t(0.01)")
 cfg.addParameter("visSpacing",             "int64_t",     "1000")
-cfg.addParameter("path",                   "std::string", '"vtk_out"')
+cfg.addParameter("vtk_out",                "std::string", '"vtk_out"')
 cfg.addParameter("sqlFile",                "std::string", '"benchmark.sqlite"')
 
+cfg.addParameter("recalculateBlockLevelsInRefresh",                "bool", "false");
+cfg.addParameter("alwaysRebalanceInRefresh",                       "bool", "true");
+cfg.addParameter("reevaluateMinTargetLevelsAfterForcedRefinement", "bool", "false");
+cfg.addParameter("allowRefreshChangingDepth",                      "bool", "false");
+
+cfg.addParameter("allowMultipleRefreshCycles",                     "bool", "false");
+cfg.addParameter("checkForEarlyOutInRefresh",                      "bool", "true");
+cfg.addParameter("checkForLateOutInRefresh",                       "bool", "true");
+
+cfg.addParameter("regridMin",              "uint_t",      'uint_c(100)')
+cfg.addParameter("regridMax",              "uint_t",      'uint_c(1000)')
+cfg.addParameter("maxBlocksPerProcess",    "int",         'int_c(1000)')
+cfg.addParameter("baseWeight",             "real_t",      'real_t(10.0)')
+cfg.addParameter("metisipc2redist",        "real_t",      'real_t(1000.0)')
+cfg.addParameter("LBAlgorithm",            "std::string", '"Hilbert"')
+cfg.addParameter("metisAlgorithm",         "std::string", '"PART_GEOM_KWAY"' );
+cfg.addParameter("metisWeightsToUse",      "std::string", '"BOTH_WEIGHTS"' );
+cfg.addParameter("metisEdgeSource",        "std::string", '"EDGES_FROM_EDGE_WEIGHTS"' );
+
 cfg.generate()

python/mesa_pd.py

@@ -32,6 +32,7 @@ if __name__ == '__main__':
    ps    = data.ParticleStorage()
    ch    = data.ContactHistory()
    lc    = data.LinkedCells()
+   slc   = data.SparseLinkedCells()
    ss    = data.ShapeStorage(ps, shapes)
    cs    = data.ContactStorage()
 
@@ -47,6 +48,9 @@ if __name__ == '__main__':
    ps.addProperty("torque",           "walberla::mesa_pd::Vec3", defValue="real_t(0)", syncMode="NEVER")
    ps.addProperty("oldTorque",        "walberla::mesa_pd::Vec3", defValue="real_t(0)", syncMode="MIGRATION")
 
+   ps.addInclude("blockforest/BlockForest.h")
+   ps.addProperty("currentBlock",     "blockforest::Block*",     defValue="nullptr",   syncMode="NEVER")
+
    ps.addProperty("type",             "uint_t",                  defValue="0",         syncMode="COPY")
 
    ps.addProperty("flags",            "walberla::mesa_pd::data::particle_flags::FlagT", defValue="", syncMode="COPY")
@@ -96,6 +100,7 @@ if __name__ == '__main__':
    kernels.append( kernel.InitContactsForHCSITS() )
    kernels.append( kernel.IntegrateParticlesHCSITS() )
    kernels.append( kernel.InsertParticleIntoLinkedCells() )
+   kernels.append( kernel.InsertParticleIntoSparseLinkedCells() )
    kernels.append( kernel.LinearSpringDashpot() )
    kernels.append( kernel.NonLinearSpringDashpot() )
    kernels.append( kernel.SingleCast(shapes) )
@@ -122,6 +127,7 @@ if __name__ == '__main__':
    ps.generate(args.path + "/src/mesa_pd/")
    ch.generate(args.path + "/src/mesa_pd/")
    lc.generate(args.path + "/src/mesa_pd/")
+   slc.generate(args.path + "/src/mesa_pd/")
    ss.generate(args.path + "/src/mesa_pd/")
    cs.generate(args.path + "/src/mesa_pd/")
 

python/mesa_pd/data/SparseLinkedCells.py

+# -*- coding: utf-8 -*-
+
+import numpy as np
+from ..utility import generateFile
+
+class SparseLinkedCells:
+   def generate(self, path):
+      generateFile(path, 'data/SparseLinkedCells.templ.h')

python/mesa_pd/data/__init__.py

@@ -5,9 +5,11 @@ from .ContactStorage import ContactStorage
 from .LinkedCells import LinkedCells
 from .ParticleStorage import ParticleStorage
 from .ShapeStorage import ShapeStorage
+from .SparseLinkedCells import SparseLinkedCells
 
 __all__ = ['ContactHistory',
            'ContactStorage',
            'GeometryStorage',
            'LinkedCells',
-           'ParticleStorage']
+           'ParticleStorage',
+           'SparseLinkedCells']

python/mesa_pd/kernel/InsertParticleIntoSparseLinkedCells.py

+# -*- coding: utf-8 -*-
+
+from mesa_pd.accessor import Accessor
+from mesa_pd.utility import generateFile
+
+class InsertParticleIntoSparseLinkedCells:
+   def __init__(self):
+      self.accessor = Accessor()
+      self.accessor.require("position",     "walberla::mesa_pd::Vec3",                        access="g")
+      self.accessor.require("flags",        "walberla::mesa_pd::data::particle_flags::FlagT", access="g")
+      self.accessor.require("nextParticle", "size_t",                                     access="gs" )
+
+   def getRequirements(self):
+      return self.accessor
+