Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found
Select Git revision

Target

Select target project
  • ravi.k.ayyala/lbmpy
  • brendan-waters/lbmpy
  • anirudh.jonnalagadda/lbmpy
  • jbadwaik/lbmpy
  • alexander.reinauer/lbmpy
  • itischler/lbmpy
  • he66coqe/lbmpy
  • ev81oxyl/lbmpy
  • Bindgen/lbmpy
  • da15siwa/lbmpy
  • holzer/lbmpy
  • RudolfWeeber/lbmpy
  • pycodegen/lbmpy
13 results
Select Git revision
Show changes
Showing
with 172 additions and 5 deletions
File moved
import numpy as np
import pytest
from lbmpy.boundaries import NoSlip, UBB, SimpleExtrapolationOutflow, ExtrapolationOutflow, \
FixedDensity, DiffusionDirichlet, NeumannByCopy, StreamInConstant, FreeSlip
from lbmpy.boundaries.boundaryhandling import LatticeBoltzmannBoundaryHandling
from lbmpy.boundaries import (NoSlip, NoSlipLinearBouzidi, QuadraticBounceBack,
UBB, SimpleExtrapolationOutflow, ExtrapolationOutflow, FixedDensity, DiffusionDirichlet,
NeumannByCopy, StreamInConstant, FreeSlip)
from lbmpy.boundaries.boundaryhandling import LatticeBoltzmannBoundaryHandling, create_lattice_boltzmann_boundary_kernel
from lbmpy.creationfunctions import create_lb_function, create_lb_method, LBMConfig
from lbmpy.enums import Stencil, Method
from lbmpy.geometry import add_box_boundary
from lbmpy.lbstep import LatticeBoltzmannStep
from lbmpy.stencils import LBStencil
import pystencils as ps
from pystencils import create_data_handling, make_slice, Target, CreateKernelConfig
from pystencils.slicing import slice_from_direction
from pystencils.stencil import inverse_direction
......@@ -26,7 +28,7 @@ def mirror_stencil(direction, mirror_axis):
def test_simple(target):
if target == Target.GPU:
import pytest
pytest.importorskip('pycuda')
pytest.importorskip('cupy')
dh = create_data_handling((4, 4), parallel=False, default_target=target)
dh.add_array('pdfs', values_per_cell=9, cpu=True, gpu=target != Target.CPU)
......@@ -45,7 +47,7 @@ def test_simple(target):
bh = LatticeBoltzmannBoundaryHandling(lb_func.method, dh, 'pdfs', target=target)
wall = NoSlip()
moving_wall = UBB((1, 0))
moving_wall = UBB((1, 0), density=1.0)
bh.set_boundary(wall, make_slice[0, :])
bh.set_boundary(wall, make_slice[-1, :])
bh.set_boundary(wall, make_slice[:, 0])
......@@ -435,3 +437,45 @@ def test_boundary_utility_functions():
assert stream == StreamInConstant(constant=1.0, name="stream")
assert not stream == StreamInConstant(constant=1.0, name="test")
assert not stream == noslip
@pytest.mark.parametrize("given_force_vector", [True, False])
@pytest.mark.parametrize("dtype", ["float32", "float64"])
def test_force_on_boundary(given_force_vector, dtype):
stencil = LBStencil(Stencil.D2Q9)
pdfs = ps.fields(f"pdfs_src({stencil.Q}): {dtype}[{stencil.D}D]", layout='fzyx')
method = create_lb_method(lbm_config=LBMConfig(stencil=stencil, method=Method.SRT, relaxation_rate=1.8))
noslip = NoSlip(name="noslip", calculate_force_on_boundary=True)
bouzidi = NoSlipLinearBouzidi(name="bouzidi", calculate_force_on_boundary=True)
qq_bounce_Back = QuadraticBounceBack(name="qqBB", relaxation_rate=1.8, calculate_force_on_boundary=True)
boundary_objects = [noslip, bouzidi, qq_bounce_Back]
for boundary in boundary_objects:
if given_force_vector:
force_vector_type = np.dtype([(f"F_{i}", dtype) for i in range(stencil.D)], align=True)
force_vector = ps.Field('forceVector', ps.FieldType.INDEXED, force_vector_type, layout=[0],
shape=(ps.TypedSymbol("forceVectorSize", "int32"), 1), strides=(1, 1))
else:
force_vector = None
index_struct_dtype = _numpy_data_type_for_boundary_object(boundary, stencil.D)
index_field = ps.Field('indexVector', ps.FieldType.INDEXED, index_struct_dtype, layout=[0],
shape=(ps.TypedSymbol("indexVectorSize", "int32"), 1), strides=(1, 1))
create_lattice_boltzmann_boundary_kernel(pdfs, index_field, method, boundary, force_vector=force_vector)
def _numpy_data_type_for_boundary_object(boundary_object, dim):
boundary_index_array_coordinate_names = ["x", "y", "z"]
direction_member_name = "dir"
default_index_array_dtype = np.int32
coordinate_names = boundary_index_array_coordinate_names[:dim]
return np.dtype(
[(name, default_index_array_dtype) for name in coordinate_names]
+ [(direction_member_name, default_index_array_dtype)]
+ [(i[0], i[1].numpy_dtype) for i in boundary_object.additional_data],
align=True,)
%% Cell type:code id: tags:
``` python
import pystencils as ps
from lbmpy.session import *
from lbmpy.boundaries.boundaries_in_kernel import update_rule_with_push_boundaries
from lbmpy.macroscopic_value_kernels import macroscopic_values_getter, macroscopic_values_setter
from collections import OrderedDict
from time import perf_counter
```
%% Cell type:markdown id: tags:
# Version 1: compile-in boundaries
%% Cell type:code id: tags:
``` python
domain_size = (32, 32, 32)
relaxation_rate = 1.8
time_steps = 100
lid_velocity = 0.05
stencil = LBStencil(Stencil.D3Q19)
```
%% Cell type:code id: tags:
``` python
dh = create_data_handling(domain_size, default_target=ps.Target.CPU)
pdfs = dh.add_array('pdfs', values_per_cell=stencil.Q)
u = dh.add_array('u', values_per_cell=stencil.D)
streaming_pattern = 'aa'
```
%% Cell type:code id: tags:
``` python
boundaries = OrderedDict((
((0, 1, 0), UBB([lid_velocity, 0, 0])),
((1, 0, 0), NoSlip()),
((-1, 0, 0), NoSlip()),
((0, -1, 0), NoSlip()),
((0, 0, 1), NoSlip()),
((0, 0, -1), NoSlip()),
))
lbm_opt = LBMOptimisation(symbolic_field=pdfs, cse_global=False, cse_pdfs=False)
lbm_config = LBMConfig(stencil=stencil, relaxation_rate=relaxation_rate, compressible=False)
cr_even = create_lb_collision_rule(lbm_config=lbm_config, lbm_optimisation=lbm_opt)
cr_odd = create_lb_collision_rule(lbm_config=lbm_config, lbm_optimisation=lbm_opt)
update_rule_aa_even = update_rule_with_push_boundaries(cr_even, pdfs, boundaries, streaming_pattern, Timestep.EVEN)
update_rule_aa_odd = update_rule_with_push_boundaries(cr_odd, pdfs, boundaries, streaming_pattern, Timestep.ODD)
getter_assignments = macroscopic_values_getter(update_rule_aa_even.method, velocity=u.center_vector,
pdfs=pdfs, density=None,
streaming_pattern=streaming_pattern,
previous_timestep=Timestep.EVEN)
config = ps.CreateKernelConfig(target=dh.default_target)
getter_kernel = ps.create_kernel(getter_assignments, config=config).compile()
even_kernel = ps.create_kernel(update_rule_aa_even, config=config).compile()
odd_kernel = ps.create_kernel(update_rule_aa_odd, config=config).compile()
```
%% Cell type:code id: tags:
``` python
def init():
dh.fill(pdfs.name, 0, ghost_layers=True)
def aa_time_loop(steps=100):
assert steps % 2 == 0, "Works only for an even number of time steps"
dh.all_to_gpu()
for i in range(steps // 2):
dh.run_kernel(odd_kernel)
dh.run_kernel(even_kernel)
dh.run_kernel(getter_kernel)
dh.all_to_cpu()
```
%% Cell type:code id: tags:
``` python
init()
aa_time_loop(time_steps)
vel_version1 = dh.gather_array(u.name, ghost_layers=False).copy()
plt.vector_field_magnitude(vel_version1[:, :, domain_size[2]//2, :])
plt.colorbar()
```
%% Output
<matplotlib.colorbar.Colorbar at 0x11b578710>
%% Cell type:markdown id: tags:
# Version 2: Normal boundary handling
%% Cell type:code id: tags:
``` python
ldc = create_lid_driven_cavity(domain_size, relaxation_rate=relaxation_rate, lid_velocity=lid_velocity)
ldc.run(time_steps)
vel_version2 = ldc.velocity[:, :, :, :]
plt.vector_field_magnitude(vel_version2[:, :, domain_size[2]//2, :])
plt.colorbar()
```
%% Output
<matplotlib.colorbar.Colorbar at 0x12e912a90>
%% Cell type:code id: tags:
``` python
np.testing.assert_almost_equal(vel_version1[1:-1, 1:-1, :], vel_version2[1:-1, 1:-1, :], decimal=2)
```