Skip to content

Commit

Permalink
Rewrite waLBerla MPI communication
Browse files Browse the repository at this point in the history
Split LB ghost communicator from LB boundaries communicator.
Always communicate bounce-back velocities from the ghost layer.
This fixes the regression introduced by
3fd1709.
  • Loading branch information
jngrad committed Feb 14, 2024
1 parent 639a556 commit d0c9154
Show file tree
Hide file tree
Showing 11 changed files with 306 additions and 32 deletions.
2 changes: 1 addition & 1 deletion samples/lb_circular_couette.py
Expand Up @@ -57,7 +57,7 @@
cyl_center = agrid * (grid_size // 2 + 0.5) * [1, 1, 0]
cylinder_in = espressomd.shapes.Cylinder(
center=cyl_center, axis=[0, 0, 1], length=3 * system.box_l[2],
radius=8.1 * agrid, direction=1)
radius=8.6 * agrid, direction=1)
cylinder_out = espressomd.shapes.Cylinder(
center=cyl_center, axis=[0, 0, 1], length=3 * system.box_l[2],
radius=14.5 * agrid, direction=-1)
Expand Down
3 changes: 1 addition & 2 deletions src/script_interface/walberla/LBFluid.cpp
Expand Up @@ -108,7 +108,6 @@ Variant LBFluid::do_call_method(std::string const &name,
}
if (name == "clear_boundaries") {
m_instance->clear_boundaries();
m_instance->ghost_communication();
::System::get_system().on_lb_boundary_conditions_change();
return {};
}
Expand Down Expand Up @@ -269,8 +268,8 @@ void LBFluid::load_checkpoint(std::string const &filename, int mode) {
};

auto const on_success = [&lb_obj]() {
lb_obj.reallocate_ubb_field();
lb_obj.ghost_communication();
lb_obj.reallocate_ubb_field();
};

load_checkpoint_common(*context(), "LB", filename, mode, read_metadata,
Expand Down
1 change: 1 addition & 0 deletions src/script_interface/walberla/LBFluidNode.cpp
Expand Up @@ -53,6 +53,7 @@ Variant LBFluidNode::do_call_method(std::string const &name,
m_lb_fluid->set_node_velocity_at_boundary(m_index, u);
}
m_lb_fluid->ghost_communication();
m_lb_fluid->reallocate_ubb_field();
return {};
}
if (name == "get_velocity_at_boundary") {
Expand Down
6 changes: 4 additions & 2 deletions src/script_interface/walberla/LBFluidSlice.cpp
Expand Up @@ -99,8 +99,10 @@ Variant LBFluidSlice::do_call_method(std::string const &name,
1. / m_conv_velocity);
}
if (name == "set_velocity_at_boundary") {
return call(&LatticeModel::set_slice_velocity_at_boundary, {1},
m_conv_velocity);
auto const retval = call(&LatticeModel::set_slice_velocity_at_boundary, {1},
m_conv_velocity);
m_lb_fluid->reallocate_ubb_field();
return retval;
}
if (name == "get_pressure_tensor") {
return call(&LatticeModel::get_slice_pressure_tensor, {3, 3},
Expand Down
147 changes: 147 additions & 0 deletions src/walberla_bridge/src/BoundaryPackInfo.hpp
@@ -0,0 +1,147 @@
/*
* Copyright (C) 2024 The ESPResSo project
*
* This file is part of ESPResSo.
*
* ESPResSo is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* ESPResSo is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/

#pragma once

#include <core/debug/Debug.h>
#include <core/mpi/RecvBuffer.h>
#include <core/mpi/SendBuffer.h>
#include <domain_decomposition/IBlock.h>
#include <field/communication/PackInfo.h>
#include <stencil/Directions.h>

#include <memory>
#include <tuple>
#include <utility>

namespace walberla {
namespace field {
namespace communication {

template <typename GhostLayerField_T, typename Boundary_T>
class BoundaryPackInfo : public PackInfo<GhostLayerField_T> {
protected:
using PackInfo<GhostLayerField_T>::bdId_;

public:
using PackInfo<GhostLayerField_T>::PackInfo;
using PackInfo<GhostLayerField_T>::numberOfGhostLayersToCommunicate;

~BoundaryPackInfo() override = default;

void setup_boundary_handle(std::shared_ptr<LatticeWalberla> lattice,
std::shared_ptr<Boundary_T> boundary) {
m_lattice = std::move(lattice);
m_boundary = std::move(boundary);
}

bool constantDataExchange() const override { return false; }
bool threadsafeReceiving() const override { return true; }

void communicateLocal(IBlock const *sender, IBlock *receiver,
stencil::Direction dir) override {
mpi::SendBuffer sBuffer;
packDataImpl(sender, dir, sBuffer);
mpi::RecvBuffer rBuffer(sBuffer);
unpackData(receiver, stencil::inverseDir[dir], rBuffer);
}

void unpackData(IBlock *receiver, stencil::Direction dir,
mpi::RecvBuffer &buffer) override {

auto *flag_field = receiver->getData<GhostLayerField_T>(bdId_);
WALBERLA_ASSERT_NOT_NULLPTR(flag_field);
WALBERLA_ASSERT_NOT_NULLPTR(m_boundary);
WALBERLA_ASSERT_NOT_NULLPTR(m_lattice);

auto const boundary_flag = flag_field->getFlag(Boundary_flag);
auto const gl = numberOfGhostLayersToCommunicate(flag_field);
auto const begin = [gl, dir](auto const *flag_field) {
return flag_field->beginGhostLayerOnly(gl, dir);
};

#ifndef NDEBUG
uint_t xSize, ySize, zSize, bSize;
buffer >> xSize >> ySize >> zSize >> bSize;
uint_t buf_size{0u};
for (auto it = begin(flag_field); it != flag_field->end(); ++it) {
if (isFlagSet(it, boundary_flag)) {
++buf_size;
}
}
WALBERLA_ASSERT_EQUAL(xSize, flag_field->xSize());
WALBERLA_ASSERT_EQUAL(ySize, flag_field->ySize());
WALBERLA_ASSERT_EQUAL(zSize, flag_field->zSize());
WALBERLA_ASSERT_EQUAL(bSize, buf_size);
#endif

auto const offset = std::get<0>(m_lattice->get_local_grid_range());
typename Boundary_T::value_type value;
for (auto it = begin(flag_field); it != flag_field->end(); ++it) {
if (isFlagSet(it, boundary_flag)) {
auto const node = offset + Utils::Vector3i{{it.x(), it.y(), it.z()}};
buffer >> value;
m_boundary->unpack_node(node, value);
}
}
}

protected:
void packDataImpl(IBlock const *sender, stencil::Direction dir,
mpi::SendBuffer &buffer) const override {

auto const *flag_field = sender->getData<GhostLayerField_T>(bdId_);
WALBERLA_ASSERT_NOT_NULLPTR(flag_field);
WALBERLA_ASSERT_NOT_NULLPTR(m_boundary);
WALBERLA_ASSERT_NOT_NULLPTR(m_lattice);

auto const boundary_flag = flag_field->getFlag(Boundary_flag);
auto const gl = numberOfGhostLayersToCommunicate(flag_field);
auto const begin = [gl, dir](auto const *flag_field) {
return flag_field->beginSliceBeforeGhostLayer(dir, gl);
};

#ifndef NDEBUG
uint_t buf_size{0u};
for (auto it = begin(flag_field); it != flag_field->end(); ++it) {
if (isFlagSet(it, boundary_flag)) {
++buf_size;
}
}
buffer << flag_field->xSize() << flag_field->ySize() << flag_field->zSize()
<< buf_size;
#endif

auto const offset = std::get<0>(m_lattice->get_local_grid_range());
for (auto it = begin(flag_field); it != flag_field->end(); ++it) {
if (isFlagSet(it, boundary_flag)) {
auto const node = offset + Utils::Vector3i{{it.x(), it.y(), it.z()}};
buffer << m_boundary->get_node_value_at_boundary(node);
}
}
}

private:
std::shared_ptr<LatticeWalberla> m_lattice;
std::shared_ptr<Boundary_T> m_boundary;
};

} // namespace communication
} // namespace field
} // namespace walberla
63 changes: 42 additions & 21 deletions src/walberla_bridge/src/lattice_boltzmann/LBWalberlaImpl.hpp
Expand Up @@ -45,6 +45,7 @@
#include <stencil/D3Q27.h>

#include "../BoundaryHandling.hpp"
#include "../BoundaryPackInfo.hpp"
#include "InterpolateAndShiftAtBoundary.hpp"
#include "ResetForce.hpp"
#include "lb_kernels.hpp"
Expand Down Expand Up @@ -235,8 +236,9 @@ class LBWalberlaImpl : public LBWalberlaBase {
typename FieldTrait<FloatType, Architecture>::template PackInfo<Field>;

// communicators
std::shared_ptr<FullCommunicator> m_full_communication;
std::shared_ptr<PDFStreamingCommunicator> m_pdf_streaming_communication;
std::shared_ptr<FullCommunicator> m_boundary_communicator;
std::shared_ptr<FullCommunicator> m_pdf_full_communicator;
std::shared_ptr<PDFStreamingCommunicator> m_pdf_streaming_communicator;

// ResetForce sweep + external force handling
std::shared_ptr<ResetForce<PdfField, VectorField>> m_reset_force;
Expand Down Expand Up @@ -350,27 +352,33 @@ class LBWalberlaImpl : public LBWalberlaBase {
reset_boundary_handling();

// Set up the communication and register fields
m_pdf_streaming_communication =
m_pdf_streaming_communicator =
std::make_shared<PDFStreamingCommunicator>(blocks);
m_pdf_streaming_communication->addPackInfo(
m_pdf_streaming_communicator->addPackInfo(
std::make_shared<PackInfo<PdfField>>(m_pdf_field_id, n_ghost_layers));
m_pdf_streaming_communication->addPackInfo(
m_pdf_streaming_communicator->addPackInfo(
std::make_shared<PackInfo<VectorField>>(m_last_applied_force_field_id,
n_ghost_layers));
m_pdf_streaming_communication->addPackInfo(
std::make_shared<field::communication::PackInfo<FlagField>>(
m_flag_field_id, n_ghost_layers));

m_full_communication = std::make_shared<FullCommunicator>(blocks);
m_full_communication->addPackInfo(
m_pdf_full_communicator = std::make_shared<FullCommunicator>(blocks);
m_pdf_full_communicator->addPackInfo(
std::make_shared<PackInfo<PdfField>>(m_pdf_field_id, n_ghost_layers));
m_full_communication->addPackInfo(std::make_shared<PackInfo<VectorField>>(
m_last_applied_force_field_id, n_ghost_layers));
m_full_communication->addPackInfo(std::make_shared<PackInfo<VectorField>>(
m_velocity_field_id, n_ghost_layers));
m_full_communication->addPackInfo(
m_pdf_full_communicator->addPackInfo(
std::make_shared<PackInfo<VectorField>>(m_last_applied_force_field_id,
n_ghost_layers));
m_pdf_full_communicator->addPackInfo(
std::make_shared<PackInfo<VectorField>>(m_velocity_field_id,
n_ghost_layers));

m_boundary_communicator = std::make_shared<FullCommunicator>(blocks);
m_boundary_communicator->addPackInfo(
std::make_shared<field::communication::PackInfo<FlagField>>(
m_flag_field_id, n_ghost_layers));
auto boundary_packinfo = std::make_shared<
field::communication::BoundaryPackInfo<FlagField, BoundaryModel>>(
m_flag_field_id, n_ghost_layers);
boundary_packinfo->setup_boundary_handle(m_lattice, m_boundary);
m_boundary_communicator->addPackInfo(boundary_packinfo);

// Instantiate the sweep responsible for force double buffering and
// external forces
Expand Down Expand Up @@ -439,13 +447,13 @@ class LBWalberlaImpl : public LBWalberlaBase {
integrate_reset_force(blocks);
// LB collide
integrate_collide(blocks);
m_pdf_streaming_communication->communicate();
m_pdf_streaming_communicator->communicate();
// Handle boundaries
integrate_boundaries(blocks);
// LB stream
integrate_stream(blocks);
// Refresh ghost layers
m_full_communication->communicate();
ghost_communication_pdfs();
}

void integrate_pull_scheme() {
Expand All @@ -458,7 +466,7 @@ class LBWalberlaImpl : public LBWalberlaBase {
// LB collide
integrate_collide(blocks);
// Refresh ghost layers
ghost_communication();
ghost_communication_pdfs();
}

protected:
Expand All @@ -474,7 +482,6 @@ class LBWalberlaImpl : public LBWalberlaBase {

public:
void integrate() override {
reallocate_ubb_field();
if (has_lees_edwards_bc()) {
integrate_pull_scheme();
} else {
Expand All @@ -485,7 +492,16 @@ class LBWalberlaImpl : public LBWalberlaBase {
}

void ghost_communication() override {
m_full_communication->communicate();
ghost_communication_boundary();
ghost_communication_pdfs();
}

void ghost_communication_boundary() {
m_boundary_communicator->communicate();
}

void ghost_communication_pdfs() {
m_pdf_full_communicator->communicate();
if (has_lees_edwards_bc()) {
auto const &blocks = get_lattice().get_blocks();
apply_lees_edwards_pdf_interpolation(blocks);
Expand Down Expand Up @@ -1097,14 +1113,19 @@ class LBWalberlaImpl : public LBWalberlaBase {

void reallocate_ubb_field() override { m_boundary->boundary_update(); }

void clear_boundaries() override { reset_boundary_handling(); }
void clear_boundaries() override {
reset_boundary_handling();
ghost_communication();
}

void
update_boundary_from_shape(std::vector<int> const &raster_flat,
std::vector<double> const &data_flat) override {
auto const grid_size = get_lattice().get_grid_dimensions();
auto const data = fill_3D_vector_array(data_flat, grid_size);
set_boundary_from_grid(*m_boundary, get_lattice(), raster_flat, data);
ghost_communication();
reallocate_ubb_field();
}

// Pressure tensor
Expand Down
1 change: 1 addition & 0 deletions src/walberla_bridge/tests/LBWalberlaImpl_unit_tests.cpp
Expand Up @@ -202,6 +202,7 @@ BOOST_DATA_TEST_CASE(update_boundary_from_shape, bdata::make(all_lbs()),
std::vector<double> vel_flat(vel_3d.data(),
vel_3d.data() + vel_3d.num_elements());
lb->update_boundary_from_shape(raster_flat, vel_flat);
lb->ghost_communication();
}

for (auto const &node : nodes) {
Expand Down
1 change: 1 addition & 0 deletions testsuite/python/CMakeLists.txt
Expand Up @@ -334,6 +334,7 @@ python_test(FILE thole.py MAX_NUM_PROC 4)
python_test(FILE lb_slice.py MAX_NUM_PROC 2)
python_test(FILE lb_boundary_velocity.py MAX_NUM_PROC 1)
# python_test(FILE lb_boundary_volume_force.py MAX_NUM_PROC 2) # TODO
python_test(FILE lb_boundary_ghost_layer.py MAX_NUM_PROC 2)
python_test(FILE lb_circular_couette.py MAX_NUM_PROC 2 GPU_SLOTS 1)
python_test(FILE lb_poiseuille.py MAX_NUM_PROC 4 GPU_SLOTS 1)
python_test(FILE lb_poiseuille_cylinder.py MAX_NUM_PROC 2 GPU_SLOTS 1)
Expand Down

0 comments on commit d0c9154

Please sign in to comment.