FCountKernel.hpp 10.6 KB
Newer Older
1 2 3
#ifndef FCOUNT_KERNEL_HPP
#define FCOUNT_KERNEL_HPP

4 5 6
#include <iostream>

#include <algorithm>
7 8
#include <cstddef>
#include <tuple>
9
#include <array>
10
#include <unordered_map>
11
#include <vector>
12

13
#include <mutex>
14 15
#include <thread>
#include <chrono>
16

17
#include "Containers/FTreeCoordinate.hpp"
18
#include "Components/FBasicParticle.hpp"
19 20

#include "Kernels/Generic/FGenericData.hpp"
21

22 23
#include "Utils/FOstreamTuple.hpp"

24 25 26 27 28 29 30
/**
 * \brief Particle with a count attribute for use with FCountKernel
 *
 * FBasicParticle additional tempalte parameters:
 *   - particle index (std::size_t)
 *   - particle interaction count (std::size_t)
 */
31
template<typename FReal, std::size_t Dim>
32 33
struct TestCountParticle : public FBasicParticle<FReal, Dim, std::size_t, std::size_t> {
    using FBase = FBasicParticle<FReal, Dim, std::size_t, std::size_t>;
34 35 36
    using FBase::FBase;
    /// Data index names
    enum {POS1, POS2, POS3};
37
    enum {IDX = Dim, COUNT};
38 39
};

40 41 42 43 44
/** \brief Node data to use with FCountKernel
 *
 * Holds two members that count the interactions that the node multipole and
 * local development represent.
 */
45
namespace test_count_data {
46 47 48 49
    /// FMM operators names
    enum {P2P, P2M, M2M, M2L, L2L, L2P, P2L, M2P, OP_COUNT};

    /// Multipole equivalent particle count
50 51 52
    struct multipole_t {
        std::size_t dummy = 0;
        std::size_t up = 0;
53
    } m_data;
54

55
    /// Local development particle interaction count
56 57 58
    struct local_expansion_t {
        std::size_t down = 0;
        std::size_t dummy = 0;
59 60
    } l_data;

61 62
};

63 64
using TestCountNodeData = FGenericData<test_count_data::multipole_t, test_count_data::local_expansion_t>;

65

66

67 68 69 70
/** \brief Test Kernel to count particles */
template<class CellClass, class ContainerClass>
struct FCountKernel {

71 72
    using particle_t = typename ContainerClass::particle_t;

73 74
    constexpr static std::size_t child_count = Fpow(2, particle_t::position_t::Dim);

75 76
    std::mutex mtx;
    bool is_main_kernel = true;
77
    int tmp = 0;
78

79
    std::unordered_map<std::string, std::size_t> call_count;
80

81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104
    FCountKernel() = default;
    FCountKernel(const FCountKernel& other) : call_count(other.call_count) {}
    FCountKernel(FCountKernel&& other) : call_count(std::move(other.call_count)) {}

    FCountKernel& operator=(const FCountKernel& other) {
        call_count = other.call_count;
    }
    FCountKernel& operator=(FCountKernel&& other) {
        call_count = std::move(other.call_count);
    }

    void fuse_results(FCountKernel& other) {
        other.is_main_kernel = false;
        for(auto value : other.call_count) {
            this->call_count[value.first] += value.second;
        }
    }

    template<class tree_t>
    void cleanup(tree_t& tree) {
        if(! this->is_main_kernel) {
            return;
        }

105
        using particle_t = typename ContainerClass::particle_t;
106 107 108 109 110 111 112 113 114 115 116 117

        std::cerr << "Checking particle counts\n";

        std::size_t particle_count = 0;
        // Count particles
        for(auto* leaf : tree.leaves()) {
            particle_count += leaf->getParticleContainer()->size();
        }

        for(auto* leaf : tree.leaves()) {
            for(auto p : *(leaf->getParticleContainer())) {
                if(std::get<particle_t::COUNT>(p) != particle_count) {
118
                    std::cerr << scalfmm::tuple_out(p) << "does not have the right particle count\n";
119 120 121 122 123 124 125 126 127 128 129 130 131
                }
            }
        }
    }


    /** \brief Particle to multipole operator
     *
     * Creates a leaf multipole expantion from its particles.
     *
     * \param leaf_multipole Leaf multipole
     * \param source_particle_container Leaf particle container
     */
132
    template<class Symb>
133
    void P2M(typename CellClass::multipole_t* const leaf_multipole,
134
             const Symb* const /* leaf_symbolic_data */,
135 136 137
             const ContainerClass* const source_particle_container)
    {
        mtx.lock();
138
        this->call_count["P2M"] += 1;
139
        mtx.unlock();
140
        leaf_multipole->dummy += 1;
141
        leaf_multipole->up += source_particle_container->size();
142
        leaf_multipole->dummy -= 1;
143 144
    }

145

146
    template<class Symb>
147
    void M2M(typename CellClass::multipole_t * const node_multipole,
148 149 150
             const Symb* const node_symbolic_data ,
             const typename CellClass::multipole_t * const * const child_multipoles,
             const Symb* const /* child_symbolic_data */ [])
151
    {
152
        mtx.lock();
153
        this->call_count["M2M"] += 1;
154
        mtx.unlock();
155
        node_multipole->dummy++;
156
        for(std::size_t idx = 0 ; idx < child_count ; ++idx) {
157 158 159 160 161
            if(child_multipoles[idx]) {
                node_multipole->up += child_multipoles[idx]->up;
                if(child_multipoles[idx]->dummy) {
                    std::cerr << "P2M/M2M parallelism overlap... lvl:" << node_symbolic_data->getLevel() << "\n";
                    std::cerr << child_multipoles[idx]->dummy << '\n';
162
                }
163
            }
164
        }
165
        node_multipole->dummy--;
166 167
    }

168

169
    template<class Symb>
170
    void M2L(typename CellClass::local_expansion_t* const node_local_expansion,
171
             const Symb* const /*node_symbolic_data*/,
172
             const typename CellClass::multipole_t* const v_item_multipoles[],
173
             const Symb* const /*v_item_symbolic_data*/[],
174
             const int /*position*/[],
175
             const int v_item_data_size)
176 177 178 179 180 181 182
    {
        mtx.lock();
        this->call_count["M2L"] += v_item_data_size;
        mtx.unlock();
        for(int idx = 0 ; idx < v_item_data_size ; ++idx) {
            if(v_item_multipoles[idx]) {
                node_local_expansion->down += v_item_multipoles[idx]->up;
183
            }
184 185 186
        }
    }

187

188
    template<class Symb>
189
    void L2L(const typename CellClass::local_expansion_t* const node_local_exp,
190
             const Symb* const /*node_symbolic_data*/,
191
             typename CellClass::local_expansion_t** const child_local_exps,
192
             const Symb* const /*child_symbolic_data*/[])
193
    {
194
        mtx.lock();
195
        this->call_count["L2L"] += 1;
196
        mtx.unlock();
197
        for(std::size_t idx = 0 ; idx < child_count ; ++idx) {
198 199
            if(child_local_exps[idx]) {
                child_local_exps[idx]->down += node_local_exp->down;
200
            }
201 202 203 204
        }

    }

205
    template<class Symb>
206
    void P2L(typename CellClass::local_expansion_t* const node_local_expansion,
207
             const Symb* const /* node_symbolic_data */,
208 209
             const ContainerClass* const source_particle_container,
             const Symb* const /* leaf_symbolic_data */)
210 211
    {
        mtx.lock();
212
        this->call_count["P2L"] += 1;
213 214
        mtx.unlock();
        node_local_expansion->down += source_particle_container->size();
215 216
    }

217

218
    template<class Symb>
219
    void L2P(const typename CellClass::local_expansion_t* const leaf_local_exp,
220
             const Symb* const /* target_symbolic_data */,
221 222 223
             ContainerClass* const target_particle_container)
    {
        mtx.lock();
224
        this->call_count["L2P"] += 1;
225
        mtx.unlock();
226
        for(auto&& particle_ref_tuple : *target_particle_container) {
227
            std::get<particle_t::COUNT>(particle_ref_tuple) += leaf_local_exp->down;
228 229 230
        }
    }

231

232
    template<class Symb>
233
    void M2P(const typename CellClass::multipole_t* const node_multipole,
234
             const Symb* const /* node_symbolic_data */,
235 236
             ContainerClass* const target_particle_container,
             const Symb* const /* leaf_symbolic_data */)
237 238
    {
        mtx.lock();
239
        this->call_count["M2P"] += 1;
240
        mtx.unlock();
241
        for(auto&& particle_ref_tuple : *target_particle_container) {
242
            std::get<particle_t::COUNT>(particle_ref_tuple) += node_multipole->up;
243 244 245 246 247 248 249 250
        }
    }

    /** \brief Particle to particle
     *
     * Direct computation of interactions between particles of a leaf and those
     * of its neighbours.
     *
251 252 253 254
     * \note This overload computes the interactions between the cell and itself in
     * addition to the given neighbours.
     *
     * \param c Unused, coordinates of the target node
255 256 257
     * \param target_particle_container Target node target particle container
     * \param source_particle_container Target node source particle container
     * \param u_item_source_particle_container Source particle container pointer array
258
     * \param positions Unused, source particle container array positions relative to current leaf
259 260
     * \param u_item_count u_item_source_particle_container size
     */
261
    void P2P(const FTreeCoordinate& c,
262 263 264
             ContainerClass* const target_particle_container,
             const ContainerClass* const source_particle_container,
             ContainerClass* const u_item_source_particle_container[],
265 266 267 268 269 270 271 272 273 274 275 276
             const int positions[],
             const int u_item_count
        )
    {
        P2P(c, target_particle_container, source_particle_container,
            u_item_source_particle_container, positions,
            u_item_count, true);
    }

    static bool NeedFinishedM2LEvent() {return false;}
    void finishedLevelM2L(int){}

277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292
    /**
     * \brief Particle to particle
     *
     * Direct computation of interactions between particles of a leaf and those
     * of its neighbours.
     *
     * \note This overload computes the interactions between the cell and its given
     * neighbours.
     *
     * \param c Unused, coordinates of the target node
     * \param target_particle_container Target node target particle container
     * \param source_particle_container Target node source particle container
     * \param u_item_source_particle_container Source particle container pointer array
     * \param positions Unused, source particle container array positions relative to current leaf
     * \param u_item_count u_item_source_particle_container size
     */
293 294 295 296
    void P2P(const FTreeCoordinate&,
             ContainerClass* const target_particle_container,
             const ContainerClass* const  source_particle_container,
             ContainerClass* const u_item_source_particle_container[],
297
             const int /*positions*/[],
298
             const int u_item_count,
299
             bool do_inner
300
        )
301
    {
302
        if(do_inner) {
303
            mtx.lock();
304
            this->call_count["P2P"] += 1;
305
            mtx.unlock();
306
        } else {
307
            mtx.lock();
308
            this->call_count["partial P2P"] += 1;
309
            mtx.unlock();
310 311
        }

312
        for(auto&& particle_ref_tuple : *target_particle_container) {
313 314 315 316
            if(do_inner) {
                std::get<particle_t::COUNT>(particle_ref_tuple) +=
                    source_particle_container->size();
            }
317 318 319 320 321 322 323 324 325 326 327 328 329 330

            for(int i = 0; i < u_item_count; ++i) {
                std::get<particle_t::COUNT>(particle_ref_tuple) +=
                    (u_item_source_particle_container[i])->size();
            }
        }
    }

};




#endif /* FCOUNT_KERNEL_HPP */