Skip to content
Permalink
Browse files

Fix (for gpu) +test lb velocity interpolation close to boundary

  • Loading branch information
RudolfWeeber committed Mar 25, 2020
1 parent 41bacc8 commit d7ad89006e2d6451f7602a676724749c11fd3854
Showing with 22 additions and 8 deletions.
  1. +3 −1 src/core/grid_based_algorithms/lbgpu_cuda.cu
  2. +19 −7 testsuite/python/lb_interpolation.py
@@ -1372,8 +1372,10 @@ __device__ __inline__ float3 node_velocity(float rho_eq, LB_nodes_gpu n_a,
auto const boundary_index = n_a.boundary[index];

if (boundary_index) {
auto const inv_lattice_speed = para->tau / para->agrid;
auto const &u = n_a.boundary_velocity[index];
return make_float3(u[0], u[1], u[2]);
return make_float3(inv_lattice_speed * u[0], inv_lattice_speed * u[1],
inv_lattice_speed * u[2]);
}

auto const rho = rho_eq + calc_mode_x_from_n(n_a, index, 0);
@@ -23,11 +23,11 @@
import espressomd.shapes
import espressomd.lb

AGRID = 1.5
AGRID = 1.5
VISC = 4.2
DENS = 1.3
FRIC = 1.4
TAU = 0.2
TAU = 0.3
BOX_L = 18.0
TIME_STEP = TAU
LB_PARAMETERS = {
@@ -57,9 +57,9 @@ class LBInterpolation:
def set_boundaries(self, velocity):
"""Place boundaries *not* exactly on a LB node."""
wall_shape1 = espressomd.shapes.Wall(
normal=[1, 0, 0], dist=0.6 * AGRID)
normal=[1, 0, 0], dist=AGRID)
wall_shape2 = espressomd.shapes.Wall(
normal=[-1, 0, 0], dist=-(BOX_L - 0.6 * AGRID))
normal=[-1, 0, 0], dist=-(BOX_L - AGRID))
self.system.lbboundaries.add(
espressomd.lbboundaries.LBBoundary(shape=wall_shape1))
self.system.lbboundaries.add(
@@ -72,9 +72,21 @@ def test_interpolated_velocity(self):
"""
self.set_boundaries([0.0, 0.0, V_BOUNDARY])
self.system.integrator.run(250)
# Shear plane for boundary 1
# for pos in itertools.product((AGRID,), np.arange(0.5 * AGRID, BOX_L, AGRID), np.arange(0.5 * AGRID, BOX_L, AGRID)):
# np.testing.assert_almost_equal(self.lbf.get_interpolated_velocity(pos)[2], 0.0)
# Check interpolated vel at upper boundary. The node position is at
# box_l[0]-agrid/2.
np.testing.assert_allclose(self.lbf.get_interpolated_velocity(
[self.system.box_l[0] - AGRID / 2, 0, 0]), np.array([0, 0, V_BOUNDARY]))

# Cehck interpolated velocity involving boundary and neighboring node
# The boundary node index is lbf.shape[0]-1, so -2 refers to the
# node in front of the boundary.
node_next_to_boundary = self.lbf[
self.lbf.shape[0] - 2, 0, 0]
# The midpoint between the boundary and
# that node is box_l - agrid.
np.testing.assert_allclose(self.lbf.get_interpolated_velocity([self.system.box_l[0] - AGRID, 0, 0]),
0.5 * (np.array([0, 0, V_BOUNDARY]) + node_next_to_boundary.velocity))

# Bulk
for pos in itertools.product(
np.arange(1.5 * AGRID, BOX_L - 1.5 * AGRID, 0.5 * AGRID),

0 comments on commit d7ad890

Please sign in to comment.
You can’t perform that action at this time.