Skip to content
Snippets Groups Projects
Commit 358ebb11 authored by Helen Schottenhamml's avatar Helen Schottenhamml
Browse files

Merge branch 'Fixdatatype' into 'master'

Fix data type for LB weight info

See merge request !109
parents 5bd88f01 b22c2846
No related branches found
No related tags found
1 merge request!109Fix data type for LB weight info
Pipeline #35761 passed
...@@ -283,7 +283,7 @@ class UBB(LbBoundary): ...@@ -283,7 +283,7 @@ class UBB(LbBoundary):
Returns: Returns:
list containing LbmWeightInfo and NeighbourOffsetArrays list containing LbmWeightInfo and NeighbourOffsetArrays
""" """
return [LbmWeightInfo(lb_method), NeighbourOffsetArrays(lb_method.stencil)] return [LbmWeightInfo(lb_method, self.data_type), NeighbourOffsetArrays(lb_method.stencil)]
@property @property
def velocity_is_callable(self): def velocity_is_callable(self):
...@@ -312,7 +312,8 @@ class UBB(LbBoundary): ...@@ -312,7 +312,8 @@ class UBB(LbBoundary):
velocity = [eq.rhs for eq in shifted_vel_eqs.new_filtered(cqc.first_order_moment_symbols).main_assignments] 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) c_s_sq = sp.Rational(1, 3)
weight_of_direction = LbmWeightInfo.weight_of_direction 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( 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) dir_symbol, lb_method)
...@@ -595,10 +596,11 @@ class DiffusionDirichlet(LbBoundary): ...@@ -595,10 +596,11 @@ class DiffusionDirichlet(LbBoundary):
name: optional name of the boundary. 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: if name is None:
name = "Diffusion Dirichlet " + str(concentration) name = "Diffusion Dirichlet " + str(concentration)
self.concentration = concentration self.concentration = concentration
self._data_type = data_type
super(DiffusionDirichlet, self).__init__(name) super(DiffusionDirichlet, self).__init__(name)
...@@ -611,10 +613,11 @@ class DiffusionDirichlet(LbBoundary): ...@@ -611,10 +613,11 @@ class DiffusionDirichlet(LbBoundary):
Returns: Returns:
list containing LbmWeightInfo list containing LbmWeightInfo
""" """
return [LbmWeightInfo(lb_method)] return [LbmWeightInfo(lb_method, self._data_type)]
def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field): def __call__(self, f_out, f_in, dir_symbol, inv_dir, lb_method, index_field):
w_dir = LbmWeightInfo.weight_of_direction(dir_symbol, lb_method) 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]), return [Assignment(f_in(inv_dir[dir_symbol]),
2 * w_dir * self.concentration - f_out(dir_symbol))] 2 * w_dir * self.concentration - f_out(dir_symbol))]
......
...@@ -38,7 +38,7 @@ class LatticeBoltzmannBoundaryHandling(BoundaryHandling): ...@@ -38,7 +38,7 @@ class LatticeBoltzmannBoundaryHandling(BoundaryHandling):
self._prev_timestep = None self._prev_timestep = None
def add_fixed_steps(self, fixed_loop, **kwargs): def add_fixed_steps(self, fixed_loop, **kwargs):
if self._inplace: # Fixed Loop can't do timestep selection if self._inplace: # Fixed Loop can't do timestep selection
raise NotImplementedError("Adding to fixed loop is currently not supported for inplace kernels") raise NotImplementedError("Adding to fixed loop is currently not supported for inplace kernels")
super(LatticeBoltzmannBoundaryHandling, self).add_fixed_steps(fixed_loop, **kwargs) super(LatticeBoltzmannBoundaryHandling, self).add_fixed_steps(fixed_loop, **kwargs)
...@@ -52,10 +52,12 @@ class LatticeBoltzmannBoundaryHandling(BoundaryHandling): ...@@ -52,10 +52,12 @@ class LatticeBoltzmannBoundaryHandling(BoundaryHandling):
if boundary_obj not in self._boundary_object_to_boundary_info: if boundary_obj not in self._boundary_object_to_boundary_info:
sym_index_field = Field.create_generic('indexField', spatial_dimensions=1, sym_index_field = Field.create_generic('indexField', spatial_dimensions=1,
dtype=numpy_data_type_for_boundary_object(boundary_obj, self.dim)) dtype=numpy_data_type_for_boundary_object(boundary_obj, self.dim))
kernels = [self._create_boundary_kernel(
self._data_handling.fields[self._field_name], sym_index_field, boundary_obj, Timestep.EVEN).compile(), ast_even = self._create_boundary_kernel(self._data_handling.fields[self._field_name], sym_index_field,
self._create_boundary_kernel( boundary_obj, Timestep.EVEN)
self._data_handling.fields[self._field_name], sym_index_field, boundary_obj, Timestep.ODD).compile()] ast_odd = self._create_boundary_kernel(self._data_handling.fields[self._field_name], sym_index_field,
boundary_obj, Timestep.ODD)
kernels = [ast_even.compile(), ast_odd.compile()]
if flag is None: if flag is None:
flag = self.flag_interface.reserve_next_flag() flag = self.flag_interface.reserve_next_flag()
boundary_info = self.InplaceStreamingBoundaryInfo(self, boundary_obj, flag, kernels) boundary_info = self.InplaceStreamingBoundaryInfo(self, boundary_obj, flag, kernels)
...@@ -84,6 +86,7 @@ class LatticeBoltzmannBoundaryHandling(BoundaryHandling): ...@@ -84,6 +86,7 @@ class LatticeBoltzmannBoundaryHandling(BoundaryHandling):
self.boundary_object = boundary_obj self.boundary_object = boundary_obj
self.flag = flag self.flag = flag
self._kernels = kernels self._kernels = kernels
# end class InplaceStreamingBoundaryInfo # end class InplaceStreamingBoundaryInfo
# ------------------------------ Force On Boundary ------------------------------------------------------------ # ------------------------------ Force On Boundary ------------------------------------------------------------
...@@ -148,29 +151,32 @@ class LatticeBoltzmannBoundaryHandling(BoundaryHandling): ...@@ -148,29 +151,32 @@ class LatticeBoltzmannBoundaryHandling(BoundaryHandling):
return dh.reduce_float_sequence(list(result), 'sum') return dh.reduce_float_sequence(list(result), 'sum')
# end class LatticeBoltzmannBoundaryHandling # end class LatticeBoltzmannBoundaryHandling
class LbmWeightInfo(CustomCodeNode): class LbmWeightInfo(CustomCodeNode):
def __init__(self, lb_method, data_type='double'):
self.weights_symbol = TypedSymbol("weights", data_type)
data_type_string = "double" if self.weights_symbol.dtype.numpy_dtype == np.float64 else "float"
# --------------------------- Functions to be used by boundaries -------------------------- weights = [str(w.evalf()) for w in lb_method.weights]
if data_type_string == "float":
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"
super(LbmWeightInfo, self).__init__(code, symbols_read=set(), symbols_defined={w_sym})
@staticmethod def weight_of_direction(self, dir_idx, lb_method=None):
def weight_of_direction(dir_idx, lb_method=None):
if isinstance(sp.sympify(dir_idx), sp.Integer): if isinstance(sp.sympify(dir_idx), sp.Integer):
return lb_method.weights[dir_idx].evalf() return lb_method.weights[dir_idx].evalf()
else: else:
return sp.IndexedBase(LbmWeightInfo.WEIGHTS_SYMBOL, shape=(1,))[dir_idx] return sp.IndexedBase(self.weights_symbol, shape=(1,))[dir_idx]
# ---------------------------------- Internal ---------------------------------------------
WEIGHTS_SYMBOL = TypedSymbol("weights", "double")
def __init__(self, lb_method):
weights = [str(w.evalf()) for w in lb_method.weights]
w_sym = LbmWeightInfo.WEIGHTS_SYMBOL
code = "const double %s [] = { %s };\n" % (w_sym.name, ",".join(weights))
super(LbmWeightInfo, self).__init__(code, symbols_read=set(), symbols_defined={w_sym})
# end class LbmWeightInfo # end class LbmWeightInfo
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment