ESPResSo
Extensible Simulation Package for Research on Soft Matter Systems
Loading...
Searching...
No Matches
walberla_bridge/include/walberla_bridge/LatticeWalberla.hpp
Go to the documentation of this file.
1/*
2 * Copyright (C) 2021-2023 The ESPResSo project
3 *
4 * This file is part of ESPResSo.
5 *
6 * ESPResSo is free software: you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License as published by
8 * the Free Software Foundation, either version 3 of the License, or
9 * (at your option) any later version.
10 *
11 * ESPResSo is distributed in the hope that it will be useful,
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 * GNU General Public License for more details.
15 *
16 * You should have received a copy of the GNU General Public License
17 * along with this program. If not, see <http://www.gnu.org/licenses/>.
18 */
19
20#pragma once
21
22#include <utils/Vector.hpp>
23
24#include <cassert>
25#include <cmath>
26#include <initializer_list>
27#include <memory>
28#include <utility>
29
30// forward declarations
32class StructuredBlockForest;
33} // namespace walberla::blockforest
35class IBlock;
36} // namespace walberla::domain_decomposition
37
38/** Class that runs and controls the BlockForest in waLBerla. */
40public:
41 using Lattice_T = walberla::blockforest::StructuredBlockForest;
42
43private:
44 Utils::Vector3i m_grid_dimensions;
45 unsigned int m_n_ghost_layers;
46
47 /** Block forest */
48 std::shared_ptr<Lattice_T> m_blocks;
49 using IBlock = walberla::domain_decomposition::IBlock;
50 std::vector<IBlock *> m_cached_blocks;
51
52public:
53 LatticeWalberla(Utils::Vector3i const &grid_dimensions,
54 Utils::Vector3i const &node_grid,
55 unsigned int n_ghost_layers);
56
57 // Grid, domain, halo
58 [[nodiscard]] auto get_ghost_layers() const { return m_n_ghost_layers; }
59 [[nodiscard]] auto get_grid_dimensions() const { return m_grid_dimensions; }
60 [[nodiscard]] auto get_blocks() const { return m_blocks; }
61 [[nodiscard]] auto const &get_cached_blocks() const {
62 return m_cached_blocks;
63 }
64 [[nodiscard]] std::pair<Utils::Vector3d, Utils::Vector3d>
65 get_local_domain() const;
66 [[nodiscard]] auto get_local_grid_range() const {
67 auto const conversion = [](Utils::Vector3d const &pos) -> Utils::Vector3i {
68 auto const dim =
69 Utils::Vector3i{{static_cast<int>(pos[0]), static_cast<int>(pos[1]),
70 static_cast<int>(pos[2])}};
71#ifndef NDEBUG
72 for (auto const i : {0u, 1u, 2u}) {
73 assert(std::abs(static_cast<double>(dim[i]) - pos[i]) < 1e-10);
74 }
75#endif
76 return dim;
77 };
78 auto const [lower_corner, upper_corner] = get_local_domain();
79 return std::make_pair(conversion(lower_corner), conversion(upper_corner));
80 }
81
82 [[nodiscard]] bool node_in_local_domain(Utils::Vector3i const &node) const;
83 [[nodiscard]] bool node_in_local_halo(Utils::Vector3i const &node) const;
84 [[nodiscard]] bool pos_in_local_domain(Utils::Vector3d const &pos) const;
85 [[nodiscard]] bool pos_in_local_halo(Utils::Vector3d const &pos) const;
86 [[nodiscard]] static Utils::Vector3i
87 calc_grid_dimensions(Utils::Vector3d const &box_size, double agrid);
88};
Vector implementation and trait types for boost qvm interoperability.
Class that runs and controls the BlockForest in waLBerla.
bool pos_in_local_halo(Utils::Vector3d const &pos) const
static Utils::Vector3i calc_grid_dimensions(Utils::Vector3d const &box_size, double agrid)
bool node_in_local_halo(Utils::Vector3i const &node) const
std::pair< Utils::Vector3d, Utils::Vector3d > get_local_domain() const
walberla::blockforest::StructuredBlockForest Lattice_T
bool node_in_local_domain(Utils::Vector3i const &node) const
bool pos_in_local_domain(Utils::Vector3d const &pos) const