diff --git a/apps/showcases/RayleighBenardConvection/CMakeLists.txt b/apps/showcases/RayleighBenardConvection/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..cccb8326e04677ce82a2d3f1b60c8b2c908b2146
--- /dev/null
+++ b/apps/showcases/RayleighBenardConvection/CMakeLists.txt
@@ -0,0 +1,4 @@
+if ( WALBERLA_BUILD_WITH_CODEGEN )
+ add_subdirectory( Codegen )
+ add_subdirectory( Native )
+endif ()
diff --git a/apps/showcases/RayleighBenardConvection/Codegen/RayleighBenardConvection_codegen.py b/apps/showcases/RayleighBenardConvection/Codegen/RayleighBenardConvection_codegen.py
index da8bd3e236b894cd4bf6b8afebd3e4cff947c9ff..b187973f5c6e314211c44d699b7265487c9d7227 100644
--- a/apps/showcases/RayleighBenardConvection/Codegen/RayleighBenardConvection_codegen.py
+++ b/apps/showcases/RayleighBenardConvection/Codegen/RayleighBenardConvection_codegen.py
@@ -95,6 +95,7 @@ lbm_config_thermal = LBMConfig(
method=Method.SRT,
relaxation_rate=omega_thermal,
compressible=True,
+ zero_centered=False,
velocity_input=velocity_field,
equilibrium_order=1,
entropic=False,
diff --git a/apps/showcases/RayleighBenardConvection/Native/CMakeLists.txt b/apps/showcases/RayleighBenardConvection/Native/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..eabc2d9c1a8ee84b5ec11245a48e3b90362f8ea1
--- /dev/null
+++ b/apps/showcases/RayleighBenardConvection/Native/CMakeLists.txt
@@ -0,0 +1,15 @@
+waLBerla_link_files_to_builddir(*.prm)
+waLBerla_link_files_to_builddir(*.py)
+waLBerla_link_files_to_builddir(*.obj)
+
+if( WALBERLA_BUILD_WITH_CODEGEN )
+ walberla_generate_target_from_python( NAME RayleighBenardConvectionLatticeModelGeneration
+ FILE RayleighBenardConvectionLatticeModelGeneration.py
+ OUT_FILES RayleighBenardConvectionLatticeModel_fluid.cpp RayleighBenardConvectionLatticeModel_fluid.h
+ RayleighBenardConvectionLatticeModel_thermal.cpp RayleighBenardConvectionLatticeModel_thermal.h )
+
+ waLBerla_add_executable(NAME RayleighBenardConvectionNative
+ FILES RayleighBenardConvectionNative.cpp #InitializerFunctions.cpp #SQLProperties.cpp
+ DEPENDS blockforest boundary core domain_decomposition field lbm postprocessing python_coupling#sqlite
+ timeloop vtk RayleighBenardConvectionLatticeModelGeneration)
+endif()
diff --git a/apps/showcases/RayleighBenardConvection/Native/RayleighBenardConvectionLatticeModelGeneration.py b/apps/showcases/RayleighBenardConvection/Native/RayleighBenardConvectionLatticeModelGeneration.py
new file mode 100644
index 0000000000000000000000000000000000000000..047d05affc4ed6a046fbbee9e83da9ef625a901f
--- /dev/null
+++ b/apps/showcases/RayleighBenardConvection/Native/RayleighBenardConvectionLatticeModelGeneration.py
@@ -0,0 +1,73 @@
+import sympy as sp
+import pystencils as ps
+from lbmpy.creationfunctions import LBMConfig, LBMOptimisation, create_lb_collision_rule
+from lbmpy.enums import ForceModel, Method, Stencil
+from lbmpy.stencils import LBStencil
+
+from pystencils_walberla import CodeGeneration
+from lbmpy_walberla import generate_lattice_model
+
+
+with CodeGeneration() as ctx:
+ stencil_thermal = LBStencil(Stencil.D3Q7)
+ stencil_fluid = LBStencil(Stencil.D3Q19)
+
+ target = ps.Target.CPU
+ layout = "fzyx"
+ zero_centered = False
+ force_model = ForceModel.SIMPLE
+
+ omega_fluid = sp.Symbol("omega_fluid")
+ omega_thermal = sp.Symbol("omega_thermal")
+ gravity_LBM = sp.Symbol("gravity") #! gravity? gravitational acceleration? ...? what exactly do I use / need here?
+ rho = sp.Symbol("rho")
+
+ # optimizations to be used by the code generator
+ lbm_opt = LBMOptimisation(cse_global=True, field_layout=layout)
+
+ velocity_field = ps.fields(
+ f"vel_field({stencil_fluid.D}): [{stencil_fluid.D}D]", layout=layout
+ )
+ temperature_field = ps.fields(f"temperature_field: [{stencil_thermal.D}D]", layout=layout)
+
+ force = sp.Matrix([0, gravity_LBM * temperature_field(0) * rho, 0])
+
+ #------------------------------ Fluid LBM ------------------------------
+ # Fluid LBM
+ lbm_config_fluid = LBMConfig(
+ stencil=stencil_fluid,
+ method=Method.SRT,
+ relaxation_rate=omega_fluid,
+ compressible=False,
+ zero_centered=zero_centered,
+ output={"velocity": velocity_field},
+ force=force,
+ force_model=ForceModel.SIMPLE,
+ ) # Simple force model (automatisch w_i, c_i, 3)
+
+ #> currently optimizations turned off!
+ collision_rule_fluid = create_lb_collision_rule(lbm_config=lbm_config_fluid)#,
+ #lbm_optimisation=lbm_opt)
+
+ generate_lattice_model(ctx, "RayleighBenardConvectionLatticeModel_fluid", collision_rule_fluid, field_layout=layout)
+
+ #------------------------------ Thermal LBM ------------------------------
+ # Thermal LBM
+ lbm_config_thermal = LBMConfig(
+ stencil=stencil_thermal,
+ method=Method.SRT,
+ relaxation_rate=omega_thermal,
+ compressible=True,
+ zero_centered=False, #? think about that -> if "defualt" delta_rho not known
+ velocity_input=velocity_field,
+ equilibrium_order=1,
+ entropic=False,
+ output={"density": temperature_field},
+ )
+ #> currently optimizations turned off!
+ collision_rule_thermal = create_lb_collision_rule(lbm_config=lbm_config_thermal)#,
+ #lbm_optimisation=lbm_opt)
+
+ generate_lattice_model(ctx, "RayleighBenardConvectionLatticeModel_thermal", collision_rule_thermal, field_layout=layout)
+
+print("\t >>> finished code generation successfully <<<")
diff --git a/apps/showcases/RayleighBenardConvection/Native/RayleighBenardConvectionNative.cpp b/apps/showcases/RayleighBenardConvection/Native/RayleighBenardConvectionNative.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..a8eca2f17b1f2f72ff3011baafa18e28ba8e4218
--- /dev/null
+++ b/apps/showcases/RayleighBenardConvection/Native/RayleighBenardConvectionNative.cpp
@@ -0,0 +1,378 @@
+//======================================================================================================================
+//
+// 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 RayleighBenardConvection.cpp
+//! \author Jonas Plewinski <jonas.plewinski@fau.de>
+//
+//======================================================================================================================
+
+#include "blockforest/Initialization.h"
+// #include "blockforest/communication/UniformBufferedScheme.h"
+#include "core/Environment.h"
+#include "core/logging/Initialization.h"
+#include "core/math/Constants.h"
+#include "core/timing/RemainingTimeLogger.h"
+
+#include "field/AddToStorage.h"
+
+#include "lbm/blockforest/communication/SimpleCommunication.h"
+// #include "field/Gather.h"
+#include "field/FlagField.h"
+#include "field/vtk/VTKWriter.h"
+
+#include "geometry/InitBoundaryHandling.h"
+#include "geometry/mesh/TriangleMeshIO.h"
+
+#include "lbm/PerformanceEvaluation.h"
+#include "lbm/boundary/all.h"
+#include "lbm/field/AddToStorage.h"
+#include "lbm/vtk/Density.h"
+
+#include "python_coupling/CreateConfig.h"
+// #include "python_coupling/PythonCallback.h"
+
+#include "timeloop/SweepTimeloop.h"
+
+//#include "InitializerFunctions.h"
+#include "RayleighBenardConvectionLatticeModel_fluid.h"
+#include "RayleighBenardConvectionLatticeModel_thermal.h"
+
+using namespace walberla;
+
+///////////////////////
+/// Typedef Aliases ///
+///////////////////////
+using ScalarField_T = GhostLayerField< real_t, 1 >;
+using VectorField_T = GhostLayerField< Vector3< real_t >, 1 >;
+using VectorFieldFlattened_T = GhostLayerField< real_t, 3 >;
+
+//using CollisionModel_T = lbm::collision_model::SRT;
+
+//--------------------------- fluid -----------------------------
+using LatticeModelFluid_T = lbm::RayleighBenardConvectionLatticeModel_fluid;
+using FluidStencil_T = LatticeModelFluid_T::Stencil;
+using PdfFieldFluid_T = lbm::PdfField< LatticeModelFluid_T >;
+using PdfFluidCommunication_T = blockforest::SimpleCommunication< FluidStencil_T >;
+
+//--------------------------- thermal ---------------------------
+using LatticeModelThermal_T = lbm::RayleighBenardConvectionLatticeModel_thermal;
+using ThermalStencil_T = LatticeModelThermal_T::Stencil;
+using PdfFieldThermal_T = lbm::PdfField< LatticeModelThermal_T >;
+using PdfThermalCommunication_T = blockforest::SimpleCommunication< ThermalStencil_T >;
+
+using flag_t = walberla::uint8_t;
+using FlagField_T = FlagField< flag_t >;
+
+// function describing the initialization profile (in global coordinates)
+real_t initializationProfile(real_t x, real_t amplitude, real_t offset, real_t wavelength)
+{
+ return amplitude * std::cos(x / wavelength * real_c(2) * math::pi + math::pi) + offset;
+}
+
+template< typename LatticeModel_T >
+void initTemperatureFieldTest(const shared_ptr< StructuredBlockStorage >& blocks,
+ BlockDataID temperature_field_ID,
+ BlockDataID pdf_field_ID,
+ real_t amplitude, Vector3< uint_t > domainSize, real_t temperatureRange)
+{
+ for (auto& block : *blocks)
+ {
+ auto temperatureField = block.getData< ScalarField_T >(temperature_field_ID);
+ auto pdfField = block.getData< lbm::PdfField< LatticeModel_T > >(pdf_field_ID);
+
+ WALBERLA_FOR_ALL_CELLS_INCLUDING_GHOST_LAYER_XYZ(temperatureField, {
+ //WALBERLA_FOR_ALL_CELLS_INCLUDING_GHOST_LAYER_XYZ(pdfField, {
+ // cell in block-local coordinates
+ Cell globalCell;
+ blocks->transformBlockLocalToGlobalCell(globalCell, block, Cell(x, y, z));
+
+ uint_t sampleSize = 100;
+ real_t stepSize = real_c(1.) / real_c(sampleSize);
+ real_t offset = real_c(domainSize[1] / 2);
+ real_t wavelength = real_c(domainSize[0]);
+
+ //--
+ Vector3 vel(real_c(0));
+ //--
+ for (uint_t xSample = uint_c(0); xSample <= sampleSize; ++xSample)
+ {
+ // value of the sine-function
+ const real_t functionValue = initializationProfile(real_c(globalCell[0]) + real_c(xSample) * stepSize, amplitude, offset, wavelength);
+ for (uint_t ySample = uint_c(0); ySample <= sampleSize; ++ySample)
+ {
+ const real_t yPoint = real_c(globalCell[1]) + real_c(ySample) * stepSize;
+ //temperatureField->get(x, y, z) = (yPoint < functionValue) ? temperatureRange/50 : -temperatureRange/50;
+ auto temp = (yPoint < functionValue) ? temperatureRange/50 : -temperatureRange/50;
+ pdfField->setDensityAndVelocity(x, y, z, vel, temp);
+ }
+ }
+ }) // WALBERLA_FOR_ALL_CELLS
+ }
+}
+
+/////////////////////
+/// Main Function ///
+/////////////////////
+int main(int argc, char** argv)
+{
+ mpi::Environment Env(argc, argv);
+ // exportDataStructuresToPython(); //> what does that do?
+
+ for (auto cfg = python_coupling::configBegin(argc, argv); cfg != python_coupling::configEnd(); ++cfg)
+ {
+ WALBERLA_MPI_WORLD_BARRIER()
+
+ auto config = *cfg;
+ logging::configureLogging(config);
+ // WALBERLA_LOG_DEVEL_VAR_ON_ROOT(*config)
+
+ ///////////////////////////
+ // ADD DOMAIN PARAMETERS //
+ ///////////////////////////
+
+ uint_t nrOfProcesses = uint_c(MPIManager::instance()->numProcesses());
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(nrOfProcesses)
+
+ auto domainSetup = config->getOneBlock("DomainSetup");
+ Vector3< uint_t > cellsPerBlock = domainSetup.getParameter< Vector3< uint_t > >("cellsPerBlock");
+ std::vector< config::Config::Block* > configDomainSetupBlock;
+ config->getWritableGlobalBlock().getWritableBlocks("DomainSetup", configDomainSetupBlock, 1, 1);
+ Vector3< uint_t > blocksPerDimension;
+ Vector3< uint_t > cellsPerBlockDummy;
+ blockforest::calculateCellDistribution(cellsPerBlock, nrOfProcesses, blocksPerDimension, cellsPerBlockDummy);
+
+ Vector3< uint_t > domainSize;
+ domainSize[0] = blocksPerDimension[0] * cellsPerBlock[0];
+ domainSize[1] = blocksPerDimension[1] * cellsPerBlock[1];
+ domainSize[2] = blocksPerDimension[2] * cellsPerBlock[2];
+ // WALBERLA_LOG_INFO_ON_ROOT("cellsPerBlock = (" << cellsPerBlock[0] << ", " << cellsPerBlock[1] << ", " <<
+ // cellsPerBlock[2] << ")") WALBERLA_LOG_INFO_ON_ROOT("blocksPerDimension = (" << blocksPerDimension[0] << ", " <<
+ // blocksPerDimension[1] << ", " << blocksPerDimension[2] << ")") WALBERLA_LOG_INFO_ON_ROOT("domain size = (" <<
+ // domainSize[0] << ", " << domainSize[1] << ", " << domainSize[2] << ")")
+
+ std::string tmp = "< " + std::to_string(blocksPerDimension[0]) + ", " + std::to_string(blocksPerDimension[1]) +
+ ", " + std::to_string(blocksPerDimension[2]) + " >";
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(tmp)
+ configDomainSetupBlock[0]->setOrAddParameter("blocks", tmp);
+
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(*config)
+
+ const shared_ptr< StructuredBlockForest > blockForest = blockforest::createUniformBlockGridFromConfig(config);
+ /*WALBERLA_LOG_DEVEL_VAR_ON_ROOT(blocks->getXSize())
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(blocks->getYSize())
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(blocks->getZSize())
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(blocks->getRootBlockXSize())
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(blocks->getRootBlockYSize())
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(blocks->getRootBlockZSize())*/
+
+ ///////////////////////////////////////
+ // ADD GENERAL SIMULATION PARAMETERS //
+ ///////////////////////////////////////
+ auto parameters = config->getOneBlock("Parameters");
+ //const std::string timeStepStrategy = parameters.getParameter< std::string >("timeStepStrategy", "normal");
+ const uint_t timesteps = parameters.getParameter< uint_t >("timesteps", uint_c(50));
+ const real_t remainingTimeLoggerFrequency =
+ parameters.getParameter< real_t >("remainingTimeLoggerFrequency", real_c(3.0));
+ // const uint_t scenario = parameters.getParameter< uint_t >("scenario", uint_c(1));
+ // WALBERLA_LOG_DEVEL_VAR_ON_ROOT(timesteps)
+ // WALBERLA_LOG_DEVEL_VAR_ON_ROOT(remainingTimeLoggerFrequency)
+
+ /////////////////////////////
+ // ADD PHYSICAL PARAMETERS //
+ /////////////////////////////
+ auto physical_parameters = config->getOneBlock("PhysicalParameters");
+ const real_t omegaFluid = physical_parameters.getParameter< real_t >("omegaFluid", real_c(1.95));
+ const real_t omegaThermal = physical_parameters.getParameter< real_t >("omegaThermal");
+ const real_t temperatureHot = physical_parameters.getParameter< real_t >("temperatureHot", real_c(0.5));
+ const real_t temperatureCold = physical_parameters.getParameter< real_t >("temperatureCold", real_c(-0.5));
+ const real_t gravity = physical_parameters.getParameter< real_t >("gravitationalAcceleration");
+ // WALBERLA_LOG_DEVEL_VAR_ON_ROOT(omegaFluid)
+ // WALBERLA_LOG_DEVEL_VAR_ON_ROOT(omegaThermal)
+ // WALBERLA_LOG_DEVEL_VAR_ON_ROOT(temperatureHot)
+ // WALBERLA_LOG_DEVEL_VAR_ON_ROOT(temperatureCold)
+
+ ///////////////////////////////////
+ // ADD INITIALIZATION PARAMETERS //
+ ///////////////////////////////////
+ auto initialization_parameters = config->getOneBlock("InitializationParameters");
+ const real_t initAmplitude = initialization_parameters.getParameter< real_t >("initAmplitude");
+ const real_t initTemperatureRange = initialization_parameters.getParameter< real_t >("initTemperatureRange");
+ // WALBERLA_LOG_DEVEL_VAR_ON_ROOT(initAmplitude)
+ // WALBERLA_LOG_DEVEL_VAR_ON_ROOT(initTemperatureRange)
+
+ ////////////////////////
+ // ADD DATA TO BLOCKS //
+ ////////////////////////
+ const BlockDataID temperatureFieldId =
+ field::addToStorage< ScalarField_T >(blockForest, "Temperature field", real_c(0), field::fzyx);
+ const BlockDataID velocityFieldId = field::addToStorage< VectorFieldFlattened_T >(
+ blockForest, "Velocity field", real_c(0), field::fzyx, uint_c(1));
+ //const BlockDataID densityFieldID = field::addToStorage< ScalarField_T >(blockForest, "density", real_c(0), field::fzyx);
+
+ // create lattice model
+ LatticeModelFluid_T latticeModelFluid = //? currently force
+ LatticeModelFluid_T(temperatureFieldId, velocityFieldId, gravity, omegaFluid); //? gravity needs a sign????
+ LatticeModelThermal_T latticeModelThermal =
+ LatticeModelThermal_T(temperatureFieldId, velocityFieldId, omegaThermal);
+
+ // add pdf fields for fluid and thermal parts
+ const BlockDataID pdfFieldFluidId = lbm::addPdfFieldToStorage(blockForest, "PDF field fluid", latticeModelFluid,
+ Vector3(real_c(0), real_c(0), real_c(0)),
+ real_t(1), uint_t(1), field::fzyx);
+ const BlockDataID pdfFieldThermalId =
+ lbm::addPdfFieldToStorage(blockForest, "PDF field thermal", latticeModelThermal, Vector3(real_c(0), real_c(0), real_c(0)),
+ real_t(0), uint_t(1),field::fzyx);
+ //? what initial density makes sense here?
+ //initTemperatureField(blockForest, temperatureFieldId, initAmplitude, domainSize, initTemperatureRange);
+ initTemperatureFieldTest<LatticeModelThermal_T>(blockForest, temperatureFieldId, pdfFieldThermalId, initAmplitude, domainSize, initTemperatureRange);
+
+ PdfThermalCommunication_T pdfThermalCommunication(blockForest, pdfFieldThermalId);
+
+ ////////////////
+ // ADD SWEEPS //
+ ////////////////
+ auto lbmSweepFluid = LatticeModelFluid_T::Sweep(pdfFieldFluidId);
+ auto lbmSweepThermal = LatticeModelThermal_T::Sweep(pdfFieldThermalId);
+
+ ///////////////////////
+ // BOUNDARY HANDLING //
+ ///////////////////////
+ const FlagUID fluidFlagUID("Fluid");
+ const FlagUID thermalFlagUID("Thermal");
+
+ // Boundaries Fluid
+ typedef lbm::DefaultBoundaryHandlingFactory< LatticeModelFluid_T, FlagField_T > BHFactoryFluid;
+
+ BlockDataID flagFieldFluidId = field::addFlagFieldToStorage< FlagField_T >(blockForest, "fluid flag field");
+
+ BlockDataID fluidBoundaryHandlingId = BHFactoryFluid::addBoundaryHandlingToStorage(
+ blockForest, "Fluid Boundary Handling", flagFieldFluidId, pdfFieldFluidId, fluidFlagUID,
+ Vector3< real_t >(real_c(0), real_c(0), real_c(0)), Vector3< real_t >(real_c(0), real_c(0), real_c(0)),
+ real_c(1), real_c(1));
+
+ auto fluidBoundariesConfig = config->getBlock("Boundaries_Fluid");
+ geometry::initBoundaryHandling< BHFactoryFluid::BoundaryHandling >(*blockForest, fluidBoundaryHandlingId,
+ fluidBoundariesConfig);
+ geometry::setNonBoundaryCellsToDomain< BHFactoryFluid::BoundaryHandling >(*blockForest, fluidBoundaryHandlingId);
+
+ // Boundaries Thermal
+ typedef lbm::DefaultBoundaryHandlingFactory< LatticeModelThermal_T, FlagField_T > BHFactoryThermal;
+
+ BlockDataID flagFieldThermalId = field::addFlagFieldToStorage< FlagField_T >(blockForest, "Thermal flag field");
+
+ BlockDataID thermalBoundaryHandlingId = BHFactoryThermal::addBoundaryHandlingToStorage(
+ blockForest, "Thermal Boundary Handling", flagFieldThermalId, pdfFieldThermalId,
+ thermalFlagUID, //? fluidFlagUID, does that work like that?
+ Vector3< real_t >(real_c(0), real_c(0), real_c(0)), Vector3< real_t >(real_c(0), real_c(0), real_c(0)),
+ temperatureCold, temperatureHot);
+
+ auto thermalBoundariesConfig = config->getBlock("Boundaries_Thermal");
+ geometry::initBoundaryHandling< BHFactoryThermal::BoundaryHandling >(*blockForest, thermalBoundaryHandlingId,
+ thermalBoundariesConfig);
+ geometry::setNonBoundaryCellsToDomain< BHFactoryThermal::BoundaryHandling >(*blockForest,
+ thermalBoundaryHandlingId);
+
+ ///////////////
+ // TIME LOOP //
+ ///////////////
+ auto benchmark_parameters = config->getOneBlock("BenchmarkParameters");
+ std::string scaling_type = benchmark_parameters.getParameter< std::string >("scalingType");
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(scaling_type)
+ WALBERLA_LOG_DEVEL_VAR_ON_ROOT(nrOfProcesses)
+ WALBERLA_LOG_INFO_ON_ROOT("#blocks = ("
+ << blockForest->getXSize() << ", " << blockForest->getYSize() << ", "
+ << blockForest->getZSize() << ") | total #blocks = "
+ << blockForest->getXSize() * blockForest->getYSize() * blockForest->getZSize())
+ WALBERLA_LOG_INFO_ON_ROOT("block size = (" << blockForest->getRootBlockXSize() << ", "
+ << blockForest->getRootBlockYSize() << ", "
+ << blockForest->getRootBlockZSize() << ")")
+ WALBERLA_LOG_INFO_ON_ROOT("domain size = (" << domainSize[0] << ", " << domainSize[1] << ", " << domainSize[2]
+ << ")")
+
+ bool weak_scaling = benchmark_parameters.getParameter< bool >("weakScaling");
+
+ SweepTimeloop timeloop(blockForest->getBlockStorage(), timesteps);
+ if (!weak_scaling)
+ {
+ timeloop.add() << Sweep(BHFactoryThermal::BoundaryHandling::getBlockSweep(thermalBoundaryHandlingId),
+ "Thermal boundary conditions")
+ << AfterFunction(PdfThermalCommunication_T(blockForest, pdfFieldThermalId),
+ "Communication of thermal PDFs");
+ timeloop.add() << Sweep(lbmSweepThermal, "Thermal LB Step");
+
+ timeloop.add() << BeforeFunction(PdfFluidCommunication_T(blockForest, pdfFieldFluidId),
+ "Communication of fluid PDFs")
+ << Sweep(BHFactoryFluid::BoundaryHandling::getBlockSweep(fluidBoundaryHandlingId),
+ "Fluid NoSlip boundary conditions");
+ timeloop.add() << Sweep(lbmSweepFluid, "Fluid LB Step");
+
+ // remaining time logger
+ timeloop.addFuncAfterTimeStep(
+ RemainingTimeLogger(timeloop.getNrOfTimeSteps(), remainingTimeLoggerFrequency),
+ "remaining time logger");
+
+ // write VTK files
+ uint_t vtkWriteFrequency = parameters.getParameter< uint_t >("vtkWriteFrequency", 0);
+ if (vtkWriteFrequency > 0)
+ {
+ auto vtkOutput = vtk::createVTKOutput_BlockData(*blockForest, "vtk", vtkWriteFrequency, 1, false, "vtk_out",
+ "simulation_step", false, true, true, false, 0);
+
+ vtkOutput->addCellDataWriter( make_shared< lbm::DensityVTKWriter < LatticeModelFluid_T, float > >( pdfFieldFluidId, "DensityFromPDF" ) );
+ vtkOutput->addCellDataWriter( make_shared< lbm::DensityVTKWriter < LatticeModelThermal_T, float > >( pdfFieldThermalId, "TemperatureFromPDF" ) );
+
+ // add velocity field as VTK output
+ auto velWriter = make_shared< field::VTKWriter< VectorFieldFlattened_T > >(velocityFieldId, "Velocity");
+ vtkOutput->addCellDataWriter(velWriter);
+
+ // add temperature field as VTK output
+ auto tempWriter = make_shared< field::VTKWriter< ScalarField_T > >(temperatureFieldId, "Temperature");
+ vtkOutput->addCellDataWriter(tempWriter);
+
+ // add thermal flag field as VTK output
+ auto thermalFlagWriter =
+ make_shared< field::VTKWriter< FlagField_T > >(flagFieldThermalId, "FlagFieldThermal");
+ vtkOutput->addCellDataWriter(thermalFlagWriter);
+
+ // add fluid flag field as VTK output
+ auto fluidFlagWriter = make_shared< field::VTKWriter< FlagField_T > >(flagFieldFluidId, "FlagFieldFluid");
+ vtkOutput->addCellDataWriter(fluidFlagWriter);
+
+ timeloop.addFuncBeforeTimeStep(writeFiles(vtkOutput), "VTK Output");
+ }
+
+ // Performance evaluation
+ lbm::PerformanceEvaluation< FlagField_T > performance(blockForest, flagFieldFluidId, fluidFlagUID);
+ WcTimingPool timeloopTiming;
+ WcTimer simTimer;
+
+ WALBERLA_LOG_INFO_ON_ROOT("Starting simulation with " << timesteps << " time steps")
+ WALBERLA_MPI_WORLD_BARRIER()
+ simTimer.start();
+ timeloop.run(timeloopTiming);
+ WALBERLA_MPI_WORLD_BARRIER()
+ simTimer.end();
+
+ auto time = real_c(simTimer.max());
+ WALBERLA_MPI_SECTION() { reduceInplace(time, walberla::mpi::MAX); }
+ performance.logResultOnRoot(timesteps, time);
+
+ const auto reducedTimeloopTiming = timeloopTiming.getReduced();
+ WALBERLA_LOG_RESULT_ON_ROOT("Time loop timing:\n" << *reducedTimeloopTiming)
+ WALBERLA_LOG_INFO_ON_ROOT("Simulation done!")
+ }
+ }
+ return EXIT_SUCCESS;
+}
\ No newline at end of file
diff --git a/apps/showcases/RayleighBenardConvection/Native/RayleighBenardConvection_config.py b/apps/showcases/RayleighBenardConvection/Native/RayleighBenardConvection_config.py
new file mode 100755
index 0000000000000000000000000000000000000000..30cb3c5040919867802229cd41c5d5f2d72a0504
--- /dev/null
+++ b/apps/showcases/RayleighBenardConvection/Native/RayleighBenardConvection_config.py
@@ -0,0 +1,134 @@
+import numpy as np
+import waLBerla as wlb
+import math
+
+
+class Scenario:
+ def __init__(self, weak_scaling=False, scaling_type=None, cells=(32, 32, 1), timesteps=100):
+ #> Domain Parameters
+ self.domain_size = (32, 32, 1)
+ self.blocks = (4, 1, 1)
+ self.periodic = (1, 0, 1)
+ self.cells = (self.domain_size[0] // self.blocks[0], self.domain_size[1] // self.blocks[1], self.domain_size[2] // self.blocks[2])
+ self.cells = cells
+ #print(f"self.cells = {self.cells}")
+ #> Standard Parameters
+ self.timesteps = timesteps
+ self.vtk_write_frequency = 200
+ self.scenario = 1
+ #> Physical Parameters
+ #! Prandtl
+ self.Prandtl = 1
+ #! Rayleigh
+ self.Rayleigh = 50000
+ #! omega_fluid
+ self.omega_fluid = 1.95
+ #! omega_thermal
+ self.length_x_SI = 2
+ self.length_conversion = self.length_x_SI / self.domain_size[0]
+ self.viscosity_LBM = 1./3 * (1. / self.omega_fluid - 1./2)
+ self.time_conversion = self.viscosity_LBM * self.length_conversion * self.length_conversion / \
+ np.sqrt(self.Prandtl / self.Rayleigh)
+ self.thermal_diffusivity_LBM = np.sqrt(1. / (self.Prandtl * self.Rayleigh)) * self.time_conversion / \
+ (self.length_conversion * self.length_conversion)
+ self.omega_thermal = 1. / (3 * self.thermal_diffusivity_LBM + 1./2)
+ #print(f"omega_fluid = {self.omega_fluid} | omega_thermal = {self.omega_thermal}")
+ #! temperature_hot
+ self.temperature_hot = 0.5
+ #! temperature_cold
+ self.temperature_cold = - self.temperature_hot
+ #print(f"temperature_hot = {self.temperature_hot} | temperature_cold = {self.temperature_cold}")
+ #! gravity_LBM
+ self.gravity_SI = 9.81
+ self.gravity_conversion = self.length_conversion / (self.time_conversion * self.time_conversion) #! look
+ self.gravity_LBM = self.gravity_SI / self.gravity_conversion
+ #> Initialization Parameters
+ self.delta_temperature = abs(self.temperature_hot) + abs(self.temperature_cold)
+ self.init_amplitude = self.delta_temperature / 20 * self.domain_size[0]
+ self.init_temperature_range = 2 * self.delta_temperature / self.domain_size[0]
+ #print(f"init_amplitude = {self.init_amplitude}")
+ #print(f"init_temperature_range = {self.init_temperature_range}")
+
+ #> Benchmark Parameters
+ self.weakScaling = weak_scaling # True
+ self.scalingType = scaling_type #"fluid", "thermal", "rbc"
+ self.benchmarkingIterations = 5
+ self.warmupSteps = 10
+
+ #?self.viscosity_SI = 1.516e-5 #? look
+ #?self.thermal_diffusivity_SI = 2.074e-5 #? look
+
+ #?self.alpha = 3.43e-3 #? look
+
+ @wlb.member_callback
+ def config(self):
+ return {
+ 'DomainSetup': {
+ #'blocks': self.blocks,
+ #'domainSize': self.domain_size,
+ #'cellsPerBlock': self.cells,
+ 'cellsPerBlock': self.cells,
+ 'blocks': self.blocks,
+ 'periodic': self.periodic,
+ 'dx': 1.0,
+ },
+ 'Parameters': {
+ 'timesteps': self.timesteps,
+ 'vtkWriteFrequency': self.vtk_write_frequency,
+ 'remainingTimeLoggerFrequency': 10.0,
+ 'scenario': self.scenario,
+ },
+ 'PhysicalParameters': { #todo check gravity!
+ 'omegaFluid': self.omega_fluid,
+ 'omegaThermal': self.omega_thermal,
+ 'temperatureHot': self.temperature_hot,
+ 'temperatureCold': self.temperature_cold,
+ 'gravitationalAcceleration': self.gravity_LBM, #? is this the right gravity I use here?
+ #? What do I have to use here?
+ #? - gravitaional acceleration
+ #? - gravitaional density
+ #? - ...
+ },
+ 'InitializationParameters': {
+ 'initAmplitude': self.init_amplitude,
+ 'initTemperatureRange': self.init_temperature_range,
+ },
+ 'BenchmarkParameters': {
+ 'weakScaling': self.weakScaling,
+ 'scalingType': self.scalingType, #"fluid", "thermal", "rbc"
+ 'benchmarkingIterations': self.benchmarkingIterations,
+ 'warmupSteps': self.warmupSteps,
+ },
+ 'Boundaries_Fluid': {
+ 'Border': [
+ {'direction': 'N', 'walldistance': -1, 'flag': 'NoSlip'},
+ {'direction': 'S', 'walldistance': -1, 'flag': 'NoSlip'},
+ ]
+ },
+ 'Boundaries_Thermal': {
+ 'Border': [
+ {'direction': 'N', 'walldistance': -1, 'flag': 'Pressure0'},
+ {'direction': 'S', 'walldistance': -1, 'flag': 'Pressure1'},
+ ]
+ },
+ }
+
+
+def runSimulation():
+ scenarios = wlb.ScenarioManager()
+ scenarios.add(Scenario(timesteps=10000))
+
+
+def weakScalingBenchmark():
+ scenarios = wlb.ScenarioManager()
+ block_sizes = [(i, i, i) for i in (8, 16, 32, 64, 128)]
+ timestepsPerBlockSize = [1024, 512, 256, 128, 64]
+
+ for scaling_type in ['rbc', 'fluid', 'thermal']:
+ for block_size, timesteps in [(block_sizes[i], timestepsPerBlockSize[i]) for i in range(len(timestepsPerBlockSize))]:
+ scenario = Scenario(weak_scaling=True, scaling_type=scaling_type, cells=block_size, timesteps=timesteps)
+ scenarios.add(scenario)
+
+
+runSimulation()
+#weakScalingBenchmark()