diff --git a/lbmpy/boundaries/boundaryconditions.py b/lbmpy/boundaries/boundaryconditions.py
index b82b745e796815d6aa1038a18a34f92c6fee64a7..29e20bc4c0ff3969da4f4ef47c1678396dcc0019 100644
--- a/lbmpy/boundaries/boundaryconditions.py
+++ b/lbmpy/boundaries/boundaryconditions.py
@@ -26,9 +26,8 @@ class LbBoundary(abc.ABC):
def __init__(self, name=None):
self._name = name
- self._data_type = None
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
"""
This function defines the boundary behavior and must therefore be implemented by all boundaries.
The boundary is defined through a list of sympy equations from which a boundary kernel is generated.
@@ -108,7 +107,7 @@ class NoSlip(LbBoundary):
"""Set an optional name here, to mark boundaries, for example for force evaluations"""
super(NoSlip, self).__init__(name)
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
return Assignment(f_in(inv_dir[dir_symbol]), f_out(dir_symbol))
@@ -213,7 +212,7 @@ class FreeSlip(LbBoundary):
else:
return []
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
if self.normal_direction:
normal_direction = self.normal_direction
mirrored_stencil_symbol = MirroredStencilDirections._mirrored_symbol(self.mirror_axis)
@@ -246,7 +245,7 @@ class UBB(LbBoundary):
name: optional name of the boundary.
"""
- def __init__(self, velocity, adapt_velocity_to_force=False, dim=None, name=None, data_type=None):
+ def __init__(self, velocity, adapt_velocity_to_force=False, dim=None, name=None, data_type='double'):
self._velocity = velocity
self._adaptVelocityToForce = adapt_velocity_to_force
if callable(self._velocity) and not dim:
@@ -254,6 +253,7 @@ class UBB(LbBoundary):
elif not callable(self._velocity):
dim = len(velocity)
self.dim = dim
+ self.data_type = data_type
super(UBB, self).__init__(name)
@@ -262,7 +262,7 @@ class UBB(LbBoundary):
""" In case of the UBB boundary additional data is a velocity vector. This vector is added to each cell to
realize velocity profiles for the inlet."""
if self.velocity_is_callable:
- return [(f'vel_{i}', create_type(self._data_type)) for i in range(self.dim)]
+ return [(f'vel_{i}', create_type(self.data_type)) for i in range(self.dim)]
else:
return []
@@ -283,7 +283,7 @@ class UBB(LbBoundary):
Returns:
list containing LbmWeightInfo and NeighbourOffsetArrays
"""
- return [LbmWeightInfo(lb_method, data_type=self._data_type), NeighbourOffsetArrays(lb_method.stencil)]
+ return [LbmWeightInfo(lb_method, self.data_type), NeighbourOffsetArrays(lb_method.stencil)]
@property
def velocity_is_callable(self):
@@ -291,11 +291,7 @@ class UBB(LbBoundary):
This is useful if the inflow velocity should have a certain profile for instance"""
return callable(self._velocity)
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
- if data_type is None:
- self._data_type = "double"
- else:
- self._data_type = data_type
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
vel_from_idx_field = callable(self._velocity)
vel = [index_field(f'vel_{i}') for i in range(self.dim)] if vel_from_idx_field else self._velocity
@@ -316,7 +312,7 @@ class UBB(LbBoundary):
velocity = [eq.rhs for eq in shifted_vel_eqs.new_filtered(cqc.first_order_moment_symbols).main_assignments]
c_s_sq = sp.Rational(1, 3)
- weight_info = LbmWeightInfo(lb_method, data_type=self._data_type)
+ weight_info = LbmWeightInfo(lb_method, data_type=self.data_type)
weight_of_direction = weight_info.weight_of_direction
vel_term = 2 / c_s_sq * sum([d_i * v_i for d_i, v_i in zip(neighbor_offset, velocity)]) * weight_of_direction(
dir_symbol, lb_method)
@@ -380,7 +376,7 @@ class SimpleExtrapolationOutflow(LbBoundary):
"""
return [NeighbourOffsetArrays(lb_method.stencil)]
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
neighbor_offset = NeighbourOffsetArrays.neighbour_offset(dir_symbol, lb_method.stencil)
tangential_offset = tuple(offset - normal for offset, normal in zip(neighbor_offset, self.normal_direction))
@@ -420,7 +416,7 @@ class ExtrapolationOutflow(LbBoundary):
def __init__(self, normal_direction, lb_method, dt=1, dx=1, name=None,
streaming_pattern='pull', zeroth_timestep=Timestep.BOTH,
- initial_density=None, initial_velocity=None, data_type=None):
+ initial_density=None, initial_velocity=None, data_type='double'):
self.lb_method = lb_method
self.stencil = lb_method.stencil
@@ -443,6 +439,8 @@ class ExtrapolationOutflow(LbBoundary):
self.initial_velocity = initial_velocity
self.equilibrium_calculation = None
+ self.data_type = data_type
+
if initial_density and initial_velocity:
equilibrium = lb_method.get_equilibrium(conserved_quantity_equations=AssignmentCollection([]))
rho = lb_method.zeroth_order_equilibrium_moment_symbol
@@ -489,7 +487,7 @@ class ExtrapolationOutflow(LbBoundary):
def additional_data(self):
"""Used internally only. For the ExtrapolationOutflow information of the previous PDF values is needed. This
information is stored in the index vector."""
- data = [('pdf', create_type(self._data_type)), ('pdf_nd', create_type(self._data_type))]
+ data = [('pdf', create_type(self.data_type)), ('pdf_nd', create_type(self.data_type))]
return data
@property
@@ -508,11 +506,7 @@ class ExtrapolationOutflow(LbBoundary):
"""
return [NeighbourOffsetArrays(lb_method.stencil)]
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
- if data_type is None:
- self._data_type = "double"
- else:
- self._data_type = data_type
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
subexpressions = []
boundary_assignments = []
dtdx = sp.Rational(self.dt, self.dx)
@@ -555,7 +549,7 @@ class FixedDensity(LbBoundary):
super(FixedDensity, self).__init__(name)
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
def remove_asymmetric_part_of_main_assignments(assignment_collection, degrees_of_freedom):
new_main_assignments = [Assignment(a.lhs, get_symmetric_part(a.rhs, degrees_of_freedom))
for a in assignment_collection.main_assignments]
@@ -602,10 +596,11 @@ class DiffusionDirichlet(LbBoundary):
name: optional name of the boundary.
"""
- def __init__(self, concentration, name=None):
+ def __init__(self, concentration, name=None, data_type='double'):
if name is None:
name = "Diffusion Dirichlet " + str(concentration)
self.concentration = concentration
+ self._data_type = data_type
super(DiffusionDirichlet, self).__init__(name)
@@ -620,15 +615,11 @@ class DiffusionDirichlet(LbBoundary):
"""
return [LbmWeightInfo(lb_method, self._data_type)]
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
- if data_type is None:
- self._data_type = "double"
- else:
- self._data_type = data_type
- weight_info = LbmWeightInfo(lb_method, data_type=self._data_type)
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
+ weight_info = LbmWeightInfo(lb_method, self._data_type)
w_dir = weight_info.weight_of_direction(dir_symbol, lb_method)
-
- return [Assignment(f_in(inv_dir[dir_symbol]), 2 * w_dir * self.concentration - f_out(dir_symbol))]
+ return [Assignment(f_in(inv_dir[dir_symbol]),
+ 2 * w_dir * self.concentration - f_out(dir_symbol))]
# end class DiffusionDirichlet
@@ -649,7 +640,7 @@ class NeumannByCopy(LbBoundary):
"""
return [NeighbourOffsetArrays(lb_method.stencil)]
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
neighbour_offset = NeighbourOffsetArrays.neighbour_offset(dir_symbol, lb_method.stencil)
return [Assignment(f_in(inv_dir[dir_symbol]), f_out(inv_dir[dir_symbol])),
Assignment(f_out[neighbour_offset](dir_symbol), f_out(dir_symbol))]
@@ -682,7 +673,7 @@ class StreamInConstant(LbBoundary):
"""
return [NeighbourOffsetArrays(lb_method.stencil)]
- def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field, data_type=None):
+ def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
neighbour_offset = NeighbourOffsetArrays.neighbour_offset(dir_symbol, lb_method.stencil)
return [Assignment(f_in(inv_dir[dir_symbol]), self.constant),
Assignment(f_out[neighbour_offset](dir_symbol), self.constant)]
diff --git a/lbmpy/boundaries/boundaryhandling.py b/lbmpy/boundaries/boundaryhandling.py
index 1ca7b972a76e3b704d2ba3d0ed7f5efea8a96913..b840c6c708c0fd0055bb05798160d4392bcf3fce 100644
--- a/lbmpy/boundaries/boundaryhandling.py
+++ b/lbmpy/boundaries/boundaryhandling.py
@@ -161,11 +161,12 @@ class LbmWeightInfo(CustomCodeNode):
weights = [str(w.evalf()) for w in lb_method.weights]
if data_type_string == "float":
- weights = ".f, ".join(weights)
+ weights = "f, ".join(weights)
+ weights += "f" # suffix for the last element
else:
weights = ", ".join(weights)
w_sym = self.weights_symbol
- code = f"const {data_type_string} {w_sym.name} [] = { weights };\n"
+ code = f"const {data_type_string} {w_sym.name} [] = {{{ weights }}};\n"
super(LbmWeightInfo, self).__init__(code, symbols_read=set(), symbols_defined={w_sym})
def weight_of_direction(self, dir_idx, lb_method=None):
@@ -181,6 +182,7 @@ class LbmWeightInfo(CustomCodeNode):
def create_lattice_boltzmann_boundary_kernel(pdf_field, index_field, lb_method, boundary_functor,
prev_timestep=Timestep.BOTH, streaming_pattern='pull',
target=Target.CPU, **kernel_creation_args):
+
indexing = BetweenTimestepsIndexing(
pdf_field, lb_method.stencil, prev_timestep, streaming_pattern, np.int64, np.int64)
@@ -188,8 +190,7 @@ def create_lattice_boltzmann_boundary_kernel(pdf_field, index_field, lb_method,
dir_symbol = indexing.dir_symbol
inv_dir = indexing.inverse_dir_symbol
- boundary_assignments = boundary_functor(f_out, f_in, dir_symbol, inv_dir, lb_method, index_field,
- data_type=pdf_field.dtype.numpy_dtype)
+ boundary_assignments = boundary_functor(f_out, f_in, dir_symbol, inv_dir, lb_method, index_field)
boundary_assignments = indexing.substitute_proxies(boundary_assignments)
# Code Elements inside the loop