Tutorial - Useful Features

This tutorial will introduce some useful features of the waLBerla framework which can make your life easier.


You can checkpoint the current state of your rigid body dynamics simulation at any point to restore it afterwards. First you have to store the current domain partitioning using blockforest::BlockForest::saveToFile().

auto forest = blockforest::createBlockForest( math::AABB(0,0,0,60,60,60),
Vector3<uint_t>(2,2,2), // number of blocks
Vector3<bool>(false, false, false)); // periodicity

Then you have to store the current simulation data using domain_decomposition::BlockStorage::saveBlockData().

forest->saveBlockData("SerializeDeserialize.dump", storageID);

This will store all non global rigid bodies to the file system.

To load everything again you start by creating the blockforest::BlockForest. This time you will use a different constructor.

auto forest = make_shared< BlockForest >( uint_c( MPIManager::instance()->rank() ), "SerializeDeserialize.sbf", true, false );

Instead of initializing the Storage BlockDatum like you normally would

auto storageID = forest->addBlockData(createStorageDataHandling<BodyTuple>(), "Storage");

you have to use domain_decomposition::BlockStorage::loadBlockData()

auto storageID = forest->loadBlockData("SerializeDeserialize.dump", createStorageDataHandling<BodyTuple>(), "Storage");

Unfortunately due to a misorder in the loading scheme you have to reload your coarse collision detection.

for (auto blockIt = forest->begin(); blockIt != forest->end(); ++blockIt)
ccd::ICCD* ccd = blockIt->getData< ccd::ICCD >( ccdID );

Hopefully this gets fixed in the future. ;)

This method does not save global bodies nor solver settings. You have to take care to restore these settings on your own.

A fully working example can be found in the SerializeDeserialize.cpp test of the pe module.

VTK Output

For VTK Output you have to create vtk::VTKOutput objects. To output the domain partitioning use vtk::createVTKOutput_DomainDecomposition.

auto vtkDomainOutput = vtk::createVTKOutput_DomainDecomposition( forest, "domain_decomposition", 1, "vtk_out", "simulation_step" );

To output all sphere particles use vtk::createVTKOutput_PointData in conjunction with SphereVtkOutput:

auto vtkSphereHelper = make_shared<SphereVtkOutput>(storageID, *forest) ;
auto vtkSphereOutput = vtk::createVTKOutput_PointData(vtkSphereHelper, "Bodies", 1, "vtk_out", "simulation_step", false, false);

Currently only spheres are supported for VTK output but you can easily write your own SphereVtkOutput and adapt it to the body you like.

To actually write something to disc call vtk::VTKOutput::write():

vtkDomainOutput->write( );
vtkSphereOutput->write( );

You can call this every time step if you want. The files will be automatically numbered so that ParaView can generate an animation.

Loading from Config

You can specify a config file as the first command line parameter. To access it you can use the Environment::config() function. You can access subblocks of the config with config::Config::getBlock().

auto cfg = env.config();
if (cfg == nullptr) WALBERLA_ABORT("No config specified!");
const Config::BlockHandle configBlock = cfg->getBlock( "LoadFromConfig" );

To get values from the config call config::Config::getParameter():

real_t radius = configBlock.getParameter<real_t>("radius", real_c(0.4) );

Certain tasks already have predefined loading functions. You can for example directly create a BlockForest from the config file.

shared_ptr<BlockForest> forest = blockforest::createBlockForestFromConfig( configBlock );

The corresponding block in the config file looks like:

simulationCorner < -15, -15, 0 >;
simulationDomain < 12, 23, 34 >;
blocks < 3, 4, 5 >;
isPeriodic < 0, 1, 0 >;

Also the HardContact solver can be configured directly from the config file:

BlockDataID blockDataID;
cr::HCSITS hcsits( globalBodyStorage, forest, blockDataID, blockDataID, blockDataID);
configure(configBlock, hcsits);

The config file looks like:

HCSITSmaxIterations 123;
HCSITSRelaxationParameter 0.123;
HCSITSErrorReductionParameter 0.123;
HCSITSRelaxationModelStr ApproximateInelasticCoulombContactByDecoupling;
globalLinearAcceleration < 1, -2, 3 >;


To get additional information where your application spends its time you can use the WcTimingTree. It will give you a hierarchical view of the time used. Usage example:

//WcTimingTree tt;
tt.start("Initial Sync");
tt.stop("Initial Sync");

Before you output the information you should collect all the information from all the processes if you are running in parallel.

auto temp = tt.getReduced( );
std::cout << temp;

Many built-in functions like solver or synchronization methods come with an additional parameter where you can specify your timing tree. They will then include detailed information in your timing tree.

SQLite Output

waLBerla also supports SQLite databases for simulation data output. This can come in handy in parallel simulations as well as in data analysis. To store information in a SQLite database you have to fill three property maps depending on the type of information you want to store.

std::map< std::string, walberla::int64_t > integerProperties;
std::map< std::string, double > realProperties;
std::map< std::string, std::string > stringProperties;

You can then dump the information to disc. timing::TimingPool and timing::TimingTree already have predefined save functions so you do not have to extract all the information yourself and save it in the property array.

auto runId = sqlite::storeRunInSqliteDB( sqlFile, integerProperties, stringProperties, realProperties );
sqlite::storeTimingPoolInSqliteDB( sqlFile, runId, *tpReduced, "Timeloop" );
sqlite::storeTimingTreeInSqliteDB( sqlFile, runId, tt, "TimingTree" );

Image Output

Using the pe::raytracing::Raytracer you can generate images of the simulation for each timestep and output them as PNG files. Setup the raytracer by reading a config object and optionally supply a shading and a visibility function, as it is done in the second pe tutorial. The shading function will be called during rendering for each body to apply user defined coloring to bodies, the visibility function to determine if a body should be visible or not.

std::function<ShadingParameters (const BodyID body)> customShadingFunction = [](const BodyID body) {
if (body->getTypeID() == Sphere::getStaticTypeID()) {
Raytracer raytracer(forest, storageID, globalBodyStorage, ccdID,

Alternatively it may also be setup entirely in code:

Lighting lighting(Vec3(-12, 12, 12),
Color(1, 1, 1),
Color(1, 1, 1),
Color(real_t(0.4), real_t(0.4), real_t(0.4)));
Raytracer raytracer(forest, storageID, globalBodyStorage, ccdID,
640, 480,
real_t(49.13), 2,
Vec3(-25, 10, 10), Vec3(-5, 10, 10), Vec3(0, 0, 1),
Color(real_t(0.1), real_t(0.1), real_t(0.1)),

After the configuration is done, images can be generated each timestep by calling Raytracer::generateImage<BodyTuple>() which will be output to the specified directory.

raytracer.generateImage<BodyTuple>(size_t(i), &tt);

To hide certain bodies during rendering, the visibility function will be called with a BodyID as its sole argument and should return true if the object is supposed to be visible, false if not:

// [...]
std::function<bool (const BodyID body)> customVisibilityFunction = [](const BodyID body) {
if (body->getTypeID() == Plane::getStaticTypeID()) {
// hide all planes
return false;
return true;
Raytracer raytracer(forest, storageID, globalBodyStorage, ccdID,

For an overview over predefined shading functions, visit the file ShadingFunctions.h. For further information see the documentation for the classes pe::raytracing::Raytracer, pe::raytracing::Lighting and pe::raytracing::Color, the pe::raytracing::ShadingParameters struct and ShadingFunctions.h file may also be useful.

shared_ptr< BlockForest > createBlockForest(const AABB &domainAABB, const uint_t numberOfXBlocks, const uint_t numberOfYBlocks, const uint_t numberOfZBlocks, const uint_t numberOfXProcesses, const uint_t numberOfYProcesses, const uint_t numberOfZProcesses, const bool xPeriodic, const bool yPeriodic, const bool zPeriodic, const bool keepGlobalBlockInformation)
Function for creating a block forest that represents a uniform block grid.
Definition: Initialization.cpp:203
uint_t storeRunInSqliteDB(const string &dbFile, const map< string, int > &integerProperties, const map< string, string > &stringProperties, const map< string, double > &realProperties, const int busyTimeout)
Definition: SQLite.cpp:457
Definition: MPIManager.h:174
RigidBody * BodyID
Handle for a rigid body.
Definition: Types.h:81
std::tuple< Sphere, Plane > BodyTuple
Definition: 02_ConfinedGasExtended.cpp:46
#define WALBERLA_ABORT(msg)
Definition: Abort.h:62
ShadingParameters & makeGlossy(real_t _shininess=30)
Makes a rendered object shiny by setting the shininess and adjusting the specularColor.
Definition: ShadingParameters.h:68
math::Vector3< real_t > Vec3
Definition: Types.h:46
std::size_t size_t
Definition: DataTypes.h:134
shared_ptr< BlockForest > createBlockForestFromConfig(const Config::BlockHandle &mainConf, const bool keepGlobalBlockInformation)
Definition: Initialization.cpp:364
HardContactSemiImplicitTimesteppingSolvers HCSITS
Definition: HCSITS.h:273
ShadingParameters processRankDependentShadingParams(const BodyID body)
Definition: ShadingFunctions.h:66
float real_t
Definition: DataTypes.h:167
shared_ptr< VTKOutput > createVTKOutput_DomainDecomposition(const BlockStorage &bs, const std::string &identifier=std::string("domain_decomposition"), const uint_t writeFrequency=1, const std::string &baseFolder=std::string("vtk_out"), const std::string &executionFolder=std::string("simulation_step"), const bool continuousNumbering=false, const bool binary=true, const bool littleEndian=true, const bool useMPIIO=true, const uint_t initialExecutionCount=0)
Definition: VTKOutput.h:523
void configure(const Config::BlockHandle &config, pe::cr::HCSITS &cr)
configures HardContactSemiImplicitTimesteppingSolvers with parameters from config file
Definition: HCSITS.impl.h:1860
void storeTimingPoolInSqliteDB(const string &dbFile, uint_t runId, const WcTimingPool &tp, const std::string &timingPoolName, const int busyTimeout)
Definition: SQLite.cpp:501
static id_t getStaticTypeID()
Returns unique type id of this type.
Definition: Plane.h:359
void storeTimingTreeInSqliteDB(const string &dbFile, uint_t runId, const WcTimingTree &tt, const std::string &timingTreeName, const int busyTimeout)
Definition: SQLite.cpp:510
shared_ptr< VTKOutput > createVTKOutput_PointData(const shared_ptr< PointDataSource > pds, const std::string &identifier=std::string("point_data"), const uint_t writeFrequency=1, const std::string &baseFolder=std::string("vtk_out"), const std::string &executionFolder=std::string("simulation_step"), const bool continuousNumbering=false, const bool binary=true, const bool littleEndian=true, const bool useMPIIO=true, const uint_t initialExecutionCount=0)
Definition: VTKOutput.h:612
GenericAABB< real_t > AABB
Definition: AABBFwd.h:33
static id_t getStaticTypeID()
Returns unique type id of this type.
Definition: Sphere.h:466
ShadingParameters defaultBodyTypeDependentShadingParams(const BodyID body)
Definition: ShadingFunctions.h:48