Une nouvelle version du portail de gestion des comptes externes sera mise en production lundi 09 août. Elle permettra d'allonger la validité d'un compte externe jusqu'à 3 ans. Pour plus de détails sur cette version consulter : https://doc-si.inria.fr/x/FCeS

FGroupSeqAlgorithm.hpp 21 KB
Newer Older
1 2 3 4 5 6 7
#ifndef FGROUPSEQALGORITHM_HPP
#define FGROUPSEQALGORITHM_HPP

#include "../Utils/FGlobal.hpp"
#include "../Core/FCoreCommon.hpp"
#include "../Utils/FQuickSort.hpp"
#include "../Containers/FTreeCoordinate.hpp"
8 9
#include "../Utils/FLog.hpp"
#include "../Utils/FTic.hpp"
10 11 12 13

#include <list>
#include <vector>

14
template <class OctreeClass, class CellContainerClass, class CellClass, class KernelClass, class ParticleGroupClass, class ParticleContainerClass>
15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42
class FGroupSeqAlgorithm {
protected:
    struct OutOfBlockInteraction{
        MortonIndex outIndex;
        MortonIndex insideIndex;
        int outPosition;

        operator long long() const{
            return static_cast<long long>(outIndex);
        }
    };

    const int MaxThreads;         //< The number of threads
    OctreeClass*const tree;       //< The Tree
    KernelClass*const kernels;    //< The kernels

public:
    FGroupSeqAlgorithm(OctreeClass*const inTree, KernelClass* inKernels) : MaxThreads(1), tree(inTree), kernels(inKernels){
        FAssertLF(tree, "tree cannot be null");
        FAssertLF(kernels, "kernels cannot be null");

        FLOG(FLog::Controller << "FGroupSeqAlgorithm (Max Thread " << MaxThreads << ")\n");
    }

    ~FGroupSeqAlgorithm(){
    }

    void execute(const unsigned operationsToProceed = FFmmNearAndFarFields){
43 44
        FLOG( FLog::Controller << "\tStart FGroupSeqAlgorithm\n" );

45 46 47 48 49 50 51 52 53 54 55 56 57
        if(operationsToProceed & FFmmP2M) bottomPass();

        if(operationsToProceed & FFmmM2M) upwardPass();

        if(operationsToProceed & FFmmM2L) transferPass();

        if(operationsToProceed & FFmmL2L) downardPass();

        if( (operationsToProceed & FFmmP2P) || (operationsToProceed & FFmmL2P) ) directPass();
    }

protected:
    void bottomPass(){
58
        FLOG( FTic timer; );
59 60
        typename std::list<ParticleGroupClass*>::iterator iterParticles = tree->leavesBegin();
        const typename std::list<ParticleGroupClass*>::iterator endParticles = tree->leavesEnd();
61

62 63
        typename std::list<CellContainerClass*>::iterator iterCells = tree->cellsBegin(tree->getHeight()-1);
        const typename std::list<CellContainerClass*>::iterator endCells = tree->cellsEnd(tree->getHeight()-1);
64 65 66 67 68 69 70 71 72

        while(iterParticles != endParticles && iterCells != endCells){
            { // Can be a task(in:iterParticles, out:iterCells)
                const MortonIndex blockStartIdx = (*iterCells)->getStartingIndex();
                const MortonIndex blockEndIdx = (*iterCells)->getEndingIndex();

                for(MortonIndex mindex = blockStartIdx ; mindex < blockEndIdx ; ++mindex){
                    CellClass* cell = (*iterCells)->getCell(mindex);
                    if(cell){
73
                        ParticleContainerClass particles = (*iterParticles)->template getLeaf<ParticleContainerClass>(mindex);
74 75 76 77 78 79 80 81 82 83 84
                        FAssertLF(particles.isAttachedToSomething());
                        kernels->P2M(cell, &particles);
                    }
                }
            }

            ++iterParticles;
            ++iterCells;
        }

        FAssertLF(iterParticles == endParticles && iterCells == endCells);
85
        FLOG( FLog::Controller << "\t\t bottomPass in " << timer.tacAndElapsed() << "s\n" );
86 87 88
    }

    void upwardPass(){
89
        FLOG( FTic timer; );
90
        for(int idxLevel = tree->getHeight()-2 ; idxLevel >= 2 ; --idxLevel){
91 92
            typename std::list<CellContainerClass*>::iterator iterCells = tree->cellsBegin(idxLevel);
            const typename std::list<CellContainerClass*>::iterator endCells = tree->cellsEnd(idxLevel);
93

94 95
            typename std::list<CellContainerClass*>::iterator iterChildCells = tree->cellsBegin(idxLevel+1);
            const typename std::list<CellContainerClass*>::iterator endChildCells = tree->cellsEnd(idxLevel+1);
96 97 98 99 100 101 102 103 104

            while(iterCells != endCells && iterChildCells != endChildCells){
                { // Can be a task(in:iterParticles, out:iterChildCells ...)
                    const MortonIndex blockStartIdx = (*iterCells)->getStartingIndex();
                    const MortonIndex blockEndIdx = (*iterCells)->getEndingIndex();

                    for(MortonIndex mindex = blockStartIdx ; mindex < blockEndIdx && iterChildCells != endChildCells; ++mindex){
                        CellClass* cell = (*iterCells)->getCell(mindex);
                        if(cell){
105
                            CellClass* child[8] = {nullptr,nullptr,nullptr,nullptr,nullptr,nullptr,nullptr,nullptr};
106 107

                            for(int idxChild = 0 ; idxChild < 8 ; ++idxChild){
BRAMAS Berenger's avatar
BRAMAS Berenger committed
108
                                while(iterChildCells != endChildCells && (*iterChildCells)->getEndingIndex() <= ((mindex<<3)+idxChild) ){
109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124
                                    ++iterChildCells;
                                }
                                if( iterChildCells == endChildCells ){
                                    break;
                                }
                                child[idxChild] = (*iterChildCells)->getCell((mindex<<3)+idxChild);
                            }

                            kernels->M2M(cell, child, idxLevel);
                        }
                    }
                }

                ++iterCells;
            }

BRAMAS Berenger's avatar
BRAMAS Berenger committed
125 126
            FAssertLF(iterCells == endCells);
            FAssertLF((iterChildCells == endChildCells || (++iterChildCells) == endChildCells));
127 128
            FAssertLF(iterCells == endCells && (iterChildCells == endChildCells || (++iterChildCells) == endChildCells));
        }
129
        FLOG( FLog::Controller << "\t\t upwardPass in " << timer.tacAndElapsed() << "s\n" );
130 131 132
    }

    void transferPass(){
133
        FLOG( FTic timer; );
134
        for(int idxLevel = tree->getHeight()-1 ; idxLevel >= 2 ; --idxLevel){
135 136
            typename std::list<CellContainerClass*>::iterator iterCells = tree->cellsBegin(idxLevel);
            const typename std::list<CellContainerClass*>::iterator endCells = tree->cellsEnd(idxLevel);
137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152

            while(iterCells != endCells){
                std::vector<OutOfBlockInteraction> outsideInteractions;

                { // Can be a task(inout:iterCells, out:outsideInteractions)
                    const MortonIndex blockStartIdx = (*iterCells)->getStartingIndex();
                    const MortonIndex blockEndIdx = (*iterCells)->getEndingIndex();

                    for(MortonIndex mindex = blockStartIdx ; mindex < blockEndIdx ; ++mindex){
                        CellClass* cell = (*iterCells)->getCell(mindex);
                        if(cell){
                            MortonIndex interactionsIndexes[189];
                            int interactionsPosition[189];
                            FTreeCoordinate coord(mindex, idxLevel);
                            int counter = coord.getInteractionNeighbors(idxLevel,interactionsIndexes,interactionsPosition);

153
                            const CellClass* interactions[343];
154 155 156 157 158 159 160 161 162 163 164
                            memset(interactions, 0, 343*sizeof(CellClass*));
                            int counterExistingCell = 0;

                            for(int idxInter = 0 ; idxInter < counter ; ++idxInter){
                                if( blockStartIdx <= interactionsIndexes[idxInter] && interactionsIndexes[idxInter] < blockEndIdx ){
                                    CellClass* interCell = (*iterCells)->getCell(interactionsIndexes[idxInter]);
                                    if(interCell){
                                        interactions[interactionsPosition[idxInter]] = interCell;
                                        counterExistingCell += 1;
                                    }
                                }
165
                                else if(interactionsIndexes[idxInter] < mindex){
166 167 168 169 170 171 172 173 174 175 176 177 178 179 180
                                    OutOfBlockInteraction property;
                                    property.insideIndex = mindex;
                                    property.outIndex    = interactionsIndexes[idxInter];
                                    property.outPosition = interactionsPosition[idxInter];
                                    outsideInteractions.push_back(property);
                                }
                            }

                            kernels->M2L( cell , interactions, counterExistingCell, idxLevel);
                        }
                    }
                }


                // Manage outofblock interaction
181
                FQuickSort<OutOfBlockInteraction, long long, int>::QsSequential(outsideInteractions.data(),int(outsideInteractions.size()));
182

183
                typename std::list<CellContainerClass*>::iterator iterLeftCells = tree->cellsBegin(idxLevel);
184
                int currentOutInteraction = 0;
185
                while(iterLeftCells != iterCells && currentOutInteraction < int(outsideInteractions.size())){
186 187 188 189 190 191 192
                    const MortonIndex blockStartIdx = (*iterLeftCells)->getStartingIndex();
                    const MortonIndex blockEndIdx = (*iterLeftCells)->getEndingIndex();

                    while(outsideInteractions[currentOutInteraction].outIndex < blockStartIdx){
                        currentOutInteraction += 1;
                    }

193
                    int lastOutInteraction = currentOutInteraction;
194
                    while(lastOutInteraction < int(outsideInteractions.size()) && outsideInteractions[lastOutInteraction].outIndex < blockEndIdx){
195 196 197 198 199 200 201 202 203
                        lastOutInteraction += 1;
                    }

                    { // Can be a task(in:currentOutInteraction, in:outsideInteractions, in:lastOutInteraction, inout:iterLeftCells, inout:iterCells)
                        for(int outInterIdx = currentOutInteraction ; outInterIdx < lastOutInteraction ; ++outInterIdx){
                            CellClass* interCell = (*iterLeftCells)->getCell(outsideInteractions[outInterIdx].outIndex);
                            if(interCell){
                                CellClass* cell = (*iterCells)->getCell(outsideInteractions[outInterIdx].insideIndex);
                                FAssertLF(cell);
204
                                const CellClass* interactions[343];
205 206 207 208 209
                                memset(interactions, 0, 343*sizeof(CellClass*));
                                interactions[outsideInteractions[outInterIdx].outPosition] = interCell;
                                const int counter = 1;
                                kernels->M2L( cell , interactions, counter, idxLevel);

210
                                interactions[outsideInteractions[outInterIdx].outPosition] = nullptr;
211 212 213 214 215 216 217 218 219 220 221 222 223 224
                                interactions[getOppositeInterIndex(outsideInteractions[outInterIdx].outPosition)] = cell;
                                kernels->M2L( interCell , interactions, counter, idxLevel);
                            }
                        }
                    }

                    currentOutInteraction = lastOutInteraction;
                    ++iterLeftCells;
                }

                ++iterCells;
            }

        }
225
        FLOG( FLog::Controller << "\t\t transferPass in " << timer.tacAndElapsed() << "s\n" );
226 227 228
    }

    void downardPass(){
229
        FLOG( FTic timer; );
230
        for(int idxLevel = 2 ; idxLevel <= tree->getHeight()-2 ; ++idxLevel){
231 232
            typename std::list<CellContainerClass*>::iterator iterCells = tree->cellsBegin(idxLevel);
            const typename std::list<CellContainerClass*>::iterator endCells = tree->cellsEnd(idxLevel);
233

234 235
            typename std::list<CellContainerClass*>::iterator iterChildCells = tree->cellsBegin(idxLevel+1);
            const typename std::list<CellContainerClass*>::iterator endChildCells = tree->cellsEnd(idxLevel+1);
236 237 238 239 240 241 242 243 244

            while(iterCells != endCells && iterChildCells != endChildCells){
                { // Can be a task(in:iterParticles, inout:iterChildCells ...)
                    const MortonIndex blockStartIdx = (*iterCells)->getStartingIndex();
                    const MortonIndex blockEndIdx = (*iterCells)->getEndingIndex();

                    for(MortonIndex mindex = blockStartIdx ; mindex < blockEndIdx && iterChildCells != endChildCells; ++mindex){
                        CellClass* cell = (*iterCells)->getCell(mindex);
                        if(cell){
245
                            CellClass* child[8] = {nullptr,nullptr,nullptr,nullptr,nullptr,nullptr,nullptr,nullptr};
246 247

                            for(int idxChild = 0 ; idxChild < 8 ; ++idxChild){
BRAMAS Berenger's avatar
BRAMAS Berenger committed
248
                                while(iterChildCells != endChildCells && (*iterChildCells)->getEndingIndex() <= ((mindex<<3)+idxChild) ){
249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266
                                    ++iterChildCells;
                                }
                                if( iterChildCells == endChildCells ){
                                    break;
                                }
                                child[idxChild] = (*iterChildCells)->getCell((mindex<<3)+idxChild);
                            }

                            kernels->L2L(cell, child, idxLevel);
                        }
                    }
                }

                ++iterCells;
            }

            FAssertLF(iterCells == endCells && (iterChildCells == endChildCells || (++iterChildCells) == endChildCells));
        }
267
        FLOG( FLog::Controller << "\t\t downardPass in " << timer.tacAndElapsed() << "s\n" );
268 269 270
    }

    void directPass(){
271
        FLOG( FTic timer; );
272
        {
273 274
            typename std::list<ParticleGroupClass*>::iterator iterParticles = tree->leavesBegin();
            const typename std::list<ParticleGroupClass*>::iterator endParticles = tree->leavesEnd();
275

276 277
            typename std::list<CellContainerClass*>::iterator iterCells = tree->cellsBegin(tree->getHeight()-1);
            const typename std::list<CellContainerClass*>::iterator endCells = tree->cellsEnd(tree->getHeight()-1);
278 279 280 281

            while(iterParticles != endParticles && iterCells != endCells){
                { // Can be a task(in:iterCells, inout:iterParticles)
                    const MortonIndex blockStartIdx = (*iterCells)->getStartingIndex();
BRAMAS Berenger's avatar
BRAMAS Berenger committed
282
                    const MortonIndex blockEndIdx = (*iterCells)->getEndingIndex();
283 284 285 286

                    for(MortonIndex mindex = blockStartIdx ; mindex < blockEndIdx ; ++mindex){
                        CellClass* cell = (*iterCells)->getCell(mindex);
                        if(cell){
287
                            ParticleContainerClass particles = (*iterParticles)->template getLeaf<ParticleContainerClass>(mindex);
288
                            FAssertLF(particles.isAttachedToSomething());
289
                            kernels->L2P(cell, &particles);
290 291 292 293 294 295 296 297 298 299 300
                        }
                    }
                }

                ++iterParticles;
                ++iterCells;
            }

            FAssertLF(iterParticles == endParticles && iterCells == endCells);
        }
        {
301 302
            typename std::list<ParticleGroupClass*>::iterator iterParticles = tree->leavesBegin();
            const typename std::list<ParticleGroupClass*>::iterator endParticles = tree->leavesEnd();
303 304 305 306 307 308 309 310 311

            while(iterParticles != endParticles){
                typename std::vector<OutOfBlockInteraction> outsideInteractions;

                { // Can be a task(inout:iterCells, out:outsideInteractions)
                    const MortonIndex blockStartIdx = (*iterParticles)->getStartingIndex();
                    const MortonIndex blockEndIdx = (*iterParticles)->getEndingIndex();

                    for(MortonIndex mindex = blockStartIdx ; mindex < blockEndIdx ; ++mindex){
312
                        ParticleContainerClass particles = (*iterParticles)->template getLeaf<ParticleContainerClass>(mindex);
313 314 315 316 317 318 319 320
                        if(particles.isAttachedToSomething()){
                            MortonIndex interactionsIndexes[26];
                            int interactionsPosition[26];
                            FTreeCoordinate coord(mindex, tree->getHeight()-1);
                            int counter = coord.getNeighborsIndexes(tree->getHeight(),interactionsIndexes,interactionsPosition);

                            ParticleContainerClass interactionsObjects[27];
                            ParticleContainerClass* interactions[27];
321
                            memset(interactions, 0, 27*sizeof(ParticleContainerClass*));
322 323 324 325
                            int counterExistingCell = 0;

                            for(int idxInter = 0 ; idxInter < counter ; ++idxInter){
                                if( blockStartIdx <= interactionsIndexes[idxInter] && interactionsIndexes[idxInter] < blockEndIdx ){
326
                                    interactionsObjects[counterExistingCell] = (*iterParticles)->template getLeaf<ParticleContainerClass>(interactionsIndexes[idxInter]);
327
                                    if(interactionsObjects[counterExistingCell].isAttachedToSomething()){
BRAMAS Berenger's avatar
BRAMAS Berenger committed
328
                                        FAssertLF(interactions[interactionsPosition[idxInter]] == nullptr);
329 330 331 332
                                        interactions[interactionsPosition[idxInter]] = &interactionsObjects[counterExistingCell];
                                        counterExistingCell += 1;
                                    }
                                }
333
                                else if(interactionsIndexes[idxInter] < mindex){
334 335 336 337 338 339 340 341 342 343 344 345 346 347 348
                                    OutOfBlockInteraction property;
                                    property.insideIndex = mindex;
                                    property.outIndex    = interactionsIndexes[idxInter];
                                    property.outPosition = interactionsPosition[idxInter];
                                    outsideInteractions.push_back(property);
                                }
                            }

                            kernels->P2P( coord, &particles, &particles , interactions, counterExistingCell);
                        }
                    }
                }


                // Manage outofblock interaction
349
                FQuickSort<OutOfBlockInteraction, long long, int>::QsSequential(outsideInteractions.data(),int(outsideInteractions.size()));
350

351
                typename std::list<ParticleGroupClass*>::iterator iterLeftParticles = tree->leavesBegin();
352
                int currentOutInteraction = 0;
353
                while(iterLeftParticles != iterParticles && currentOutInteraction < int(outsideInteractions.size())){
354 355 356 357 358 359 360
                    const MortonIndex blockStartIdx = (*iterLeftParticles)->getStartingIndex();
                    const MortonIndex blockEndIdx = (*iterLeftParticles)->getEndingIndex();

                    while(outsideInteractions[currentOutInteraction].outIndex < blockStartIdx){
                        currentOutInteraction += 1;
                    }

361
                    int lastOutInteraction = currentOutInteraction;
362
                    while(lastOutInteraction < int(outsideInteractions.size()) && outsideInteractions[lastOutInteraction].outIndex < blockEndIdx){
363 364 365 366 367
                        lastOutInteraction += 1;
                    }

                    { // Can be a task(in:currentOutInteraction, in:outsideInteractions, in:lastOutInteraction, inout:iterLeftParticles, inout:iterParticles)
                        for(int outInterIdx = currentOutInteraction ; outInterIdx < lastOutInteraction ; ++outInterIdx){
368
                            ParticleContainerClass interParticles = (*iterLeftParticles)->template getLeaf<ParticleContainerClass>(outsideInteractions[outInterIdx].outIndex);
369
                            if(interParticles.isAttachedToSomething()){
370
                                ParticleContainerClass particles = (*iterParticles)->template getLeaf<ParticleContainerClass>(outsideInteractions[outInterIdx].insideIndex);
371
                                FAssertLF(particles.isAttachedToSomething());
372
                                ParticleContainerClass* interactions[27];
373
                                memset(interactions, 0, 27*sizeof(ParticleContainerClass*));
374 375 376 377
                                interactions[outsideInteractions[outInterIdx].outPosition] = &interParticles;
                                const int counter = 1;
                                kernels->P2PRemote( FTreeCoordinate(outsideInteractions[outInterIdx].insideIndex, tree->getHeight()-1), &particles, &particles , interactions, counter);

378
                                interactions[outsideInteractions[outInterIdx].outPosition] = nullptr;
379 380 381 382 383 384 385 386 387 388 389 390 391
                                interactions[getOppositeNeighIndex(outsideInteractions[outInterIdx].outPosition)] = &particles;
                                kernels->P2PRemote( FTreeCoordinate(outsideInteractions[outInterIdx].outIndex, tree->getHeight()-1), &interParticles, &interParticles , interactions, counter);
                            }
                        }
                    }

                    currentOutInteraction = lastOutInteraction;
                    ++iterLeftParticles;
                }

                ++iterParticles;
            }
        }
392
        FLOG( FLog::Controller << "\t\t directPass in " << timer.tacAndElapsed() << "s\n" );
393 394 395 396
    }

    int getOppositeNeighIndex(const int index) const {
        // ((idxX+1)*3 + (idxY+1)) * 3 + (idxZ+1)
397
        return 27-index-1;
398 399 400 401
    }

    int getOppositeInterIndex(const int index) const {
        // ((( (xdiff+3) * 7) + (ydiff+3))) * 7 + zdiff + 3
402
        return 343-index-1;
403 404 405 406 407
    }
};


#endif // FGROUPSEQALGORITHM_HPP