ESPResSo
Extensible Simulation Package for Research on Soft Matter Systems
Loading...
Searching...
No Matches
LatticeSlice.impl.hpp
Go to the documentation of this file.
1/*
2 * Copyright (C) 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>
24
25#include <boost/mpi/collectives/gather.hpp>
26#include <boost/mpi/collectives/scatter.hpp>
27#include <boost/mpi/communicator.hpp>
28#include <boost/multi_array.hpp>
29#include <boost/serialization/vector.hpp>
30
31#include <algorithm>
32#include <functional>
33#include <optional>
34#include <stdexcept>
35#include <tuple>
36#include <type_traits>
37#include <vector>
38
40
41namespace detail {
42
43// boundary types are always std::optional types
44template <class> struct is_optional_type : public std::false_type {};
45template <class T>
46struct is_optional_type<std::optional<T>> : public std::true_type {};
47
48template <class ArrayView, typename T>
49void unflatten_grid(ArrayView &view, std::vector<T> const &values) {
50 using array_type = boost::multi_array<T, 4>;
51 auto it = values.begin();
52 auto const dim_i = static_cast<typename array_type::index>(view.shape()[0]);
53 auto const dim_j = static_cast<typename array_type::index>(view.shape()[1]);
54 auto const dim_k = static_cast<typename array_type::index>(view.shape()[2]);
55 auto const dim_t = static_cast<typename array_type::index>(view.shape()[3]);
56 for (typename array_type::index i = 0; i != dim_i; ++i) {
57 for (typename array_type::index j = 0; j != dim_j; ++j) {
58 for (typename array_type::index k = 0; k != dim_k; ++k) {
59 for (typename array_type::index t = 0; t != dim_t; ++t) {
60 view[i][j][k][t] = *it;
61 ++it;
62 }
63 }
64 }
65 }
66}
67
68template <class FieldSerializer, class ArrayView, typename T>
69void flatten_grid(ArrayView const &view, std::vector<T> &out,
70 double units_conversion) {
71 using array_type = boost::multi_array<T, 4>;
72 out.reserve(view.num_elements());
73 auto const dim_i = static_cast<typename array_type::index>(view.shape()[0]);
74 auto const dim_j = static_cast<typename array_type::index>(view.shape()[1]);
75 auto const dim_k = static_cast<typename array_type::index>(view.shape()[2]);
76 auto const dim_t = static_cast<typename array_type::index>(view.shape()[3]);
77 for (typename array_type::index i = 0; i != dim_i; ++i) {
78 for (typename array_type::index j = 0; j != dim_j; ++j) {
79 for (typename array_type::index k = 0; k != dim_k; ++k) {
80 for (typename array_type::index t = 0; t != dim_t; ++t) {
81 if constexpr (std::is_floating_point_v<T>) {
82 out.emplace_back(view[i][j][k][t] * units_conversion);
83 } else if constexpr (is_optional_type<T>{}) {
84 if (view[i][j][k][t]) {
85 out.emplace_back(*(view[i][j][k][t]) * units_conversion);
86 } else {
87 out.emplace_back(std::nullopt);
88 }
89 } else {
90 out.emplace_back(view[i][j][k][t]);
91 }
92 }
93 }
94 }
95 }
96}
97
98} // namespace detail
99
100inline auto gather_slices_topology(boost::mpi::communicator const &comm,
103 std::vector<Utils::Vector3i> nodes_lower_corners;
104 std::vector<Utils::Vector3i> nodes_upper_corners;
105 boost::mpi::gather(comm, local_lower_corner, nodes_lower_corners, 0);
106 boost::mpi::gather(comm, local_upper_corner, nodes_upper_corners, 0);
107 return std::make_tuple(nodes_lower_corners, nodes_upper_corners);
108}
109
110template <class FieldSerializer>
111template <class LatticeModel, typename T>
113 VariantMap const &params, std::vector<int> const &data_dims,
115 std::vector<T> (LatticeModel::*getter)(Utils::Vector3i const &,
116 Utils::Vector3i const &) const,
118 auto const &comm = context()->get_comm();
120 local_upper_corner] = get_slices_bounding_boxes();
123 auto const data_size = std::accumulate(data_dims.cbegin(), data_dims.cend(),
124 1, std::multiplies<>());
125 auto const local_values =
127 std::vector<std::vector<T>> nodes_values;
128 boost::mpi::gather(comm, local_values, nodes_values, 0);
129 if (comm.rank() == 0) {
130 auto const dims = slice_upper_corner - slice_lower_corner;
131 using index_range = boost::multi_array_types::index_range;
132 using array_type = boost::multi_array<T, 4>;
133 array_type array(boost::extents[dims[0]][dims[1]][dims[2]][data_size]);
134 // populate the 3D array with data from each node
135 for (std::size_t rank = 0; rank < nodes_values.size(); ++rank) {
136 if (nodes_values[rank].empty()) {
137 continue;
138 }
139 auto const range_lower_corner =
141 auto const range_upper_corner =
143 auto const local_range = [&](int j) {
145 };
146 typename array_type::template array_view<4>::type view =
147 array[boost::indices[local_range(0)][local_range(1)][local_range(2)]
148 [index_range()]];
149 detail::unflatten_grid(view, nodes_values[rank]);
150 }
151 // create the output flat array
152 std::vector<T> out;
153 detail::flatten_grid<FieldSerializer>(array, out, units_conversion);
154 std::vector<int> shape = {m_shape.begin(), m_shape.end()};
155 if (data_dims.size() != 1ul or data_dims[0] != 1) {
156 shape.insert(shape.end(), data_dims.begin(), data_dims.end());
157 }
158 auto const variant = FieldSerializer::serialize(out);
159 return {std::vector<Variant>{{variant, Variant{shape}}}};
160 }
161 return {};
162}
163
164template <class FieldSerializer>
165template <class LatticeModel, typename T>
167 VariantMap const &params, std::vector<int> const &data_dims,
168 LatticeModel &lattice_model,
169 void (LatticeModel::*setter)(Utils::Vector3i const &,
170 Utils::Vector3i const &,
171 std::vector<T> const &),
172 double units_conversion) {
173 auto const &comm = context()->get_comm();
175 local_upper_corner] = get_slices_bounding_boxes();
178 auto const data_size = std::accumulate(data_dims.cbegin(), data_dims.cend(),
179 1, std::multiplies<>());
180 auto const sentinel = get_sentinel_index(lattice_model.get_lattice());
181 std::vector<std::vector<T>> nodes_values(comm.size());
182 if (comm.rank() == 0) {
183 auto const values =
184 FieldSerializer::template deserialize<T>(params.at("values"));
185 auto const dims = slice_upper_corner - slice_lower_corner;
186 using index_range = boost::multi_array_types::index_range;
187 using array_type = boost::multi_array<T, 4>;
188 array_type array(boost::extents[dims[0]][dims[1]][dims[2]][data_size]);
189 // populate the 3D array from the input flat array
190 detail::unflatten_grid(array, values);
191 // partition the 3D array into individual flat arrays for each MPI rank
192 for (std::size_t rank = 0; rank < nodes_lower_corners.size(); ++rank) {
196 continue;
197 }
198 auto const local_range = [&](int j) {
200 };
201 typename array_type::template array_view<4>::type view =
202 array[boost::indices[local_range(0)][local_range(1)][local_range(2)]
203 [index_range()]];
204 detail::flatten_grid<FieldSerializer>(view, nodes_values[rank],
206 }
207 }
208 std::vector<T> local_values;
209 boost::mpi::scatter(comm, nodes_values, local_values, 0);
211 lattice_model.ghost_communication();
212}
213
214} // namespace ScriptInterface::walberla
Vector implementation and trait types for boost qvm interoperability.
virtual boost::mpi::communicator const & get_comm() const =0
Context * context() const
Responsible context.
Variant gather_3d(VariantMap const &params, std::vector< int > const &data_dims, LatticeModel const &lattice_model, std::vector< T >(LatticeModel::*getter)(Utils::Vector3i const &, Utils::Vector3i const &) const, double units_conversion=1.) const
static DEVICE_QUALIFIER constexpr Vector< T, N > broadcast(typename Base::value_type const &value)
Create a vector that has all entries set to the same value.
Definition Vector.hpp:110
auto gather_slices_topology(boost::mpi::communicator const &comm, Utils::Vector3i const &local_lower_corner, Utils::Vector3i const &local_upper_corner)
T get_value(Variant const &v)
Extract value of specific type T from a Variant.
std::unordered_map< std::string, Variant > VariantMap
Definition Variant.hpp:69
boost::make_recursive_variant< None, bool, int, std::size_t, double, std::string, ObjectRef, Utils::Vector3b, Utils::Vector3i, Utils::Vector2d, Utils::Vector3d, Utils::Vector4d, std::vector< int >, std::vector< double >, std::vector< boost::recursive_variant_ >, std::unordered_map< int, boost::recursive_variant_ >, std::unordered_map< std::string, boost::recursive_variant_ > >::type Variant
Possible types for parameters.
Definition Variant.hpp:67
static SteepestDescentParameters params
Currently active steepest descent instance.