diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
index 7340d091065e0a1de1cf7510914f0021827c6a93..d501e83f223767a186a1fac6c5beccf5e144a4d0 100644
--- a/.gitlab-ci.yml
+++ b/.gitlab-ci.yml
@@ -4,6 +4,23 @@ stages:
- Deploy
+testsuite-clang18:
+ image: i10git.cs.fau.de:5005/walberla/buildenvs/clang-18:latest
+ tags:
+ - docker
+ stage: "Test"
+ needs: []
+ variables:
+ CXX: clang++
+ CC: clang
+ before_script:
+ - cd tests
+ - cmake -G Ninja -S . -B build
+ - cd build
+ - cmake --build . --target SfgTests
+ script:
+ - ctest
+
build-examples:
image: i10git.cs.fau.de:5005/walberla/buildenvs/clang-18:latest
tags:
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 16854db8f42545bcf0434a71bf60939b9eb91b2a..f9e705726123ca02a17edeab80e8f9ebe52c6ecb 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -10,7 +10,10 @@ include( PrepareSFG )
add_library( sfg_walberla INTERFACE )
target_sources( sfg_walberla
- INTERFACE include/sfg/SparseIteration.hpp include/sfg/GenericHbbBoundary.hpp
+ INTERFACE
+ include/sfg/SparseIteration.hpp
+ include/sfg/GenericHbbBoundary.hpp
+ include/sfg/IrregularFreeSlip.hpp
)
target_include_directories( sfg_walberla INTERFACE ${CMAKE_CURRENT_SOURCE_DIR}/include )
diff --git a/include/sfg/IrregularFreeSlip.hpp b/include/sfg/IrregularFreeSlip.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..c887b8e7bd55d437b4a0d83160c47b7fe66682a7
--- /dev/null
+++ b/include/sfg/IrregularFreeSlip.hpp
@@ -0,0 +1,185 @@
+#pragma once
+
+#include "core/DataTypes.h"
+#include "core/cell/Cell.h"
+#include "domain_decomposition/BlockDataID.h"
+#include "domain_decomposition/IBlock.h"
+#include "field/GhostLayerField.h"
+#include "field/FlagField.h"
+#include "sfg/SparseIteration.hpp"
+#include "stencil/Directions.h"
+#include <tuple>
+
+#include "SparseIteration.hpp"
+
+namespace walberla::sfg
+{
+ struct IrregularFreeSlipLinkInfo
+ {
+ int32_t x;
+ int32_t y;
+ int32_t z;
+ int32_t dir;
+ int32_t source_offset_x;
+ int32_t source_offset_y;
+ int32_t source_offset_z;
+ int32_t source_dir;
+ IrregularFreeSlipLinkInfo(
+ walberla::Cell fluidCell, walberla::stencil::Direction linkDir,
+ walberla::Cell sourceCell, walberla::stencil::Direction sourceDir)
+ : x{fluidCell.x()}, y{fluidCell.y()}, z{fluidCell.z()},
+ dir{int32_t(linkDir)}, source_offset_x{sourceCell.x() - fluidCell.x()},
+ source_offset_y{sourceCell.y() - fluidCell.y()}, source_offset_z{sourceCell.z() - fluidCell.z()},
+ source_dir{int32_t(sourceDir)} {}
+
+ bool operator==(const IrregularFreeSlipLinkInfo &other) const
+ {
+ return std::tie(x, y, z, dir, source_offset_x, source_offset_y, source_offset_z, source_dir) ==
+ std::tie(other.x, other.y, other.z, other.dir, other.source_offset_x,
+ other.source_offset_y, other.source_offset_z, other.source_dir);
+ }
+ };
+
+ namespace detail
+ {
+ template <typename Stencil, typename FlagField_T>
+ class FreeSlipLinksFromFlagField
+ {
+ public:
+ using IndexList = SparseIndexList<IrregularFreeSlipLinkInfo>;
+ using flag_t = typename FlagField_T::flag_t;
+
+ FreeSlipLinksFromFlagField(
+ StructuredBlockForest &sbfs,
+ BlockDataID flagFieldID,
+ field::FlagUID boundaryFlagUID,
+ field::FlagUID fluidFlagUID) : sbfs_{sbfs}, flagFieldID_(flagFieldID), boundaryFlagUID_(boundaryFlagUID), fluidFlagUID_(fluidFlagUID) {}
+
+ IndexList collectLinks()
+ {
+ IndexList indexList{sbfs_};
+
+ for (auto &block : sbfs_)
+ {
+ FlagField_T *flagField = block.getData<FlagField_T>(flagFieldID_);
+ flag_t boundaryFlag{flagField->getFlag(boundaryFlagUID_)};
+ flag_t fluidFlag{flagField->getFlag(fluidFlagUID_)};
+
+ auto &idxVector = indexList.get(block);
+
+ for (auto it = flagField->beginXYZ(); it != flagField->end(); ++it)
+ {
+ Cell c{it.cell()};
+ if (!isFlagSet(it, fluidFlag))
+ {
+ continue;
+ }
+
+ for (auto dIt = Stencil::beginNoCenter(); dIt != Stencil::end(); ++dIt)
+ {
+ stencil::Direction dir{*dIt};
+ Cell neighbor{c + dir};
+ if (flagField->isFlagSet(neighbor, boundaryFlag))
+ {
+ idxVector.emplace_back(createLink(flagField, c, dir));
+ }
+ }
+ }
+ }
+
+ return indexList;
+ }
+
+ private:
+ IrregularFreeSlipLinkInfo createLink(FlagField_T *flagField, const Cell &fluidCell, const stencil::Direction dir)
+ {
+ const flag_t fluidFlag{flagField->getFlag(fluidFlagUID_)};
+ const Cell wallCell{fluidCell + dir};
+
+ // inverse direction of 'dir' as lattice vector
+
+ const cell_idx_t ix = stencil::cx[stencil::inverseDir[dir]];
+ const cell_idx_t iy = stencil::cy[stencil::inverseDir[dir]];
+ const cell_idx_t iz = stencil::cz[stencil::inverseDir[dir]];
+
+ stencil::Direction sourceDir = stencil::inverseDir[dir]; // compute reflected (mirrored) of inverse direction of 'dir'
+
+ cell_idx_t wnx = 0; // compute "normal" vector of free slip wall
+ cell_idx_t wny = 0;
+ cell_idx_t wnz = 0;
+
+ if (flagField->isFlagSet(wallCell.x() + ix, wallCell.y(), wallCell.z(), fluidFlag))
+ {
+ wnx = ix;
+ sourceDir = stencil::mirrorX[sourceDir];
+ }
+ if (flagField->isFlagSet(wallCell.x(), wallCell.y() + iy, wallCell.z(), fluidFlag))
+ {
+ wny = iy;
+ sourceDir = stencil::mirrorY[sourceDir];
+ }
+ if (flagField->isFlagSet(wallCell.x(), wallCell.y(), wallCell.z() + iz, fluidFlag))
+ {
+ wnz = iz;
+ sourceDir = stencil::mirrorZ[sourceDir];
+ }
+
+ // concave corner (neighbors are non-fluid)
+ if (wnx == 0 && wny == 0 && wnz == 0)
+ {
+ wnx = ix;
+ wny = iy;
+ wnz = iz;
+ sourceDir = dir;
+ }
+
+ const Cell sourceCell{
+ wallCell.x() + wnx,
+ wallCell.y() + wny,
+ wallCell.z() + wnz};
+
+ return {
+ fluidCell,
+ dir,
+ sourceCell,
+ sourceDir};
+ }
+
+ StructuredBlockForest &sbfs_;
+ const BlockDataID flagFieldID_;
+ const field::FlagUID boundaryFlagUID_;
+ const field::FlagUID fluidFlagUID_;
+ };
+ }
+
+ template <typename Impl>
+ class IrregularFreeSlipFactory
+ {
+ public:
+ IrregularFreeSlipFactory(const shared_ptr<StructuredBlockForest> blocks, BlockDataID pdfFieldID)
+ : blocks_{blocks}, pdfFieldID_{pdfFieldID} {}
+
+ template <typename FlagField_T>
+ auto fromFlagField(BlockDataID flagFieldID, field::FlagUID boundaryFlagUID, field::FlagUID fluidFlagUID)
+ {
+ detail::FreeSlipLinksFromFlagField<typename Impl::Stencil, FlagField_T> linksFromFlagField{*blocks_, flagFieldID, boundaryFlagUID, fluidFlagUID};
+ auto indexVector = linksFromFlagField.collectLinks();
+ return impl().irregularFromIndexVector(indexVector);
+ }
+
+ protected:
+ shared_ptr<StructuredBlockForest> blocks_;
+ BlockDataID pdfFieldID_;
+
+ private:
+ Impl &impl()
+ {
+ return *static_cast<Impl *>(this);
+ }
+
+ const Impl &impl() const
+ {
+ return *static_cast<const Impl *>(this);
+ }
+ };
+}
diff --git a/src/sfg_walberla/api.py b/src/sfg_walberla/api.py
index 0be30b8dc161a51ba386b036ec897fbc019f5b5c..207c433230b72df39c1927664e7878c010cb896f 100644
--- a/src/sfg_walberla/api.py
+++ b/src/sfg_walberla/api.py
@@ -19,9 +19,10 @@ from pystencilssfg.lang import (
Ref,
ExprLike,
cpptype,
- CppClass
+ CppClass,
)
from pystencilssfg.lang.types import CppTypeFactory, CppType
+from pystencilssfg.lang.cpp import std
real_t = PsCustomType("walberla::real_t")
@@ -84,14 +85,41 @@ class AABB(_PlainCppClass):
return Vector3(real_t).bind("{}.max()", self)
+class Cell(_PlainCppClass):
+ _type = cpptype("walberla::Cell", "core/cell/Cell.h")
+
+ def x(self) -> AugExpr:
+ return AugExpr.format("{}.x()", self)
+
+ def y(self) -> AugExpr:
+ return AugExpr.format("{}.y()", self)
+
+ def z(self) -> AugExpr:
+ return AugExpr.format("{}.z()", self)
+
+
class CellInterval(_PlainCppClass):
_type = cpptype("walberla::CellInterval", "core/cell/CellInterval.h")
+class Direction(_PlainCppClass):
+ _type = cpptype("walberla::stencil::Direction", "stencil/Directions.h")
+
+
class BlockDataID(_PlainCppClass):
_type = cpptype("walberla::BlockDataID", "domain_decomposition/BlockDataID.h")
+class IBlock(_PlainCppClass):
+ _type = cpptype("walberla::IBlock", "domain_decomposition/IBlock.h")
+
+ def getData(self, dtype: str | PsType, id: BlockDataID) -> AugExpr:
+ return AugExpr.format("{}.template getData< {} >({})", self, dtype, id)
+
+ def getAABB(self) -> AABB:
+ return AABB().bind("{}.getAABB()", self)
+
+
class IBlockPtr(_PlainCppClass):
_type = cpptype("walberla::IBlock *", "domain_decomposition/IBlock.h")
@@ -110,7 +138,9 @@ class SharedPtr(CppClass):
class StructuredBlockForest(AugExpr):
- typename = cpptype("walberla::StructuredBlockForest", "blockforest/StructuredBlockForest.h")
+ typename = cpptype(
+ "walberla::StructuredBlockForest", "blockforest/StructuredBlockForest.h"
+ )
def __init__(self, ref: bool = True, const: bool = False):
dtype = self.typename(const=const, ref=ref)
@@ -303,9 +333,29 @@ def glfield(field: Field, ci: str | AugExpr | None = None):
return GhostLayerFieldExtraction(field_ptr, ci)
+class SparseIndexList(AugExpr):
+ _template = cpptype(
+ "walberla::sfg::SparseIndexList< {IndexStruct} >", "sfg/SparseIteration.hpp"
+ )
+
+ def __init__(
+ self, idx_struct: PsStructType, const: bool = False, ref: bool = False
+ ):
+ self._idx_struct = idx_struct
+ dtype = self._template(IndexStruct=idx_struct, const=const, ref=ref)
+ super().__init__(dtype)
+
+ def get(self, block: IBlock) -> std.vector:
+ return std.vector(self._idx_struct, ref=True).bind("{}.get({})", self, block)
+
+ def bufferId(self) -> BlockDataID:
+ return BlockDataID().bind("{}.bufferId()", self)
+
+
class IndexListBufferPtr(SrcField):
_template = cpptype(
- "walberla::sfg::internal::IndexListBuffer< {IndexStruct} >", "sfg/SparseIteration.hpp"
+ "walberla::sfg::internal::IndexListBuffer< {IndexStruct} >",
+ "sfg/SparseIteration.hpp",
)
def __init__(self, idx_struct: PsStructType):
@@ -350,7 +400,7 @@ CellIdx = PsStructType(
(
("x", create_type("int64_t")),
("y", create_type("int64_t")),
- ("z", create_type("int64_t"))
+ ("z", create_type("int64_t")),
),
- "walberla::sfg::CellIdx"
+ "walberla::sfg::CellIdx",
)
diff --git a/src/sfg_walberla/boundaries/__init__.py b/src/sfg_walberla/boundaries/__init__.py
index 9a732bf3f6131d110ec614b49d0996bbcf5919f5..baeed122b07f2f2e9f07e36b6d63c5392688134f 100644
--- a/src/sfg_walberla/boundaries/__init__.py
+++ b/src/sfg_walberla/boundaries/__init__.py
@@ -1,3 +1,4 @@
from .hbb import SimpleHbbBoundary
+from .freeslip import FreeSlip
-__all__ = ["SimpleHbbBoundary"]
+__all__ = ["SimpleHbbBoundary", "FreeSlip"]
diff --git a/src/sfg_walberla/boundaries/boundary_utils.py b/src/sfg_walberla/boundaries/boundary_utils.py
new file mode 100644
index 0000000000000000000000000000000000000000..6618ddf1a50ceb2e5f00d90cb832130316c7287e
--- /dev/null
+++ b/src/sfg_walberla/boundaries/boundary_utils.py
@@ -0,0 +1,90 @@
+from functools import cache
+
+from pystencils import Field, FieldType, TypedSymbol, Assignment
+from pystencils.types import PsStructType, create_type
+
+from lbmpy import LBStencil
+from lbmpy.methods import AbstractLbMethod
+from lbmpy.boundaries.boundaryconditions import LbBoundary
+from lbmpy.advanced_streaming import Timestep
+
+from pystencilssfg import SfgComposer
+from pystencilssfg.composer.custom import CustomGenerator
+
+BoundaryIndexType = create_type("int32")
+
+HbbLinkType = PsStructType(
+ [
+ ("x", BoundaryIndexType),
+ ("y", BoundaryIndexType),
+ ("z", BoundaryIndexType),
+ ("dir", BoundaryIndexType),
+ ],
+ "walberla::sfg::HbbLink",
+)
+
+
+class WalberlaLbmBoundary:
+ idx_struct_type: PsStructType
+
+ @classmethod
+ @cache
+ def get_index_field(cls):
+ return Field(
+ "indexVector",
+ FieldType.INDEXED,
+ cls.idx_struct_type,
+ (0,),
+ (TypedSymbol("indexVectorLength", BoundaryIndexType), 1),
+ (1, 1),
+ )
+
+ def get_assignments(
+ self,
+ lb_method: AbstractLbMethod,
+ pdf_field: Field,
+ prev_timestep: Timestep = Timestep.BOTH,
+ streaming_pattern="pull",
+ ) -> list[Assignment]:
+ assert isinstance(self, LbBoundary)
+
+ index_field = self.get_index_field()
+
+ from lbmpy.advanced_streaming.indexing import BetweenTimestepsIndexing
+
+ indexing = BetweenTimestepsIndexing(
+ pdf_field,
+ lb_method.stencil,
+ prev_timestep,
+ streaming_pattern,
+ BoundaryIndexType,
+ BoundaryIndexType,
+ )
+
+ f_out, f_in = indexing.proxy_fields
+ dir_symbol = indexing.dir_symbol
+ inv_dir = indexing.inverse_dir_symbol
+
+ boundary_assignments = self(
+ f_out,
+ f_in,
+ dir_symbol,
+ inv_dir,
+ lb_method,
+ index_field,
+ None, # TODO: Fix force vector
+ )
+ boundary_assignments = indexing.substitute_proxies(boundary_assignments)
+
+ elements: list[Assignment] = []
+
+ index_arrs_node = indexing.create_code_node()
+ elements += index_arrs_node.get_array_declarations()
+
+ for node in self.get_additional_code_nodes(lb_method)[::-1]:
+ elements += node.get_array_declarations()
+
+ elements += [Assignment(dir_symbol, index_field[0]("dir"))]
+ elements += boundary_assignments.all_assignments
+
+ return elements
diff --git a/src/sfg_walberla/boundaries/freeslip.py b/src/sfg_walberla/boundaries/freeslip.py
new file mode 100644
index 0000000000000000000000000000000000000000..863d49eb13c812721d96b5541d14654f7247888e
--- /dev/null
+++ b/src/sfg_walberla/boundaries/freeslip.py
@@ -0,0 +1,129 @@
+from pystencils import Field, Assignment, CreateKernelConfig
+from pystencils.types import PsStructType
+
+from lbmpy.methods import AbstractLbMethod
+from lbmpy.boundaries.boundaryconditions import LbBoundary
+
+from pystencilssfg import SfgComposer
+from pystencilssfg.composer.custom import CustomGenerator
+
+from .hbb import BoundaryIndexType
+from .boundary_utils import WalberlaLbmBoundary
+from ..sweep import Sweep
+from ..api import SparseIndexList
+
+__all__ = ["FreeSlip"]
+
+
+class _IrregularSentinel:
+ def __repr__(self) -> str:
+ return "IRREGULAR"
+
+
+class FreeSlip(CustomGenerator):
+ """Free-Slip boundary condition"""
+
+ IRREGULAR = _IrregularSentinel()
+
+ def __init__(
+ self,
+ name: str,
+ lb_method: AbstractLbMethod,
+ pdf_field: Field,
+ wall_normal: tuple[int, int, int] | _IrregularSentinel,
+ ):
+ self._name = name
+ self._method = lb_method
+ self._pdf_field = pdf_field
+ self._wall_normal = wall_normal
+
+ def generate(self, sfg: SfgComposer) -> None:
+ if self._wall_normal == self.IRREGULAR:
+ self._generate_irregular(sfg)
+ else:
+ self._generate_regular(sfg)
+
+ def _generate_irregular(self, sfg: SfgComposer):
+ sfg.include("sfg/IrregularFreeSlip.hpp")
+
+ # Get waLBerla build config
+ bc_obj = WalberlaIrregularFreeSlip()
+
+ # Get assignments for bc
+ bc_asm = bc_obj.get_assignments(self._method, self._pdf_field)
+
+ # Build generator config
+ bc_cfg = CreateKernelConfig()
+ bc_cfg.index_dtype = BoundaryIndexType
+ index_field = bc_obj.get_index_field()
+ bc_cfg.index_field = index_field
+
+ # Prepare sweep
+ bc_sweep = Sweep(self._name, bc_asm, bc_cfg)
+
+ # Emit code
+ sfg.generate(bc_sweep)
+
+ # Build factory
+ factory_name = f"{self._name}Factory"
+ factory_crtp_base = f"walberla::sfg::IrregularFreeSlipFactory< {factory_name} >"
+ index_vector = SparseIndexList(
+ WalberlaIrregularFreeSlip.idx_struct_type, ref=True
+ ).var("indexVector")
+
+ sweep_type = bc_sweep.generated_class()
+ sweep_ctor_args = {
+ f"{self._pdf_field.name}Id": "this->pdfFieldID_",
+ f"{index_field.name}Id": index_vector.bufferId(),
+ }
+
+ stencil_name = self._method.stencil.name
+ sfg.include(f"stencil/{stencil_name}.h")
+
+ sfg.klass(factory_name, bases=[f"public {factory_crtp_base}"])(
+ sfg.public(
+ f"using Base = {factory_crtp_base};",
+ f"friend class {factory_crtp_base};",
+ f"using Stencil = walberla::stencil::{stencil_name};",
+ f"using Sweep = {sweep_type.get_dtype().c_string()};",
+ "using Base::Base;",
+ ),
+ sfg.private(
+ sfg.method(
+ "irregularFromIndexVector",
+ returns=sweep_type.get_dtype(),
+ inline=True,
+ )(sfg.expr("return {};", sweep_type.ctor(**sweep_ctor_args))),
+ ),
+ )
+
+ def _generate_regular(self, sfg: SfgComposer):
+ raise NotImplementedError("Regular-geometry free slip is not implemented yet")
+
+
+class WalberlaIrregularFreeSlip(LbBoundary, WalberlaLbmBoundary):
+ idx_struct_type = PsStructType(
+ (
+ ("x", BoundaryIndexType),
+ ("y", BoundaryIndexType),
+ ("z", BoundaryIndexType),
+ ("dir", BoundaryIndexType),
+ ("source_offset_x", BoundaryIndexType),
+ ("source_offset_y", BoundaryIndexType),
+ ("source_offset_z", BoundaryIndexType),
+ ("source_dir", BoundaryIndexType),
+ ),
+ "walberla::sfg::IrregularFreeSlipLinkInfo",
+ )
+
+ def __call__(
+ self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, force_vector
+ ):
+ source_cell = (
+ index_field("source_offset_x"),
+ index_field("source_offset_y"),
+ index_field("source_offset_z"),
+ )
+ source_dir = index_field("source_dir")
+
+ return Assignment(f_in(inv_dir[dir_symbol]), f_out[source_cell](source_dir))
diff --git a/src/sfg_walberla/boundaries/hbb.py b/src/sfg_walberla/boundaries/hbb.py
index fdc1683ccbfe3c1b010e5d47de0fd5ec16dfe8b4..d041f1b19a67df5da2107d4563c5b3b560b2d6f5 100644
--- a/src/sfg_walberla/boundaries/hbb.py
+++ b/src/sfg_walberla/boundaries/hbb.py
@@ -9,6 +9,7 @@ from pystencilssfg.composer.class_composer import SfgClassComposer
from pystencilssfg.composer.custom import CustomGenerator
from pystencilssfg.lang import AugExpr
+from lbmpy import LBStencil
from lbmpy.methods import AbstractLbMethod
from lbmpy.boundaries.boundaryconditions import LbBoundary
from lbmpy.advanced_streaming.indexing import Timestep
@@ -20,18 +21,7 @@ from ..sweep import (
combine_vectors,
)
from ..api import GhostLayerFieldPtr, BlockDataID, IBlockPtr, StructuredBlockForest, IndexListBufferPtr
-
-BoundaryIndexType = create_type("int32")
-
-HbbLinkType = PsStructType(
- [
- ("x", BoundaryIndexType),
- ("y", BoundaryIndexType),
- ("z", BoundaryIndexType),
- ("dir", BoundaryIndexType),
- ],
- "walberla::sfg::HbbLink",
-)
+from .boundary_utils import HbbLinkType, BoundaryIndexType
class HbbBoundaryProperties(SweepClassProperties):
diff --git a/src/sfg_walberla/reflection.py b/src/sfg_walberla/reflection.py
new file mode 100644
index 0000000000000000000000000000000000000000..13e422d94cf2329d339a23077a51a13553928652
--- /dev/null
+++ b/src/sfg_walberla/reflection.py
@@ -0,0 +1,29 @@
+from pystencilssfg.lang import CppClass, AugExpr
+from pystencilssfg.lang.types import cpptype
+from pystencilssfg.ir import SfgClass
+
+
+class GeneratedClassWrapperBase(CppClass):
+ _class: SfgClass
+ _namespace: str | None
+
+ def __init_subclass__(cls) -> None:
+ typename = (
+ f"{cls._namespace}::{cls._class.class_name}"
+ if cls._namespace is not None
+ else cls._class.class_name
+ )
+ cls.template = cpptype(typename)
+
+ def ctor(self, **kwargs) -> AugExpr:
+ for candidate_ctor in self._class.constructors():
+ ctor_argnames = [param.name for param in candidate_ctor.parameters]
+ if set(ctor_argnames) == set(kwargs.keys()):
+ break
+ else:
+ raise Exception(
+ f"No constructor of class {self._class.class_name} matches the argument names {kwargs.keys()}"
+ )
+
+ ctor_args = [kwargs[name] for name in ctor_argnames]
+ return self.ctor_bind(*ctor_args)
diff --git a/src/sfg_walberla/sweep.py b/src/sfg_walberla/sweep.py
index 5c16a448fef6d95daf36fa4d5097756365a348b0..a79d133f4381deea94d06e96ff57c33531e73e26 100644
--- a/src/sfg_walberla/sweep.py
+++ b/src/sfg_walberla/sweep.py
@@ -31,6 +31,8 @@ from pystencilssfg.lang import (
SrcVector,
strip_ptr_ref,
)
+from pystencilssfg.lang.types import CppTypeFactory, cpptype
+from .reflection import GeneratedClassWrapperBase
from .api import (
StructuredBlockForest,
GenericWalberlaField,
@@ -422,12 +424,18 @@ class Sweep(CustomGenerator):
config: CreateKernelConfig | None = None,
):
if config is not None:
- if config.ghost_layers is not None:
+ config = config.copy()
+
+ if config.get_option("ghost_layers") is not None:
raise ValueError(
"Specifying `ghost_layers` in your codegen config is invalid when generating a waLBerla sweep."
)
- elif config.iteration_slice is None:
- config = replace(config, ghost_layers=0)
+
+ if (
+ config.get_option("iteration_slice") is None
+ and config.get_option("index_field") is None
+ ):
+ config.ghost_layers = 0
else:
config = CreateKernelConfig(ghost_layers=0)
@@ -436,7 +444,7 @@ class Sweep(CustomGenerator):
if isinstance(assignments, AssignmentCollection):
self._assignments = assignments
else:
- self._assignments = AssignmentCollection(assignments)
+ self._assignments = AssignmentCollection(assignments) # type: ignore
self._gen_config = config
if self._gen_config.get_target() == Target.CUDA:
@@ -451,13 +459,21 @@ class Sweep(CustomGenerator):
# Map from shadow field to shadowed field
self._shadow_fields: dict[Field, Field] = dict()
+ # RESULTS - unset at this point
+ self._generated_class: type[GeneratedClassWrapperBase] | None = None
+
+ # CONFIGURATION
+
@property
def sparse(self) -> bool:
- return self._gen_config.index_field is not None
+ return self._gen_config.get_option("index_field") is not None
@sparse.setter
def sparse(self, sparse_iteration: bool):
if sparse_iteration:
+ if self._gen_config.get_option("index_field") is not None:
+ return
+
if self._gen_config.get_option("index_dtype") != create_type("int64"):
raise ValueError(
"Sparse sweeps require `int64_t` as their index data type. Check your code generator config."
@@ -471,8 +487,6 @@ class Sweep(CustomGenerator):
self._gen_config.index_field = None
self._gen_config.ghost_layers = 0
- # CONFIGURATION
-
def swap_fields(self, field: Field, shadow_field: Field):
if field in self._shadow_fields:
raise ValueError(f"Field swap for {field} was already registered.")
@@ -503,6 +517,8 @@ class Sweep(CustomGenerator):
knamespace = sfg.kernel_namespace(f"{self._name}_kernels")
assignments = BlockforestParameters.process(self._assignments)
+ # TODO: Get default config from waLBerla build system and override its entries
+ # from the user-provided config
khandle = knamespace.create(assignments, self._name, self._gen_config)
all_fields: dict[str, FieldInfo] = {
@@ -601,3 +617,24 @@ class Sweep(CustomGenerator):
sfg.public(*methods),
*shadows_cache.render(sfg),
)
+
+ gen_class = sfg.context.get_class(self._name)
+ assert gen_class is not None
+
+ class GenClassWrapper(GeneratedClassWrapperBase):
+ _class = gen_class
+ _namespace = sfg.context.fully_qualified_namespace
+
+ self._generated_class = GenClassWrapper
+
+ # CODE-GENERATION RESULTS
+ # These properties are only available after the sweep was generated
+
+ @property
+ def generated_class(self) -> type[GeneratedClassWrapperBase]:
+ if self._generated_class is None:
+ raise AttributeError(
+ "Generated class is unavailable - code generation was not run yet."
+ )
+
+ return self._generated_class
diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt
index 59a075a72448d40f6743ba8e7bfbcd15e8e71bc8..72acde9c5dc7ca88b51038b3798c79810a812883 100644
--- a/tests/CMakeLists.txt
+++ b/tests/CMakeLists.txt
@@ -1,6 +1,11 @@
cmake_minimum_required( VERSION 3.24 )
project( sfg-walberla-testsuite )
+set(WALBERLA_BUILD_TESTS OFF CACHE BOOL "")
+set(WALBERLA_BUILD_BENCHMARKS OFF CACHE BOOL "")
+set(WALBERLA_BUILD_SHOWCASES OFF CACHE BOOL "")
+set(WALBERLA_BUILD_TUTORIALS OFF CACHE BOOL "")
+
include(FetchContent)
FetchContent_Declare(
@@ -15,3 +20,7 @@ add_subdirectory(${CMAKE_SOURCE_DIR}/.. ${CMAKE_BINARY_DIR}/sfg-walberla)
# Test Directories
include(CTest)
+
+add_custom_target( SfgTests )
+
+add_subdirectory( FreeSlip )
diff --git a/tests/CMakePresets.json b/tests/CMakePresets.json
index 11b6c5797a4a848db2ef82c3e601e403848fb15a..48449c2d856200adfa38b7019aa2bcff44144615 100644
--- a/tests/CMakePresets.json
+++ b/tests/CMakePresets.json
@@ -11,6 +11,7 @@
"binaryDir": "${sourceDir}/build/${presetName}",
"generator": "Ninja",
"cacheVariables": {
+ "CMAKE_BUILD_TYPE": "Debug",
"WALBERLA_BUILD_TESTS": false
}
}
diff --git a/tests/FreeSlip/CMakeLists.txt b/tests/FreeSlip/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..e47cae2e1bef502e62a8b64dd2dea7f9fa396151
--- /dev/null
+++ b/tests/FreeSlip/CMakeLists.txt
@@ -0,0 +1,6 @@
+add_executable( TestIrregularFreeSlip TestIrregularFreeSlip.cpp )
+walberla_generate_sources( TestIrregularFreeSlip SCRIPTS IrregularFreeSlip.py )
+target_link_libraries( TestIrregularFreeSlip core blockforest field geometry sfg_walberla )
+add_test( NAME TestIrregularFreeSlip COMMAND TestIrregularFreeSlip )
+
+add_dependencies( SfgTests TestIrregularFreeSlip )
diff --git a/tests/FreeSlip/IrregularFreeSlip.py b/tests/FreeSlip/IrregularFreeSlip.py
new file mode 100644
index 0000000000000000000000000000000000000000..3fa69cf2042e054618db15ec1779945b9198214a
--- /dev/null
+++ b/tests/FreeSlip/IrregularFreeSlip.py
@@ -0,0 +1,54 @@
+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.macroscopic_value_kernels import macroscopic_values_setter
+
+
+from sfg_walberla import Sweep
+from sfg_walberla.boundaries import FreeSlip
+
+with SourceFileGenerator() as sfg:
+ sfg.namespace("TestIrregularFreeSlip::gen")
+
+ stencil = LBStencil(Stencil.D3Q19)
+ d, q = stencil.D, stencil.Q
+ f: Field
+ f_tmp: Field
+ f, f_tmp, rho, u = fields(
+ f"f({q}), f_tmp({q}), rho, u({d}): double[{d}D]", layout="fzyx"
+ ) # type: ignore
+
+ stencil_name = stencil.name
+ sfg.include(f"stencil/{stencil_name}.h")
+ sfg.code(f"using LbStencil = walberla::stencil::{stencil_name};")
+
+ lbm_config = LBMConfig(
+ stencil=stencil,
+ output={"density": rho, "velocity": u},
+ )
+ lbm_opt = LBMOptimisation(
+ symbolic_field=f,
+ symbolic_temporary_field=f_tmp,
+ )
+
+ lb_update = create_lb_update_rule(lbm_config=lbm_config, lbm_optimisation=lbm_opt)
+ assert lb_update is not None
+
+ lb_update_sweep = Sweep("LbStreamCollide", lb_update)
+ lb_update_sweep.swap_fields(f, f_tmp)
+ sfg.generate(lb_update_sweep)
+
+ lb_init = macroscopic_values_setter(
+ lb_update.method,
+ density=sp.Symbol("density"),
+ velocity=u.center_vector,
+ pdfs=f,
+ )
+ lb_init_sweep = Sweep("LbInit", lb_init)
+ sfg.generate(lb_init_sweep)
+
+ irreg_freeslip = FreeSlip("FreeSlipIrregular", lb_update.method, f, wall_normal=FreeSlip.IRREGULAR)
+ sfg.generate(irreg_freeslip)
+
diff --git a/tests/FreeSlip/TestIrregularFreeSlip.cpp b/tests/FreeSlip/TestIrregularFreeSlip.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..c08df34308cc4e8c64d02717cfe1c091c0ecb5a4
--- /dev/null
+++ b/tests/FreeSlip/TestIrregularFreeSlip.cpp
@@ -0,0 +1,147 @@
+
+#include "core/all.h"
+#include "blockforest/all.h"
+
+#include "field/all.h"
+#include "field/communication/PackInfo.h"
+#include "field/communication/StencilRestrictedPackInfo.h"
+
+#include "geometry/all.h"
+
+#include "vtk/all.h"
+
+#include <array>
+
+#include "gen/IrregularFreeSlip.hpp"
+
+namespace TestIrregularFreeSlip
+{
+ using namespace walberla;
+
+ using PdfField_T = field::GhostLayerField<real_t, gen::LbStencil::Q>;
+ using ScalarField_T = field::GhostLayerField<real_t, 1>;
+ using VectorField_T = field::GhostLayerField<real_t, gen::LbStencil::D>;
+ using FlagField_T = FlagField<uint8_t>;
+ using CommScheme = blockforest::communication::UniformBufferedScheme<gen::LbStencil>;
+ using PdfsPackInfo = field::communication::StencilRestrictedPackInfo<PdfField_T, gen::LbStencil>;
+
+ /**
+ * Pipe flow with circular cross-section and free-slip walls.
+ * The pipe flow is initialized with a constant velocity in x-direction.
+ * As free-slip walls do not impose any friction, velocity should remain constant in time.
+ */
+ void freeSlipPipe(Environment &env)
+ {
+ std::array<uint_t, 3> blocks{1, 1, 1};
+ std::array<uint_t, 3> cpb{4, 32, 32};
+ std::array<bool, 3> periodic{true, false, false};
+
+ auto sbfs = blockforest::createUniformBlockGrid(
+ blocks[0], blocks[1], blocks[2],
+ cpb[0], cpb[1], cpb[2],
+ 1.0,
+ true,
+ periodic[0], periodic[1], periodic[2]
+ );
+
+ BlockDataID pdfsId = field::addToStorage<PdfField_T>(sbfs, "f", real_c(0.0));
+ BlockDataID rhoId = field::addToStorage<ScalarField_T>(sbfs, "rho", real_c(0.0));
+ BlockDataID uId = field::addToStorage<VectorField_T>(sbfs, "u", real_c(0.0));
+
+ const BlockDataID flagFieldId = field::addFlagFieldToStorage<FlagField_T>(sbfs, "flagField");
+ const FlagUID fluidFlagUid{"Fluid"};
+ const FlagUID freeSlipFlagUID{"FreeSlip"};
+
+ const CellInterval allCells{{0, 0, 0}, {sbfs->getNumberOfXCellsPerBlock() - 1, sbfs->getNumberOfYCellsPerBlock() - 1, sbfs->getNumberOfZCellsPerBlock() - 1}};
+ const CellInterval allCellsWithGl{{-1, 0, 0}, {sbfs->getNumberOfXCellsPerBlock(), sbfs->getNumberOfYCellsPerBlock() - 1, sbfs->getNumberOfZCellsPerBlock() - 1}};
+
+ const real_t pipeRadius{14.0};
+ const Vector3<real_t> pipeAnchor{0.0, 16.0, 16.0};
+ const real_t maxVelocity{0.02};
+
+ for (auto &block : *sbfs)
+ {
+ FlagField_T &flagField = *block.getData<FlagField_T>(flagFieldId);
+ const uint8_t freeSlipFlag{flagField.getOrRegisterFlag(freeSlipFlagUID)};
+
+ VectorField_T &velField = *block.getData<VectorField_T>(uId);
+
+ for (Cell c : allCellsWithGl)
+ {
+ Cell globalCell{c};
+ sbfs->transformBlockLocalToGlobalCell(globalCell, block);
+ Vector3<real_t> cellCenter{sbfs->getCellCenter(globalCell)};
+ cellCenter[0] = 0.0;
+
+ Vector3<real_t> initVelocity;
+ real_t radialDistance = (cellCenter - pipeAnchor).length();
+ if (radialDistance > pipeRadius)
+ {
+ flagField.addFlag(c, freeSlipFlag);
+ initVelocity = Vector3<real_t>{NAN};
+ }
+ else
+ {
+ initVelocity = Vector3<real_t>{maxVelocity, 0., 0.};
+ }
+
+ velField.get(c, 0) = initVelocity[0];
+ velField.get(c, 1) = initVelocity[1];
+ velField.get(c, 2) = initVelocity[2];
+ }
+ }
+
+ geometry::setNonBoundaryCellsToDomain<FlagField_T>(*sbfs, flagFieldId, fluidFlagUid);
+
+ auto flagsOutput = field::createVTKOutput<FlagField_T>(flagFieldId, *sbfs, "flags");
+ flagsOutput();
+
+ gen::LbInit initialize { pdfsId, uId, 1.0 };
+
+ for(auto& b : *sbfs){
+ initialize(&b);
+ }
+
+ gen::FreeSlipIrregular freeSlip{
+ gen::FreeSlipIrregularFactory(sbfs, pdfsId).fromFlagField<FlagField_T>(flagFieldId, freeSlipFlagUID, fluidFlagUid)
+ };
+
+ gen::LbStreamCollide streamCollide { pdfsId, rhoId, uId, 1.0 };
+
+ CommScheme comm{sbfs};
+ auto pdfsPackInfo = std::make_shared<PdfsPackInfo>(pdfsId);
+ comm.addPackInfo(pdfsPackInfo);
+
+ auto velOutput = field::createVTKOutput<VectorField_T>(uId, *sbfs, "vel");
+
+ for(uint_t i = 0; i < 10; ++i){
+ comm();
+ for (auto &block : *sbfs)
+ {
+ velOutput();
+ freeSlip(&block);
+ streamCollide(&block);
+
+ VectorField_T & velField = *block.getData< VectorField_T >(uId);
+ FlagField_T & flagField = *block.getData< FlagField_T >(flagFieldId);
+ uint8_t fluidFlag = flagField.getFlag(fluidFlagUid);
+
+ for(Cell c : allCells){
+ if( flagField.isFlagSet(c, fluidFlag) ){
+ WALBERLA_CHECK_FLOAT_EQUAL( velField.get(c, 0), maxVelocity );
+ WALBERLA_CHECK_FLOAT_EQUAL( velField.get(c, 1), 0. );
+ WALBERLA_CHECK_FLOAT_EQUAL( velField.get(c, 2), 0. );
+ }
+ }
+ }
+ }
+
+ velOutput();
+ }
+}
+
+int main(int argc, char **argv)
+{
+ walberla::Environment env{argc, argv};
+ TestIrregularFreeSlip::freeSlipPipe(env);
+}
diff --git a/user_manual/conf.py b/user_manual/conf.py
index 2bc06da4f2d273c00770914d3e7fd0e7555755e0..6dfd366ee11f1cc57897e825c1ac61bcd85f3bd7 100644
--- a/user_manual/conf.py
+++ b/user_manual/conf.py
@@ -21,7 +21,7 @@ extensions = [
]
# templates_path = ['_templates']
-exclude_patterns = ["build"]
+exclude_patterns = ["*/build"]
myst_enable_extensions = [
"colon_fence",
diff --git a/user_manual/examples/CMakeLists.txt b/user_manual/examples/CMakeLists.txt
index f498cc3603a44d597f611da0530ef297a238f97a..499f1f083bf4419e3f4c8d72c7ccb55123500ad4 100644
--- a/user_manual/examples/CMakeLists.txt
+++ b/user_manual/examples/CMakeLists.txt
@@ -16,6 +16,7 @@ add_subdirectory(${CMAKE_SOURCE_DIR}/../.. ${CMAKE_BINARY_DIR}/sfg-walberla)
add_subdirectory( GeneratorScriptBasics )
add_subdirectory( ForceDrivenChannel )
add_subdirectory( SparseSpiral )
+add_subdirectory( FreeSlipBubble )
add_custom_target( Examples )
add_dependencies(
@@ -23,4 +24,5 @@ add_dependencies(
Ex_GeneratorScriptBasics
Ex_ForceDrivenChannel
SparseSpiral
+ FreeSlipBubble
)
diff --git a/user_manual/examples/FreeSlipBubble/CMakeLists.txt b/user_manual/examples/FreeSlipBubble/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..cd3e7af346973a31ce122c55531821489e74ce14
--- /dev/null
+++ b/user_manual/examples/FreeSlipBubble/CMakeLists.txt
@@ -0,0 +1,6 @@
+walberla_link_files_to_builddir( *.prm )
+
+add_executable( FreeSlipBubble FreeSlipBubble.cpp )
+
+walberla_generate_sources( FreeSlipBubble SCRIPTS LbmAlgorithms.py )
+target_link_libraries( FreeSlipBubble core blockforest field geometry stencil sfg_walberla )
diff --git a/user_manual/examples/FreeSlipBubble/FreeSlipBubble.cpp b/user_manual/examples/FreeSlipBubble/FreeSlipBubble.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..be5f93da51d38ebd445127201ba36a31d74f1900
--- /dev/null
+++ b/user_manual/examples/FreeSlipBubble/FreeSlipBubble.cpp
@@ -0,0 +1,168 @@
+#include "blockforest/all.h"
+#include "blockforest/communication/UniformBufferedScheme.h"
+
+#include "core/all.h"
+
+#include "domain_decomposition/SharedSweep.h"
+
+#include "field/all.h"
+#include "field/communication/PackInfo.h"
+#include "field/communication/StencilRestrictedPackInfo.h"
+
+#include "geometry/InitBoundaryHandling.h"
+
+#include "stencil/all.h"
+
+#include "timeloop/all.h"
+
+#include "vtk/all.h"
+
+#include "gen/LbmAlgorithms.hpp"
+
+namespace FreeSlipBubble
+{
+ using namespace walberla;
+ using namespace blockforest;
+
+ using std::make_unique;
+ using std::shared_ptr;
+
+ using ScalarField_T = field::GhostLayerField<real_t, 1>;
+ using VectorField_T = field::GhostLayerField<real_t, 3>;
+
+ using LbStencil = stencil::D3Q19;
+ using PdfField_T = field::GhostLayerField<real_t, LbStencil::Q>;
+
+ using CommScheme = blockforest::communication::UniformBufferedScheme<LbStencil>;
+ using PdfsPackInfo = field::communication::StencilRestrictedPackInfo<PdfField_T, LbStencil>;
+
+ using FlagField_T = FlagField<uint8_t>;
+
+ void run(const shared_ptr<Config> &config)
+ {
+ auto blocks = createUniformBlockGridFromConfig(config);
+
+ Config::BlockHandle simParams = config->getBlock("Parameters");
+ const real_t omega{simParams.getParameter<real_t>("omega")};
+ const real_t density{simParams.getParameter<real_t>("density", 1.0)};
+ const Vector3<real_t> velocity = simParams.getParameter<Vector3<real_t>>("velocity", Vector3<real_t>{0.05, 0.0, 0.0});
+
+ BlockDataID pdfsId = field::addToStorage<PdfField_T>(blocks, "pdfs", real_c(0.0), field::fzyx, 1);
+ BlockDataID rhoId = field::addToStorage<ScalarField_T>(blocks, "rho", real_c(1.0), field::fzyx, 0);
+ BlockDataID uId = field::addToStorage<VectorField_T>(blocks, "u", real_c(0.0), field::fzyx, 0);
+
+ gen::LbInit lbInit{pdfsId, density, velocity};
+
+ for (auto &block : *blocks)
+ {
+ lbInit(&block);
+ }
+
+ // Stream-Collide
+
+ auto streamCollide = make_shared<gen::LbStreamCollide>(pdfsId, rhoId, uId, omega);
+
+ // Communication Setup
+
+ CommScheme comm{blocks};
+ auto pdfsPackInfo = std::make_shared<PdfsPackInfo>(pdfsId);
+ comm.addPackInfo(pdfsPackInfo);
+
+ // Boundary Conditions
+
+ const BlockDataID flagFieldId = field::addFlagFieldToStorage<FlagField_T>(blocks, "flagField");
+ const FlagUID fluidFlagUid{"Fluid"};
+ const FlagUID boundaryFlagUID{"FreeSlip"};
+
+ const CellInterval allCells{{0, 0, 0}, {blocks->getNumberOfXCellsPerBlock() - 1, blocks->getNumberOfYCellsPerBlock() - 1, blocks->getNumberOfZCellsPerBlock() - 1}};
+ const Vector3< real_t > sphereCenter { 16., 16., 16. };
+ const real_t sphereRadius { 8. };
+
+ for (auto &block : *blocks)
+ {
+ FlagField_T &flagField = *block.getData<FlagField_T>(flagFieldId);
+ const uint8_t freeSlipFlag{flagField.getOrRegisterFlag(boundaryFlagUID)};
+
+ PdfField_T & pdfField = *block.getData< PdfField_T >(pdfsId);
+
+ for (Cell c: allCells){
+ Cell globalCell { c };
+ blocks->transformBlockLocalToGlobalCell( globalCell, block );
+ Vector3< real_t > cellCenter { blocks->getCellCenter(globalCell) };
+
+ if((cellCenter - sphereCenter).sqrLength() < sphereRadius * sphereRadius){
+ flagField.addFlag(c, freeSlipFlag);
+ for( uint_t q = 0; q < 19; ++q ){
+ pdfField.get(c, q) = NAN;
+ }
+ }
+ }
+ }
+
+ geometry::setNonBoundaryCellsToDomain<FlagField_T>(*blocks, flagFieldId, fluidFlagUid);
+
+ auto flagFieldOutput = field::createVTKOutput<FlagField_T>(flagFieldId, *blocks, "flagField", 1, 0);
+ flagFieldOutput();
+
+ // begin freeSlipFactory
+ gen::FreeSlipFactory freeSlipFactory{ blocks, pdfsId };
+ // end freeSlipFactory
+ // begin freeSlipFromFlagField
+ gen::FreeSlip freeSlip {
+ freeSlipFactory.fromFlagField<FlagField_T>(
+ flagFieldId,
+ boundaryFlagUID,
+ fluidFlagUid
+ )
+ };
+ // end freeSlipFromFlagField
+
+ // Timeloop
+ const uint_t numTimesteps{simParams.getParameter<uint_t>("timesteps")};
+ SweepTimeloop loop{blocks->getBlockStorage(), numTimesteps};
+
+ loop.add() << BeforeFunction(comm) << Sweep(freeSlip);
+ loop.add() << Sweep(makeSharedSweep(streamCollide));
+
+ RemainingTimeLogger logger{numTimesteps};
+ loop.addFuncAfterTimeStep(logger);
+
+ // VTK Output
+
+ Config::BlockHandle outputParams = config->getBlock("Output");
+
+ const uint_t vtkWriteFrequency = outputParams.getParameter<uint_t>("vtkWriteFrequency", 0);
+ if (vtkWriteFrequency > 0)
+ {
+ auto vtkOutput = vtk::createVTKOutput_BlockData(*blocks, "vtk", vtkWriteFrequency, 0, false, "vtk_out",
+ "simulation_step", false, true, true, false, 0);
+
+ auto densityWriter = make_shared<field::VTKWriter<ScalarField_T, float32>>(rhoId, "density");
+ vtkOutput->addCellDataWriter(densityWriter);
+
+ auto velWriter = make_shared<field::VTKWriter<VectorField_T, float32>>(uId, "velocity");
+ vtkOutput->addCellDataWriter(velWriter);
+
+ const cell_idx_t xCells = cell_idx_c(blocks->getNumberOfXCells());
+ const cell_idx_t yCells = cell_idx_c(blocks->getNumberOfYCells());
+ const cell_idx_t zCells = cell_idx_c(blocks->getNumberOfZCells());
+
+ loop.addFuncAfterTimeStep(vtk::writeFiles(vtkOutput), "VTK Output");
+ }
+
+ // Run the Simulation
+
+ WALBERLA_LOG_INFO_ON_ROOT("Commencing simulation with " << numTimesteps << " timesteps")
+
+ loop.run();
+ }
+}
+
+int main(int argc, char **argv)
+{
+ walberla::Environment env{argc, argv};
+
+ FreeSlipBubble::run(env.config());
+
+ return EXIT_SUCCESS;
+}
diff --git a/user_manual/examples/FreeSlipBubble/LbmAlgorithms.py b/user_manual/examples/FreeSlipBubble/LbmAlgorithms.py
new file mode 100644
index 0000000000000000000000000000000000000000..a2975f2617a73d9b129dc7bd55b1a308e9a37b2a
--- /dev/null
+++ b/user_manual/examples/FreeSlipBubble/LbmAlgorithms.py
@@ -0,0 +1,59 @@
+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.macroscopic_value_kernels import macroscopic_values_setter
+
+
+from sfg_walberla import Sweep
+
+with SourceFileGenerator() as sfg:
+ sfg.namespace("FreeSlipBubble::gen")
+
+ stencil = LBStencil(Stencil.D3Q19)
+ d, q = stencil.D, stencil.Q
+ pdfs: Field
+ pdfs_tmp: Field
+ pdfs, pdfs_tmp, rho, u = fields(
+ f"f({q}), f_tmp({q}), rho, u({d}): double[{d}D]", layout="fzyx"
+ ) # type: ignore
+
+ lbm_config = LBMConfig(
+ stencil=Stencil.D3Q19,
+ output={"density": rho, "velocity": u},
+ )
+ lbm_opt = LBMOptimisation(
+ symbolic_field=pdfs,
+ symbolic_temporary_field=pdfs_tmp,
+ )
+
+ lb_update = create_lb_update_rule(lbm_config=lbm_config, lbm_optimisation=lbm_opt)
+ lb_method = lb_update.method
+ assert lb_update is not None
+
+ lb_update_sweep = Sweep("LbStreamCollide", lb_update)
+ lb_update_sweep.swap_fields(pdfs, pdfs_tmp)
+ sfg.generate(lb_update_sweep)
+
+ lb_init = macroscopic_values_setter(
+ lb_update.method,
+ density=sp.Symbol("density"),
+ velocity=sp.symbols(f"velocity_:{d}"),
+ pdfs=pdfs,
+ )
+ lb_init_sweep = Sweep("LbInit", lb_init)
+ sfg.generate(lb_init_sweep)
+
+ # begin irregular-freeslip
+ from sfg_walberla.boundaries import FreeSlip
+
+ freeslip = FreeSlip(
+ "FreeSlip",
+ lb_method,
+ pdfs,
+ wall_normal=FreeSlip.IRREGULAR,
+ )
+
+ sfg.generate(freeslip)
+ # end irregular-freeslip
diff --git a/user_manual/examples/FreeSlipBubble/Periodic.prm b/user_manual/examples/FreeSlipBubble/Periodic.prm
new file mode 100644
index 0000000000000000000000000000000000000000..2fd9fbd63f7836cce2ffb43723485ca4e01fdbc9
--- /dev/null
+++ b/user_manual/examples/FreeSlipBubble/Periodic.prm
@@ -0,0 +1,17 @@
+DomainSetup
+{
+ blocks < 1, 1, 1 >;
+ cellsPerBlock < 32, 32, 32 >;
+ periodic < 1, 1, 1 >;
+}
+
+Parameters
+{
+ omega 1.0;
+ timesteps 2000;
+}
+
+Output
+{
+ vtkWriteFrequency 100;
+}
diff --git a/user_manual/guides/BoundaryConditions.md b/user_manual/guides/BoundaryConditions.md
new file mode 100644
index 0000000000000000000000000000000000000000..7480c64f5648fd34511320785d1e450dea36bdb2
--- /dev/null
+++ b/user_manual/guides/BoundaryConditions.md
@@ -0,0 +1,56 @@
+(guide-boundary-conditions)=
+# Boundary Conditions
+
+## Sparse Boundary Sweeps for Irregular Geometries
+
+To cover the boundary conditions at surfaces with irregular geometries
+(i.e. surfaces not aligned with the cell grid),
+we need to use a sparse implementation of the boundary sweep.
+
+### Generating Sparse Link-Wise LBM Boundary Sweeps
+
+For link-wise lattice Boltzmann boundary conditions, sparsity is enabled by
+setting the `wall_normal` parameter to `IRREGULAR`:
+
+```{literalinclude} ../examples/FreeSlipBubble/LbmAlgorithms.py
+:language: Python
+:start-after: "# begin irregular-freeslip"
+:end-before: "# end irregular-freeslip"
+:linenos:
+:dedent: 4
+:emphasize-lines: 7
+```
+
+The above code will generate two classes:
+ - The `FreeSlip` sweep class, which is a [sparse sweep](#guide-sparse-sweeps) that applies
+ the free-slip boundary condition at each lattice link crossing the boundary surface;
+ - and the `FreeSlipFactory` class; which produces instances of the `FreeSlip` sweep and populates its
+ index vector from user-provided geometry information.
+
+### Populating the Index Vector
+
+Sparse boundary sweeps are created, and their index vectors populated, by their associated factory.
+Create the factory object by passing it any required parameters:
+
+```{literalinclude} ../examples/FreeSlipBubble/FreeSlipBubble.cpp
+:language: C++
+:start-after: "// begin freeSlipFactory"
+:end-before: "// end freeSlipFactory"
+:dedent: 8
+```
+
+#### From a Flag Field
+
+One possibility to populate the sparse sweep is to extract boundary information stored in a *flag field*.
+This requires that each *solid* cell considered by the boundary condition is marked with a *boundary flag*,
+and each valid *fluid* cell is marked with the *fluid* flag.
+Assuming that this information is stored in a flag field of type `FlagField_T` identified by `flagFieldId`,
+and `boundaryFlagUID` and `fluidFlagUID` respectively identify the boundary and fluid flags,
+we can create a sparse boundary sweep using the `fromFlagField` method of the factory:
+
+```{literalinclude} ../examples/FreeSlipBubble/FreeSlipBubble.cpp
+:language: C++
+:start-after: "// begin freeSlipFromFlagField"
+:end-before: "// end freeSlipFromFlagField"
+:dedent: 8
+```
diff --git a/user_manual/guides/SparseSweeps.md b/user_manual/guides/SparseSweeps.md
index 159313f0c8542839d117a651bacd9d8165128bdf..822d92fbc3894a07e09def2e2797b95d0f349e46 100644
--- a/user_manual/guides/SparseSweeps.md
+++ b/user_manual/guides/SparseSweeps.md
@@ -1,3 +1,4 @@
+(guide-sparse-sweeps)=
# Sparse Sweeps
:::{admonition} Example Code
diff --git a/user_manual/index.md b/user_manual/index.md
index a48e3fe493164d216d069af31192de6f64a37b51..29564ac46045ef4171c75bd23f1ee63d1e19985c 100644
--- a/user_manual/index.md
+++ b/user_manual/index.md
@@ -24,6 +24,7 @@ examples/GeneratorScriptBasics/GeneratorScriptBasics
:maxdepth: 1
guides/SparseSweeps
+guides/BoundaryConditions
:::
:::{toctree}