Worked on waLBerla kernel generation

creationfunctions.py

@@ -53,6 +53,8 @@ General:
 - ``output={}``: a dictionary mapping macroscopic quantites e.g. the strings 'density' and 'velocity' to pystencils
                 fields. In each timestep the corresponding quantities are written to the given fields.
 - ``velocityInput``: symbolic field where the velocities are read from (for advection diffusion LBM)
+- ``kernelType``: supported values: 'streamPullCollide' (default), 'collideOnly' 
+
 
 Entropic methods:
 
@@ -153,7 +155,8 @@ from lbmpy.methods.relaxationrates import relaxationRateFromMagicNumber
 from lbmpy.stencils import getStencil
 import lbmpy.forcemodels as forceModels
 from lbmpy.simplificationfactory import createSimplificationStrategy
-from lbmpy.updatekernels import createStreamPullKernel, createPdfArray
+from lbmpy.updatekernels import createStreamPullCollideKernel, createPdfArray, createCollideOnlyKernel, \
+    createStreamPullOnlyKernel
 from pystencils.equationcollection.equationcollection import EquationCollection
 
 
@@ -166,7 +169,7 @@ def updateWithDefaultParameters(params, optParams, failOnUnknownParameter=True):
         'equilibriumAccuracyOrder': 2,
         'c_s_sq': sp.Rational(1, 3),
 
-        'forceModel': 'none',  # can be 'simple', 'luo' or 'guo'
+        'forceModel': 'none',
         'force': (0, 0, 0),
         'useContinuousMaxwellianEquilibrium': True,
         'cumulant': False,
@@ -178,6 +181,11 @@ def updateWithDefaultParameters(params, optParams, failOnUnknownParameter=True):
 
         'output': {},
         'velocityInput': None,
+
+        'kernelType': 'streamPullCollide',
+
+        'fieldName': 'src',
+        'secondFieldName': 'dst',
     }
 
     defaultOptimizationDescription = {
@@ -340,23 +348,31 @@ def createLatticeBoltzmannUpdateRule(lbMethod=None, optimizationParams={}, **kwa
         else:
             collisionRule = addEntropyCondition(collisionRule, omegaOutputField=params['omegaOutputField'])
 
+    kernelCreateArgs = {}
+
     if 'fieldSize' in optParams and optParams['fieldSize']:
         npField = createPdfArray(optParams['fieldSize'], len(stencil), layout=optParams['fieldLayout'])
-        updateRule = createStreamPullKernel(collisionRule, numpyField=npField,
-                                            builtinPeriodicity=optParams['builtinPeriodicity'])
+        kernelCreateArgs['numpyField'] = npField
     else:
         if 'pdfArr' in optParams:
-            updateRule = createStreamPullKernel(collisionRule, numpyField=optParams['pdfArr'],
-                                                builtinPeriodicity=optParams['builtinPeriodicity'])
+            kernelCreateArgs['numpyField'] = optParams['pdfArr']
         else:
             layoutName = optParams['fieldLayout']
             if layoutName == 'fzyx' or 'zyxf':
                 dim = len(stencil[0])
                 layoutName = tuple(reversed(range(dim)))
-            updateRule = createStreamPullKernel(collisionRule, genericLayout=layoutName,
-                                                builtinPeriodicity=optParams['builtinPeriodicity'])
-
-    return updateRule
+            kernelCreateArgs['genericLayout'] = layoutName
+
+    if params['kernelType'] == 'streamPullCollide':
+        kernelCreateArgs['srcFieldName'] = params['fieldName']
+        kernelCreateArgs['dstFieldName'] = params['secondFieldName']
+        kernelCreateArgs['builtinPeriodicity'] = optParams['builtinPeriodicity']
+        return createStreamPullCollideKernel(collisionRule, **kernelCreateArgs)
+    elif params['kernelType'] == 'collideOnly':
+        kernelCreateArgs['fieldName'] = params['fieldName']
+        return createCollideOnlyKernel(collisionRule, **kernelCreateArgs)
+    else:
+        raise ValueError("Invalid value of parameter 'kernelType'", params['kernelType'])
 
 
 def createLatticeBoltzmannMethod(**params):

macroscopic_value_kernels.py

@@ -184,9 +184,9 @@ def compileAdvancedVelocitySetter(method, velocityArray, velocityRelaxationRate=
     :param optimizationParams: dictionary with optimization hints
     :return: stream-collide update function
     """
-    from lbmpy.updatekernels import createStreamPullKernel
+    from lbmpy.updatekernels import createStreamPullCollideKernel
     from lbmpy.creationfunctions import createLatticeBoltzmannAst, createLatticeBoltzmannFunction
     newCollisionRule = createAdvancedVelocitySetterCollisionRule(method, velocityArray, velocityRelaxationRate)
-    updateRule = createStreamPullKernel(newCollisionRule, pdfArr, genericLayout=fieldLayout)
+    updateRule = createStreamPullCollideKernel(newCollisionRule, pdfArr, genericLayout=fieldLayout)
     ast = createLatticeBoltzmannAst(updateRule, optimizationParams)
     return createLatticeBoltzmannFunction(ast, optimizationParams)

plot2d.py

@@ -44,6 +44,15 @@ def scalarField(field, **kwargs):
     return res
 
 
+def multipleScalarFields(field, **kwargs):
+    subPlots = field.shape[-1]
+    for i in range(subPlots):
+        subplot(1, subPlots, i + 1)
+        title(str(i))
+        scalarField(field[..., i])
+        colorbar()
+
+
 def plotBoundaryHandling(boundaryHandling, boundaryNameToColor=None):
     """
     Shows boundary cells

updatekernels.py

 import numpy as np
 import sympy as sp
 from pystencils import Field
+from pystencils.equationcollection.equationcollection import EquationCollection
 from pystencils.field import createNumpyArrayWithLayout, layoutStringToTuple
 from pystencils.sympyextensions import fastSubs
-from lbmpy.fieldaccess import StreamPullTwoFieldsAccessor, Pseudo2DTwoFieldsAccessor, PeriodicTwoFieldsAccessor
+from lbmpy.fieldaccess import StreamPullTwoFieldsAccessor, Pseudo2DTwoFieldsAccessor, PeriodicTwoFieldsAccessor, \
+    CollideOnlyInplaceAccessor
 
 
 # -------------------------------------------- LBM Kernel Creation -----------------------------------------------------
@@ -43,9 +45,9 @@ def createLBMKernel(collisionRule, inputField, outputField, accessor):
     return result
 
 
-def createStreamPullKernel(collisionRule, numpyField=None, srcFieldName="src", dstFieldName="dst",
-                           genericLayout='numpy', genericFieldType=np.float64,
-                           builtinPeriodicity=(False, False, False)):
+def createStreamPullCollideKernel(collisionRule, numpyField=None, srcFieldName="src", dstFieldName="dst",
+                                  genericLayout='numpy', genericFieldType=np.float64,
+                                  builtinPeriodicity=(False, False, False)):
     """
     Implements a stream-pull scheme, where values are read from source and written to destination field
     :param collisionRule: a collision rule created by lbm method
@@ -78,6 +80,46 @@ def createStreamPullKernel(collisionRule, numpyField=None, srcFieldName="src", d
     return createLBMKernel(collisionRule, src, dst, accessor)
 
 
+def createCollideOnlyKernel(collisionRule, numpyField=None, fieldName="src",
+                            genericLayout='numpy', genericFieldType=np.float64):
+    """
+    Implements a collision only (no neighbor access) LBM kernel.
+    For parameters see function ``createStreamPullCollideKernel``
+    """
+    dim = collisionRule.method.dim
+    if numpyField is not None:
+        assert len(numpyField.shape) == dim + 1, "Field dimension mismatch: dimension is %s, should be %d" % \
+                                                 (len(numpyField.shape), dim + 1)
+
+    if numpyField is None:
+        field = Field.createGeneric(fieldName, dim, indexDimensions=1, layout=genericLayout, dtype=genericFieldType)
+    else:
+        field = Field.createFromNumpyArray(fieldName, numpyField, indexDimensions=1)
+
+    return createLBMKernel(collisionRule, field, field, CollideOnlyInplaceAccessor)
+
+
+def createStreamPullOnlyKernel(stencil, numpyField=None, srcFieldName="src", dstFieldName="dst",
+                               genericLayout='numpy', genericFieldType=np.float64):
+    """
+    Creates a stream-pull kernel, without collision
+    For parameters see function ``createStreamPullCollideKernel``
+    """
+
+    dim = len(stencil[0])
+    if numpyField is None:
+        src = Field.createGeneric(srcFieldName, dim, indexDimensions=1, layout=genericLayout, dtype=genericFieldType)
+        dst = Field.createGeneric(dstFieldName, dim, indexDimensions=1, layout=genericLayout, dtype=genericFieldType)
+    else:
+        src = Field.createFromNumpyArray(srcFieldName, numpyField, indexDimensions=1)
+        dst = Field.createFromNumpyArray(dstFieldName, numpyField, indexDimensions=1)
+
+    accessor = StreamPullTwoFieldsAccessor()
+    eqs = [sp.Eq(a, b) for a, b in zip(accessor.write(dst, stencil), accessor.read(src, stencil))
+           if sp.Eq(a, b) != True]
+    return EquationCollection(eqs, [])
+
+
 # ---------------------------------- Pdf array creation for various layouts --------------------------------------------