ESPResSo
Extensible Simulation Package for Research on Soft Matter Systems
Loading...
Searching...
No Matches
HybridDecomposition.hpp
Go to the documentation of this file.
1/*
2 * Copyright (C) 2010-2022 The ESPResSo project
3 * Copyright (C) 2002,2003,2004,2005,2006,2007,2008,2009,2010
4 * Max-Planck-Institute for Polymer Research, Theory Group
5 *
6 * This file is part of ESPResSo.
7 *
8 * ESPResSo is free software: you can redistribute it and/or modify
9 * it under the terms of the GNU General Public License as published by
10 * the Free Software Foundation, either version 3 of the License, or
11 * (at your option) any later version.
12 *
13 * ESPResSo is distributed in the hope that it will be useful,
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 * GNU General Public License for more details.
17 *
18 * You should have received a copy of the GNU General Public License
19 * along with this program. If not, see <http://www.gnu.org/licenses/>.
20 */
21
22#pragma once
23
27
28#include "cell_system/Cell.hpp"
29
30#include "BoxGeometry.hpp"
31#include "LocalBox.hpp"
32#include "Particle.hpp"
33#include "ghosts.hpp"
34
35#include <utils/Vector.hpp>
36
37#include <boost/mpi/communicator.hpp>
38
39#include <cstddef>
40#include <functional>
41#include <optional>
42#include <set>
43#include <span>
44#include <utility>
45#include <vector>
46
47/**
48 * @brief Hybrid decomposition cell system.
49 *
50 * Store particles with short-range interactions
51 * in a @ref RegularDecomposition cell system and
52 * particles with long-range interactions
53 * in a @ref AtomDecomposition (N-square) cell system.
54 * All regular cells are coupled to the N-square cells.
55 */
57 boost::mpi::communicator m_comm;
58 BoxGeometry const &m_box;
59 double m_cutoff_regular;
60 std::vector<Cell *> m_local_cells;
61 std::vector<Cell *> m_ghost_cells;
62
63 GhostCommunicator m_exchange_ghosts_comm;
64 GhostCommunicator m_collect_ghost_force_comm;
65
66 /** RegularDecomposition to hold the small particles */
67 RegularDecomposition m_regular_decomposition;
68 /** N-Square Decomposition to hold large particles */
69 AtomDecomposition m_n_square;
70 /** Set containing the types that should be handled using n_square */
71 std::set<int> const m_n_square_types;
72
73 std::function<bool()> m_get_global_ghost_flags;
74
75 bool is_n_square_type(int type_id) const {
76 return (m_n_square_types.find(type_id) != m_n_square_types.end());
77 }
78
79public:
80 HybridDecomposition(boost::mpi::communicator comm, double cutoff_regular,
81 double skin, std::function<bool()> get_ghost_flags,
82 BoxGeometry const &box_geo, LocalBox const &local_box,
83 std::set<int> n_square_types);
84
85 auto get_cell_grid() const { return m_regular_decomposition.cell_grid; }
86
87 auto get_cell_size() const { return m_regular_decomposition.cell_size; }
88
89 auto get_n_square_types() const { return m_n_square_types; }
90
91 void resort(bool global, std::vector<ParticleChange> &diff) override;
92
93 auto get_cutoff_regular() const { return m_cutoff_regular; }
94
95 GhostCommunicator const &exchange_ghosts_comm() const override {
96 return m_exchange_ghosts_comm;
97 }
98
100 return m_collect_ghost_force_comm;
101 }
102
103 std::span<Cell *const> local_cells() const override { return m_local_cells; }
104 std::span<Cell *const> ghost_cells() const override { return m_ghost_cells; }
105
106 Cell *particle_to_cell(Particle const &p) override {
107 if (is_n_square_type(p.type())) {
108 return m_n_square.particle_to_cell(p);
109 }
110 return m_regular_decomposition.particle_to_cell(p);
111 }
112
113 Cell const *particle_to_cell(Particle const &p) const override {
114 if (is_n_square_type(p.type())) {
115 return m_n_square.particle_to_cell(p);
116 }
117 return m_regular_decomposition.particle_to_cell(p);
118 }
119
120 Utils::Vector3d max_cutoff() const override {
121 return m_n_square.max_cutoff();
122 }
123
124 Utils::Vector3d max_range() const override { return m_n_square.max_range(); }
125
126 std::optional<BoxGeometry> minimum_image_distance() const override {
127 return m_box;
128 }
129
130 BoxGeometry const &box() const override { return m_box; }
131
132 /** @brief Count particles in child regular decompositions. */
133 std::size_t count_particles_in_regular() const {
134 return count_particles(m_regular_decomposition.get_local_cells());
135 }
136
137 /** @brief Count particles in child N-square decompositions. */
138 std::size_t count_particles_in_n_square() const {
139 return count_particles(m_n_square.get_local_cells());
140 }
141
142private:
143 std::size_t count_particles(std::vector<Cell *> const &local_cells) const;
144};
Vector implementation and trait types for boost qvm interoperability.
Atom decomposition cell system.
Utils::Vector3d max_range() const override
Utils::Vector3d max_cutoff() const override
auto const & get_local_cells() const
Cell * particle_to_cell(Particle const &p) override
Determine which cell a particle id belongs to.
Definition Cell.hpp:97
Hybrid decomposition cell system.
std::span< Cell *const > local_cells() const override
GhostCommunicator const & exchange_ghosts_comm() const override
std::size_t count_particles_in_regular() const
Count particles in child regular decompositions.
Utils::Vector3d max_range() const override
std::size_t count_particles_in_n_square() const
Count particles in child N-square decompositions.
Cell * particle_to_cell(Particle const &p) override
Cell const * particle_to_cell(Particle const &p) const override
void resort(bool global, std::vector< ParticleChange > &diff) override
BoxGeometry const & box() const override
Utils::Vector3d max_cutoff() const override
GhostCommunicator const & collect_ghost_force_comm() const override
std::optional< BoxGeometry > minimum_image_distance() const override
std::span< Cell *const > ghost_cells() const override
A distributed particle decomposition.
Ghost particles and particle exchange.
Properties for a ghost communication.
Definition ghosts.hpp:157
Struct holding all information for one particle.
Definition Particle.hpp:395
auto const & type() const
Definition Particle.hpp:418
Regular decomposition cell system.
Cell * particle_to_cell(Particle const &p) override
Utils::Vector3d cell_size
Cell size.
auto const & get_local_cells() const
Utils::Vector3i cell_grid
Grid dimensions per node.