Commit 606cd79d authored by Frederik Hennig's avatar Frederik Hennig
Browse files

In-Place Streaming Patterns in Lb Codegen Test

parent 42cc89ac
Pipeline #29327 failed with stages
in 19 minutes and 23 seconds
from .boundary import generate_boundary
from .walberla_lbm_generation import RefinementScaling, generate_lattice_model
from .packinfo import generate_lb_pack_info
__all__ = ['generate_lattice_model', 'RefinementScaling', 'generate_boundary']
__all__ = ['generate_lattice_model', 'RefinementScaling', 'generate_boundary', 'generate_lb_pack_info']
......@@ -75,7 +75,7 @@ class OutflowAdditionalDataHandler(AdditionalDataHandler):
@property
def additional_field_data(self):
return f"auto {self._field_name} = block->getData< field::GhostLayerField<double, "\
f"{len(self._stencil)}> >({self._field_name}ID); "
f"{len(self._stencil)}> >({self._field_name}ID); "
@property
def data_initialisation(self):
......
from typing import OrderedDict
import pystencils_walberla.boundary
from lbmpy.boundaries.boundaryhandling import create_lattice_boltzmann_boundary_kernel
from lbmpy.advanced_streaming import Timestep, get_timesteps, is_inplace
from lbmpy.advanced_streaming import Timestep, is_inplace
from lbmpy.boundaries import ExtrapolationOutflow, UBB
from lbmpy_walberla.additional_data_handler import UBBAdditionalDataHandler, OutflowAdditionalDataHandler
......@@ -25,7 +24,7 @@ def generate_boundary(generation_context,
**kwargs),
'OddSweep': create_lattice_boltzmann_boundary_kernel(*pargs,
streaming_pattern=streaming_pattern,
prev_timestep=Timestep.EVEN,
prev_timestep=Timestep.ODD,
**kwargs)
}
else:
......
......@@ -6,14 +6,14 @@ from pystencils_walberla.codegen import comm_directions, generate_pack_info
from pystencils import Assignment, Field
def generate_pack_infos_for_lbm_kernel(generation_context,
class_name_prefix: str,
stencil,
pdf_field,
streaming_pattern='pull',
lb_collision_rule=None,
always_generate_seperate_classes=False,
**create_kernel_params):
def generate_lb_pack_info(generation_context,
class_name_prefix: str,
stencil,
pdf_field,
streaming_pattern='pull',
lb_collision_rule=None,
always_generate_seperate_classes=False,
**create_kernel_params):
"""Generates waLBerla MPI PackInfos for an LBM kernel, based on a given method
and streaming pattern. For in-place streaming patterns, two PackInfos are generated;
one for the even and another for the odd time steps.
......@@ -27,14 +27,14 @@ def generate_pack_infos_for_lbm_kernel(generation_context,
include_non_pdf_fields: Whether to scan for and include non-local accesses to fields
apart from the PDF field in the PackInfo. If set to True, the
lb_collision_rule parameter also needs to be specified.
lb_collision_rule: Optional. The collision rule defining the lattice boltzmann kernel, as returned
lb_collision_rule: Optional. The collision rule defining the lattice boltzmann kernel, as returned
by `create_lb_collision_rule`. If specified, it will be scanned for non-local
accesses to other fields other than the PDF fields (as might be required for
computing gradients in coupled simulations), whose communication will then
be included in the PackInfo.
always_generate_seperate_classes: If True, generate a pair of Even/Odd PackInfos even for a two-
fields kernel (i.e. the pull/push patterns). Otherwise, for two-fields
kernels, only one PackInfo class will be generated without a
kernels, only one PackInfo class will be generated without a
suffix to its name.
**create_kernel_params: remaining keyword arguments are passed to `pystencils.create_kernel`
"""
......@@ -70,13 +70,13 @@ def generate_pack_infos_for_lbm_kernel(generation_context,
d = stencil.index(streaming_dir)
fa = write_accesses[d]
spec[(comm_dir,)].add(fa)
if t == Timestep.EVEN:
class_name_suffix = 'Even'
elif t == Timestep.ODD:
class_name_suffix = 'Odd'
else:
class_name_suffix = ''
class_name = class_name_prefix + class_name_suffix
generate_pack_info(generation_context, class_name, spec, **create_kernel_params)
generate_pack_info(generation_context, class_name, spec, namespace='lbm', **create_kernel_params)
......@@ -68,7 +68,7 @@ def generate_boundary(generation_context,
kernel.function_name = "boundary_" + boundary_object.name
kernel.assumed_inner_stride_one = False
kernel_info = KernelInfo(kernel)
sweep_to_kernel_info_dict = { '' : kernel_info }
sweep_to_kernel_info_dict = {'': kernel_info}
dummy_kernel_info = kernel_info
# waLBerla is a 3D framework. Therefore, a zero for the z index has to be added if we work in 2D
......@@ -94,7 +94,7 @@ def generate_boundary(generation_context,
'class_name': boundary_object.name,
'sweep_classes': sweep_to_kernel_info_dict,
'multi_sweep': multi_sweep,
'dummy_kernel_info' : dummy_kernel_info,
'dummy_kernel_info': dummy_kernel_info,
'StructName': struct_name,
'StructDeclaration': struct_from_numpy_dtype(struct_name, index_struct_dtype),
'stencil_info': stencil_info,
......
......@@ -401,6 +401,7 @@ def generate_destructor(kernel_info, class_name):
contents += delete_loop.format(original_field_name=field_name[:-len('_tmp')])
return temporary_constructor.format(contents=contents, class_name=class_name)
@jinja2.contextfilter
def nested_class_method_definition_prefix(ctx, nested_class_name):
outer_class = ctx['class_name']
......
......@@ -80,13 +80,14 @@ add_subdirectory(codegen)
waLBerla_generate_target_from_python(NAME ChannelFlowCodeGenGenerated
FILE codegen/ChannelFlowCodeGen.py
OUT_FILES ChannelFlowCodeGen_Sweep.cpp ChannelFlowCodeGen_Sweep.h
OUT_FILES ChannelFlowCodeGen_EvenSweep.cpp ChannelFlowCodeGen_EvenSweep.h
ChannelFlowCodeGen_OddSweep.cpp ChannelFlowCodeGen_OddSweep.h
ChannelFlowCodeGen_MacroSetter.cpp ChannelFlowCodeGen_MacroSetter.h
ChannelFlowCodeGen_MacroGetter.cpp ChannelFlowCodeGen_MacroGetter.h
ChannelFlowCodeGen_UBB.cpp ChannelFlowCodeGen_UBB.h
ChannelFlowCodeGen_NoSlip.cpp ChannelFlowCodeGen_NoSlip.h
ChannelFlowCodeGen_Outflow.cpp ChannelFlowCodeGen_Outflow.h
ChannelFlowCodeGen_PackInfo.cpp ChannelFlowCodeGen_PackInfo.h
ChannelFlowCodeGen_PackInfoEven.cpp ChannelFlowCodeGen_PackInfoEven.h
ChannelFlowCodeGen_PackInfoOdd.cpp ChannelFlowCodeGen_PackInfoOdd.h
ChannelFlowCodeGen_InfoHeader.h)
waLBerla_compile_test( FILES codegen/ChannelFlowCodeGen.cpp DEPENDS ChannelFlowCodeGenGenerated)
set_target_properties( ChannelFlowCodeGen PROPERTIES CXX_VISIBILITY_PRESET hidden)
......
......@@ -29,15 +29,18 @@
// CodeGen includes
#include "ChannelFlowCodeGen_InfoHeader.h"
#include "ChannelFlowCodeGen_MacroGetter.h"
// #include "ChannelFlowCodeGen_MacroGetter.h"
#include "ChannelFlowCodeGen_MacroSetter.h"
#include "ChannelFlowCodeGen_NoSlip.h"
#include "ChannelFlowCodeGen_Outflow.h"
#include "ChannelFlowCodeGen_PackInfo.h"
#include "ChannelFlowCodeGen_Sweep.h"
#include "ChannelFlowCodeGen_PackInfoEven.h"
#include "ChannelFlowCodeGen_PackInfoOdd.h"
#include "ChannelFlowCodeGen_EvenSweep.h"
#include "ChannelFlowCodeGen_OddSweep.h"
#include "ChannelFlowCodeGen_UBB.h"
typedef pystencils::ChannelFlowCodeGen_PackInfo PackInfo_T;
typedef lbm::ChannelFlowCodeGen_PackInfoEven PackInfoEven_T;
typedef lbm::ChannelFlowCodeGen_PackInfoOdd PackInfoOdd_T;
typedef walberla::uint8_t flag_t;
typedef FlagField< flag_t > FlagField_T;
......@@ -69,6 +72,57 @@ auto VelocityCallback = [](const Cell &pos, const shared_ptr<StructuredBlockFore
return result;
};
class TimestepModulusTracker{
private:
uint_t modulus_;
public:
TimestepModulusTracker(uint_t initialTimestep) : modulus_(initialTimestep & 1) {};
void setTimestep(uint_t timestep) { modulus_ = timestep & 1; }
std::function<void()> advancementFunction() {
return [this] () {
this->modulus_ = (this->modulus_ + 1) & 1;
};
}
uint_t modulus() const { return modulus_; }
bool evenStep() const { return static_cast<bool>(modulus_ ^ 1); }
bool oddStep() const { return static_cast<bool>(modulus_ & 1); }
};
class AlternatingSweep{
public:
typedef std::function< void (IBlock *) > SweepFunction;
AlternatingSweep(SweepFunction evenSweep, SweepFunction oddSweep, std::shared_ptr<TimestepModulusTracker> tracker)
: tracker_(tracker), sweeps_{ evenSweep, oddSweep } {};
void operator() (IBlock * block) {
sweeps_[tracker_->modulus()](block);
}
private:
std::shared_ptr<TimestepModulusTracker> tracker_;
std::vector< SweepFunction > sweeps_;
};
class AlternatingBeforeFunction{
public:
typedef std::function< void () > BeforeFunction;
AlternatingBeforeFunction(BeforeFunction evenFunc, BeforeFunction oddFunc, std::shared_ptr<TimestepModulusTracker> &tracker)
: tracker_(tracker), funcs_{ evenFunc, oddFunc } {};
void operator() () {
funcs_[tracker_->modulus()]();
}
private:
std::shared_ptr<TimestepModulusTracker> tracker_;
std::vector< BeforeFunction > funcs_;
};
int main(int argc, char** argv)
{
......@@ -106,10 +160,12 @@ int main(int argc, char** argv)
for (auto& block : *blocks)
setterSweep(&block);
// setter sweep only initializes interior of domain - for push schemes to work a first communication is required here
blockforest::communication::UniformBufferedScheme< Stencil_T > initialComm(blocks);
initialComm.addPackInfo(make_shared< field::communication::PackInfo< PdfField_T > >(pdfFieldID));
initialComm();
// Create communication
blockforest::communication::UniformBufferedScheme< Stencil_T > evenComm(blocks);
evenComm.addPackInfo(make_shared< PackInfoEven_T >(pdfFieldID));
blockforest::communication::UniformBufferedScheme< Stencil_T > oddComm(blocks);
oddComm.addPackInfo(make_shared< PackInfoOdd_T >(pdfFieldID));
// create and initialize boundary handling
const FlagUID fluidFlagUID("Fluid");
......@@ -133,16 +189,24 @@ int main(int argc, char** argv)
// create time loop
SweepTimeloop timeloop(blocks->getBlockStorage(), timesteps);
// create communication for PdfField
blockforest::communication::UniformBufferedScheme< Stencil_T > communication(blocks);
communication.addPackInfo(make_shared< PackInfo_T >(pdfFieldID));
pystencils::ChannelFlowCodeGen_EvenSweep LBEvenSweep(densityFieldID, pdfFieldID, velFieldID, omega);
pystencils::ChannelFlowCodeGen_OddSweep LBOddSweep(densityFieldID, pdfFieldID, velFieldID, omega);
// All the sweeps
auto tracker = make_shared<TimestepModulusTracker>(0);
AlternatingSweep LBSweep(LBEvenSweep, LBOddSweep, tracker);
AlternatingSweep noSlipSweep(noSlip.getEvenSweep(), noSlip.getOddSweep(), tracker);
AlternatingSweep outflowSweep(outflow.getEvenSweep(), outflow.getOddSweep(), tracker);
AlternatingSweep ubbSweep(ubb.getEvenSweep(), ubb.getOddSweep(), tracker);
AlternatingBeforeFunction communication(evenComm, oddComm, tracker);
pystencils::ChannelFlowCodeGen_Sweep LBSweep(pdfFieldID, omega);
// add LBM sweep and communication to time loop
timeloop.add() << Sweep(noSlip, "noSlip boundary");
timeloop.add() << Sweep(outflow, "outflow boundary");
timeloop.add() << Sweep(ubb, "ubb boundary");
timeloop.add() << Sweep(noSlipSweep, "noSlip boundary");
timeloop.add() << Sweep(outflowSweep, "outflow boundary");
timeloop.add() << Sweep(ubbSweep, "ubb boundary");
timeloop.add() << BeforeFunction(communication, "communication")
<< BeforeFunction(tracker->advancementFunction(), "Timestep Advancement")
<< Sweep(LBSweep, "LB update rule");
// LBM stability check
......@@ -156,7 +220,7 @@ int main(int argc, char** argv)
"remaining time logger");
// add VTK output to time loop
pystencils::ChannelFlowCodeGen_MacroGetter getterSweep(densityFieldID, pdfFieldID, velFieldID);
// pystencils::ChannelFlowCodeGen_MacroGetter getterSweep(densityFieldID, pdfFieldID, velFieldID);
// VTK
uint_t vtkWriteFrequency = parameters.getParameter< uint_t >("vtkWriteFrequency", 0);
if (vtkWriteFrequency > 0)
......@@ -169,10 +233,10 @@ int main(int argc, char** argv)
vtkOutput->addCellDataWriter(velWriter);
vtkOutput->addCellDataWriter(densityWriter);
vtkOutput->addBeforeFunction([&]() {
for (auto& block : *blocks)
getterSweep(&block);
});
// vtkOutput->addBeforeFunction([&]() {
// for (auto& block : *blocks)
// getterSweep(&block);
// });
timeloop.addFuncAfterTimeStep(vtk::writeFiles(vtkOutput), "VTK Output");
}
......
from lbmpy.advanced_streaming.utility import get_timesteps, Timestep
from pystencils.field import fields
from lbmpy.macroscopic_value_kernels import macroscopic_values_setter, macroscopic_values_getter
from lbmpy.stencils import get_stencil
from lbmpy.creationfunctions import create_lb_method, create_lb_update_rule
from lbmpy.creationfunctions import create_lb_collision_rule, create_lb_method, create_lb_update_rule
from lbmpy.boundaries import NoSlip, UBB, ExtrapolationOutflow
from pystencils_walberla import CodeGeneration, generate_sweep, generate_pack_info_from_kernel
from lbmpy_walberla import RefinementScaling, generate_boundary, generate_lattice_model
from pystencils_walberla import CodeGeneration, generate_sweep
from lbmpy_walberla import RefinementScaling, generate_boundary, generate_lb_pack_info
import sympy as sp
......@@ -12,17 +13,26 @@ import sympy as sp
stencil = get_stencil("D3Q27")
q = len(stencil)
dim = len(stencil[0])
streaming_pattern = 'aa'
timesteps = get_timesteps(streaming_pattern)
pdfs, velocity_field, density_field = fields(f"pdfs({q}), velocity({dim}), density(1) : double[{dim}D]", layout='fzyx')
omega = sp.Symbol("omega")
u_max = sp.Symbol("u_max")
output = {
'density': density_field,
'velocity': velocity_field
}
options = {'method': 'cumulant',
'stencil': stencil,
'relaxation_rate': omega,
'galilean_correction': True,
'field_name': 'pdfs',
'temporary_field_name': 'pdfs_tmp',
# 'temporary_field_name': 'pdfs_tmp',
'streaming_pattern': streaming_pattern,
'output': output,
'optimization': {'symbolic_field': pdfs,
'cse_global': False,
'cse_pdfs': False}}
......@@ -31,13 +41,17 @@ method = create_lb_method(**options)
# getter & setter
setter_assignments = macroscopic_values_setter(method, velocity=velocity_field.center_vector,
pdfs=pdfs.center_vector, density=1)
getter_assignments = macroscopic_values_getter(method, velocity=velocity_field.center_vector,
pdfs=pdfs.center_vector, density=density_field)
pdfs=pdfs, density=1,
streaming_pattern=streaming_pattern,
previous_timestep=timesteps[0])
# getter_assignments = macroscopic_values_getter(method, velocity=velocity_field.center_vector,
# pdfs=pdfs.center_vector, density=density_field)
# opt = {'instruction_set': 'sse', 'assume_aligned': True, 'nontemporal': False, 'assume_inner_stride_one': True}
update_rule = create_lb_update_rule(lb_method=method, **options)
collision_rule = create_lb_collision_rule(lb_method=method, **options)
update_rule_even = create_lb_update_rule(collision_rule=collision_rule, timestep=Timestep.EVEN, **options)
update_rule_odd = create_lb_update_rule(collision_rule=collision_rule, timestep=Timestep.ODD, **options)
info_header = f"""
using namespace walberla;
......@@ -52,18 +66,23 @@ with CodeGeneration() as ctx:
target = 'cpu'
# sweeps
generate_sweep(ctx, 'ChannelFlowCodeGen_Sweep', update_rule, field_swaps=[('pdfs', 'pdfs_tmp')], target=target)
generate_sweep(ctx, 'ChannelFlowCodeGen_EvenSweep', update_rule_even, target=target)
generate_sweep(ctx, 'ChannelFlowCodeGen_OddSweep', update_rule_odd, target=target)
generate_sweep(ctx, 'ChannelFlowCodeGen_MacroSetter', setter_assignments, target=target)
generate_sweep(ctx, 'ChannelFlowCodeGen_MacroGetter', getter_assignments, target=target)
# generate_sweep(ctx, 'ChannelFlowCodeGen_MacroGetter', getter_assignments, target=target)
# boundaries
generate_boundary(ctx, 'ChannelFlowCodeGen_UBB', UBB(lambda *args: None, dim=dim), method, target=target)
generate_boundary(ctx, 'ChannelFlowCodeGen_NoSlip', NoSlip(), method, target=target)
generate_boundary(ctx, 'ChannelFlowCodeGen_UBB', UBB(lambda *args: None, dim=dim), method,
target=target, streaming_pattern=streaming_pattern, always_generate_seperate_sweeps=True)
generate_boundary(ctx, 'ChannelFlowCodeGen_NoSlip', NoSlip(), method, target=target,
streaming_pattern=streaming_pattern, always_generate_seperate_sweeps=True)
outflow = ExtrapolationOutflow(stencil[4], method)
generate_boundary(ctx, 'ChannelFlowCodeGen_Outflow', outflow, method, target=target)
generate_boundary(ctx, 'ChannelFlowCodeGen_Outflow', outflow, method, target=target,
streaming_pattern=streaming_pattern, always_generate_seperate_sweeps=True)
# communication
generate_pack_info_from_kernel(ctx, 'ChannelFlowCodeGen_PackInfo', update_rule, target=target)
generate_lb_pack_info(ctx, 'ChannelFlowCodeGen_PackInfo', stencil, pdfs,
streaming_pattern=streaming_pattern, always_generate_seperate_classes=True)
# Info header containing correct template definitions for stencil and field
ctx.write_file("ChannelFlowCodeGen_InfoHeader.h", info_header)
......@@ -9,6 +9,7 @@ class Scenario:
self.cells = (128, 32, 32)
self.blocks = (2, 2, 2)
# self.blocks = (1, 1, 1)
self.periodic = (0, 0, 0)
self.diameter_sphere = 30
......@@ -43,7 +44,7 @@ class Scenario:
{'direction': 'N', 'walldistance': -1, 'flag': 'NoSlip'},
{'direction': 'S', 'walldistance': -1, 'flag': 'NoSlip'},
{'direction': 'W', 'walldistance': -1, 'flag': 'UBB'},
{'direction': 'E', 'walldistance': -1, 'flag': 'Outflow'},
{'direction': 'E', 'walldistance': 0, 'flag': 'Outflow'},
{'direction': 'T', 'walldistance': -1, 'flag': 'NoSlip'},
{'direction': 'B', 'walldistance': -1, 'flag': 'NoSlip'},
],
......
......@@ -2,7 +2,7 @@ from lbmpy.creationfunctions import create_lb_collision_rule, create_lb_update_r
from lbmpy.advanced_streaming import Timestep
from lbmpy.stencils import get_stencil
from pystencils_walberla import CodeGeneration, generate_pack_info_from_kernel
from lbmpy_walberla.packinfo import generate_pack_infos_for_lbm_kernel
from lbmpy_walberla.packinfo import generate_lb_pack_info
from pystencils.field import Field
with CodeGeneration() as ctx:
......@@ -19,7 +19,7 @@ with CodeGeneration() as ctx:
}
# Generate PackInfo specifically for streaming pattern
generate_pack_infos_for_lbm_kernel(ctx, 'AccessorBasedPackInfo', stencil, pdf_field,
generate_lb_pack_info(ctx, 'AccessorBasedPackInfo', stencil, pdf_field,
streaming_pattern=streaming_pattern, target=target)
# Generate reference using the alternating pull/push approach
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment