Bugfixes in lbmpy boundary handling

boundaries/boundaryhandling.py

@@ -3,14 +3,14 @@ from lbmpy.boundaries.handlinginterface import GenericBoundaryHandling, FlagFiel
 from lbmpy.boundaries.periodicityhandling import PeriodicityHandling
 
 
-class BoundaryHandling(GenericBoundaryHandling, PeriodicityHandling):
+class BoundaryHandling(PeriodicityHandling, GenericBoundaryHandling):  # important: first periodicity, then boundary
 
     def __init__(self, pdfField, domainShape, lbMethod, ghostLayers=1, target='cpu', openMP=True, flagDtype=np.uint32):
         shapeWithGl = [a + 2 * ghostLayers for a in domainShape]
         flagFieldInterface = NumpyFlagFieldInterface(shapeWithGl, flagDtype)
 
         GenericBoundaryHandling.__init__(self, flagFieldInterface, pdfField, lbMethod, ghostLayers, target, openMP)
-        PeriodicityHandling.__init__(self, domainShape)
+        PeriodicityHandling.__init__(self, list(domainShape) + [len(lbMethod.stencil)])
 
     def setBoundary(self, boundaryObject, indexExpr=None, maskCallback=None, sliceNormalizationGhostLayers=1):
         mask = None
@@ -65,3 +65,4 @@ class NumpyFlagFieldInterface(FlagFieldInterface):
         self._boundaryObjectToName = {}
         self._boundaryObjectToFlag = {'fluid': np.dtype(self._array.dtype)(2**0)}
         self._nextFreeExponent = 1
+

boundaries/periodicityhandling.py

@@ -6,7 +6,7 @@ class PeriodicityHandling(object):
         self._spatialShape = fieldShape[:-1]
         self._indexShape = fieldShape[-1]
         self._periodicity = list(periodicity)
-        self._dirty = False
+        self._periodicityDirty = False
         self._periodicityKernels = []
 
     @property
@@ -20,13 +20,16 @@ class PeriodicityHandling(object):
             assert isinstance(d, bool)
 
         self._periodicity = [x, y, z]
-        self._dirty = True
+        self._periodicityDirty = True
 
     def __call__(self, **kwargs):
         for k in self._periodicityKernels:
             k(**kwargs)
 
     def prepare(self):
+        if not self._periodicityDirty:
+            return
+
         self._periodicityKernels = []
         dim = len(self.flagField.shape)
         if dim == 2:
@@ -60,5 +63,5 @@ class PeriodicityHandling(object):
             else:
                 assert False
 
-        self._dirty = False
+        self._periodicityDirty = False
 

scenarios.py

@@ -141,12 +141,12 @@ def createForceDrivenChannel(force=1e-6, domainSize=None, dim=2, radius=None, le
             boundaryHandling.setBoundary(wallBoundary, sliceFromDirection(direction, dim))
     elif dim == 3:
         if roundChannel:
-            wallIdx = boundaryHandling.addBoundary(wallBoundary)
-            ff = boundaryHandling.flagField
-            yMid = ff.shape[1] // 2
-            zMid = ff.shape[2] // 2
-            y, z = np.meshgrid(range(ff.shape[1]), range(ff.shape[2]))
-            ff[(y - yMid) ** 2 + (z - zMid) ** 2 > radius ** 2] = wallIdx
+            def circleMaskCb(x, y, z):
+                yMid = np.max(y) // 2
+                zMid = np.max(z) // 2
+                return (y - yMid) ** 2 + (z - zMid) ** 2 > radius ** 2
+
+            boundaryHandling.setBoundary(wallBoundary, maskCallback=circleMaskCb)
         else:
             for direction in ('N', 'S', 'T', 'B'):
                 boundaryHandling.setBoundary(wallBoundary, sliceFromDirection(direction, dim))