diff --git a/tests/BasicLbmScenarios/LbmAlgorithms.py b/tests/BasicLbmScenarios/LbmAlgorithms.py
index 66c5add1a8c6b7eb797034e7cd66be5c2b4bf3a0..fab81e8bc4c242282668dccda60b883dedd33180 100644
--- a/tests/BasicLbmScenarios/LbmAlgorithms.py
+++ b/tests/BasicLbmScenarios/LbmAlgorithms.py
@@ -2,7 +2,14 @@ import sympy as sp
from pystencilssfg import SourceFileGenerator
from pystencils import fields, Field
-from lbmpy import LBStencil, Stencil, LBMConfig, LBMOptimisation, create_lb_update_rule
+from lbmpy import (
+ LBStencil,
+ Stencil,
+ LBMConfig,
+ LBMOptimisation,
+ create_lb_update_rule,
+ ForceModel,
+)
from lbmpy.macroscopic_value_kernels import macroscopic_values_setter
from walberla.codegen import Sweep
@@ -19,6 +26,7 @@ with SourceFileGenerator() as sfg:
f, f_tmp, rho, u = fields(
f"f({q}), f_tmp({q}), rho, u({d}): double[{d}D]", layout="fzyx"
) # type: ignore
+ force = sp.symbols(f"F_:{d}")
stencil_name = stencil.name
sfg.include(f"stencil/{stencil_name}.h")
@@ -27,6 +35,8 @@ with SourceFileGenerator() as sfg:
lbm_config = LBMConfig(
stencil=stencil,
output={"density": rho, "velocity": u},
+ force=force,
+ force_model=ForceModel.HE,
)
lbm_opt = LBMOptimisation(
symbolic_field=f,
@@ -43,13 +53,22 @@ with SourceFileGenerator() as sfg:
lb_init = macroscopic_values_setter(
lb_update.method,
- density=sp.Symbol("density"),
+ density=rho.center,
velocity=u.center_vector,
pdfs=f,
)
- lb_init_sweep = Sweep("LbInit", lb_init)
+ lb_init_sweep = Sweep("LbInitFromFields", lb_init)
sfg.generate(lb_init_sweep)
+ lb_init_constant = macroscopic_values_setter(
+ lb_update.method,
+ density=sp.Symbol("density"),
+ velocity=sp.symbols(f"velocity_:{d}"),
+ pdfs=f,
+ )
+ lb_init_constant_sweep = Sweep("LbInitConstant", lb_init_constant)
+ sfg.generate(lb_init_constant_sweep)
+
with sfg.namespace("bc_sparse"):
irreg_freeslip = FreeSlip(
"FreeSlipIrregular", lb_update.method, f, wall_normal=FreeSlip.IRREGULAR
diff --git a/tests/BasicLbmScenarios/SimDomain.hpp b/tests/BasicLbmScenarios/SimDomain.hpp
index b3693e5c86ec419d5b626c144ff12d1810469f23..f26337473cf760d1014ea17dbdb59df1fdf199ab 100644
--- a/tests/BasicLbmScenarios/SimDomain.hpp
+++ b/tests/BasicLbmScenarios/SimDomain.hpp
@@ -5,6 +5,7 @@
#include "core/all.h"
#include "field/all.h"
+#include "field/communication/StencilRestrictedPackInfo.h"
#include <memory>
@@ -44,6 +45,49 @@ struct SimDomain
const BlockDataID uId;
} gpuFields;
+#else
+
+ void initFromFields(const Vector3< real_t > force)
+ {
+ gen::bulk::LbInitFromFields initialize{ cpuFields.pdfsId, cpuFields.rhoId, cpuFields.uId, force };
+
+ for (auto& b : *blocks)
+ {
+ initialize(&b);
+ }
+ }
+
+ void initConstant(const real_t rho, const Vector3< real_t > u, const Vector3< real_t > force)
+ {
+ gen::bulk::LbInitConstant initialize{ cpuFields.pdfsId, force, rho, u };
+
+ for (auto& b : *blocks)
+ {
+ initialize(&b);
+ }
+ }
+
+ gen::bulk::LbStreamCollide streamCollideSweep(const real_t omega, const Vector3< real_t > force)
+ {
+ return { cpuFields.pdfsId, cpuFields.rhoId, cpuFields.uId, force, omega };
+ }
+
+ void sync() { /* NOP */ }
+
+ void forAllBlocks(std::function< void(IBlock&) > func)
+ {
+ for (auto& block : *blocks)
+ func(block);
+ }
+
+ void forAllCells(std::function< void(Cell) > func)
+ {
+ for (uint_t z = 0; z < blocks->getNumberOfZCellsPerBlock(); ++z)
+ for (uint_t y = 0; y < blocks->getNumberOfYCellsPerBlock(); ++y)
+ for (uint_t x = 0; x < blocks->getNumberOfXCellsPerBlock(); ++x)
+ func({ x, y, z });
+ }
+
#endif
};
diff --git a/tests/BasicLbmScenarios/TestBasicLbmScenarios.cpp b/tests/BasicLbmScenarios/TestBasicLbmScenarios.cpp
index 04413de31eb8a6f0f271842f919e64860f4e095a..55fdf16b0354abec90c50a68a498572785b935a2 100644
--- a/tests/BasicLbmScenarios/TestBasicLbmScenarios.cpp
+++ b/tests/BasicLbmScenarios/TestBasicLbmScenarios.cpp
@@ -3,10 +3,6 @@
#include "core/all.h"
-#include "field/all.h"
-#include "field/communication/PackInfo.h"
-#include "field/communication/StencilRestrictedPackInfo.h"
-
#include "geometry/all.h"
#include "vtk/all.h"
@@ -19,6 +15,39 @@ namespace BasicLbmScenarios
{
using namespace walberla;
+void fullyPeriodic(Environment& env)
+{
+ SimDomain dom{ SimDomainBuilder{
+ .blocks = { 1, 1, 1 }, .cellsPerBlock = { 32, 32, 32 }, .periodic = { true, true, true } }
+ .build() };
+
+ const Vector3< real_t > force { 0.005, 0., 0. };
+
+ dom.initConstant(1.0, Vector3< real_t > { 0.0 }, force);
+
+ auto streamCollide = dom.streamCollideSweep( 1.0, force );
+
+ for(uint_t t = 0; t < 10; ++t){
+ dom.comm();
+
+ dom.forAllBlocks([&](IBlock & b){
+ streamCollide(&b);
+ });
+
+ dom.sync();
+
+ dom.forAllBlocks([&](auto & block){
+ const VectorField_T & velField = *block.template getData< VectorField_T >(dom.cpuFields.uId);
+
+ dom.forAllCells([&](Cell c){
+ real_t expected { real_c(t) * force[0] };
+ real_t actual { velField.get(c, 0) };
+ WALBERLA_CHECK_FLOAT_EQUAL(expected, actual);
+ });
+ });
+ }
+}
+
/**
* Pipe flow with circular cross-section and free-slip walls.
* The pipe flow is initialized with a constant velocity in x-direction.
@@ -45,12 +74,14 @@ void freeSlipPipe(Environment& env)
const real_t pipeRadius{ 14.0 };
const Vector3< real_t > pipeAnchor{ 0.0, 16.0, 16.0 };
const real_t maxVelocity{ 0.02 };
+ const Vector3< real_t > force{ 0., 0., 0. };
for (auto& block : *dom.blocks)
{
FlagField_T& flagField = *block.getData< FlagField_T >(dom.cpuFields.flagFieldId);
const uint8_t freeSlipFlag{ flagField.getOrRegisterFlag(freeSlipFlagUID) };
+ ScalarField_T& densField = *block.getData< ScalarField_T >(dom.cpuFields.rhoId);
VectorField_T& velField = *block.getData< VectorField_T >(dom.cpuFields.uId);
for (Cell c : allCellsWithGl)
@@ -69,6 +100,7 @@ void freeSlipPipe(Environment& env)
}
else { initVelocity = Vector3< real_t >{ maxVelocity, 0., 0. }; }
+ densField.get(c, 0) = 1.0;
velField.get(c, 0) = initVelocity[0];
velField.get(c, 1) = initVelocity[1];
velField.get(c, 2) = initVelocity[2];
@@ -80,19 +112,14 @@ void freeSlipPipe(Environment& env)
auto flagsOutput = field::createVTKOutput< FlagField_T >(dom.cpuFields.flagFieldId, *dom.blocks, "flags");
flagsOutput();
- gen::bulk::LbInit initialize{ dom.cpuFields.pdfsId, dom.cpuFields.uId, 1.0 };
-
- for (auto& b : *dom.blocks)
- {
- initialize(&b);
- }
+ dom.initFromFields(force);
gen::bc_sparse::FreeSlipIrregular freeSlip{
gen::bc_sparse::FreeSlipIrregularFactory(dom.blocks, dom.cpuFields.pdfsId)
.fromFlagField< FlagField_T >(dom.cpuFields.flagFieldId, freeSlipFlagUID, fluidFlagUid)
};
- gen::bulk::LbStreamCollide streamCollide{ dom.cpuFields.pdfsId, dom.cpuFields.rhoId, dom.cpuFields.uId, 1.0 };
+ auto streamCollide = dom.streamCollideSweep(1.0, force);
auto velOutput = field::createVTKOutput< VectorField_T >(dom.cpuFields.uId, *dom.blocks, "vel");
@@ -128,5 +155,6 @@ void freeSlipPipe(Environment& env)
int main(int argc, char** argv)
{
walberla::Environment env{ argc, argv };
+ BasicLbmScenarios::fullyPeriodic(env);
BasicLbmScenarios::freeSlipPipe(env);
}