From fd13c4cc4a1fc1f8e32efd403a7fbbd9509ef91a Mon Sep 17 00:00:00 2001
From: Daniel Bauer <daniel.j.bauer@fau.de>
Date: Fri, 4 Mar 2022 16:08:06 +0100
Subject: [PATCH] Non uniform grids with generated boundary conditions

---
 .../FlowAroundSphereCodeGen.py                |  9 +++++---
 .../additional_data_handler.py                |  4 ++++
 python/pystencils_walberla/boundary.py        |  8 +++++--
 .../templates/Boundary.tmpl.h                 | 23 ++++++++++++++++++-
 4 files changed, 38 insertions(+), 6 deletions(-)

diff --git a/apps/benchmarks/FlowAroundSphereCodeGen/FlowAroundSphereCodeGen.py b/apps/benchmarks/FlowAroundSphereCodeGen/FlowAroundSphereCodeGen.py
index 675a3766e..8d3b556bb 100644
--- a/apps/benchmarks/FlowAroundSphereCodeGen/FlowAroundSphereCodeGen.py
+++ b/apps/benchmarks/FlowAroundSphereCodeGen/FlowAroundSphereCodeGen.py
@@ -73,14 +73,17 @@ with CodeGeneration() as ctx:
 
     generate_alternating_lbm_boundary(ctx, 'FlowAroundSphereCodeGen_UBB', ubb, lb_method,
                                       target=target, streaming_pattern=streaming_pattern,
-                                      additional_data_handler=ubb_data_handler)
+                                      additional_data_handler=ubb_data_handler,
+                                      layout='fzyx')
 
     generate_alternating_lbm_boundary(ctx, 'FlowAroundSphereCodeGen_NoSlip', NoSlip(), lb_method,
-                                      target=target, streaming_pattern=streaming_pattern)
+                                      target=target, streaming_pattern=streaming_pattern,
+                                      layout='fzyx')
 
     generate_alternating_lbm_boundary(ctx, 'FlowAroundSphereCodeGen_Outflow', outflow, lb_method,
                                       target=target, streaming_pattern=streaming_pattern,
-                                      additional_data_handler=outflow_data_handler)
+                                      additional_data_handler=outflow_data_handler,
+                                      layout='fzyx')
 
     # communication
     generate_lb_pack_info(ctx, 'FlowAroundSphereCodeGen_PackInfo', stencil, pdfs,
diff --git a/python/pystencils_walberla/additional_data_handler.py b/python/pystencils_walberla/additional_data_handler.py
index 4d4d82d52..2c4efcc65 100644
--- a/python/pystencils_walberla/additional_data_handler.py
+++ b/python/pystencils_walberla/additional_data_handler.py
@@ -47,6 +47,10 @@ class AdditionalDataHandler:
     def stencil_info(self):
         return [(i, d, ", ".join([str(e) for e in d])) for i, d in enumerate(self._walberla_stencil)]
 
+    @property
+    def walberla_stencil(self):
+        return self._walberla_stencil
+
     @property
     def inverse_directions(self):
         inv_dirs = []
diff --git a/python/pystencils_walberla/boundary.py b/python/pystencils_walberla/boundary.py
index c8bf214de..5cf8f1fe5 100644
--- a/python/pystencils_walberla/boundary.py
+++ b/python/pystencils_walberla/boundary.py
@@ -6,6 +6,8 @@ from pystencils.boundaries.createindexlist import (
     boundary_index_array_coordinate_names, direction_member_name,
     numpy_data_type_for_boundary_object)
 from pystencils.data_types import TypedSymbol, create_type
+from pystencils.stencil import inverse_direction
+
 from pystencils_walberla.codegen import config_from_context
 from pystencils_walberla.jinja_filters import add_pystencils_filters_to_jinja_env
 from pystencils_walberla.additional_data_handler import AdditionalDataHandler
@@ -29,6 +31,7 @@ def generate_boundary(generation_context,
                       additional_data_handler=None,
                       interface_mappings=(),
                       generate_functor=True,
+                      layout='fzyx',
                       **create_kernel_params):
 
     if boundary_object.additional_data and additional_data_handler is None:
@@ -50,7 +53,7 @@ def generate_boundary(generation_context,
 
     field = Field.create_generic(field_name, dim,
                                  field_data_type,
-                                 index_dimensions=len(index_shape), layout='fzyx', index_shape=index_shape,
+                                 index_dimensions=len(index_shape), layout=layout, index_shape=index_shape,
                                  field_type=field_type)
 
     index_field = Field('indexVector', FieldType.INDEXED, index_struct_dtype, layout=[0],
@@ -88,7 +91,8 @@ def generate_boundary(generation_context,
         'namespace': namespace,
         'inner_or_boundary': boundary_object.inner_or_boundary,
         'single_link': boundary_object.single_link,
-        'additional_data_handler': additional_data_handler
+        'additional_data_handler': additional_data_handler,
+        'layout': layout
     }
 
     env = Environment(loader=PackageLoader('pystencils_walberla'), undefined=StrictUndefined)
diff --git a/python/pystencils_walberla/templates/Boundary.tmpl.h b/python/pystencils_walberla/templates/Boundary.tmpl.h
index cc5351d51..daf3e0c96 100644
--- a/python/pystencils_walberla/templates/Boundary.tmpl.h
+++ b/python/pystencils_walberla/templates/Boundary.tmpl.h
@@ -198,7 +198,27 @@ public:
         indexVectorOuter.clear();
 
         {% if inner_or_boundary -%}
-        for( auto it = flagField->begin(); it != flagField->end(); ++it )
+        {% if layout == "fzyx" -%}
+        {%- for dirIdx, dirVec, offset in additional_data_handler.stencil_info %}
+        for( auto it = flagField->beginWithGhostLayerXYZ( cell_idx_c( flagField->nrOfGhostLayers() - 1 ) ); it != flagField->end(); ++it )
+        {
+           if( ! isFlagSet(it, domainFlag) )
+              continue;
+
+           if ( isFlagSet( it.neighbor({{offset}} {%if dim == 3%}, 0 {%endif %}), boundaryFlag ) )
+           {
+              auto element = {{StructName}}(it.x(), it.y(), {%if dim == 3%} it.z(), {%endif %} {{dirIdx}} );
+              {{additional_data_handler.data_initialisation(dirIdx)|indent(16)}}
+              indexVectorAll.push_back( element );
+              if( inner.contains( it.x(), it.y(), it.z() ) )
+                 indexVectorInner.push_back( element );
+              else
+                 indexVectorOuter.push_back( element );
+           }
+        }
+        {% endfor %}
+        {%else%}
+        for( auto it = flagField->beginWithGhostLayerXYZ( cell_idx_c( flagField->nrOfGhostLayers() - 1 ) ); it != flagField->end(); ++it )
         {
             if( ! isFlagSet(it, domainFlag) )
                 continue;
@@ -215,6 +235,7 @@ public:
             }
             {% endfor %}
         }
+        {%endif%}
         {%else%}
         auto flagWithGLayers = flagField->xyzSizeWithGhostLayer();
         real_t dot = 0.0; real_t maxn = 0.0;
-- 
GitLab