milupHPC documentation
  • include
  • subdomain_key_tree
tree.cuh
Go to the documentation of this file.
1
12#ifndef MILUPHPC_TREE_CUH
13#define MILUPHPC_TREE_CUH
14
15#include "../cuda_utils/cuda_utilities.cuh"
16#include "../parameter.h"
17#include "../particles.cuh"
18
19#include <iostream>
20#include <stdio.h>
21#include <cuda.h>
22#include <assert.h>
23#include <cmath>
24
25namespace KeyNS {
26
30#if DIM == 1
31 CUDA_CALLABLE_MEMBER const unsigned char DirTable[1][1] = {{1}}; //TODO: 1D DirTable?
32#elif DIM == 2
33 CUDA_CALLABLE_MEMBER const unsigned char DirTable[4][4] =
34 {{1,2,0,0},
35 {0,1,3,1},
36 {2,0,2,3},
37 {3,3,1,2}};
38#else DIM == 3
39#ifndef __CUDACC__
40 const unsigned char DirTable[12][8] =
41 {{8, 10, 3, 3, 4, 5, 4, 5},
42 {2, 2, 11, 9, 4, 5, 4, 5},
43 {7, 6, 7, 6, 8, 10, 1, 1},
44 {7, 6, 7, 6, 0, 0, 11, 9},
45 {0, 8, 1, 11, 6, 8, 6, 11},
46 {10, 0, 9, 1, 10, 7, 9, 7},
47 {10, 4, 9, 4, 10, 2, 9, 3},
48 {5, 8, 5, 11, 2, 8, 3, 11},
49 {4, 9, 0, 0, 7, 9, 2, 2},
50 {1, 1, 8, 5, 3, 3, 8, 6},
51 {11, 5, 0, 0, 11, 6, 2, 2},
52 {1, 1, 4, 10, 3, 3, 7, 10}};
53#else
54 __device__ const unsigned char DirTable[12][8] =
55 {{8, 10, 3, 3, 4, 5, 4, 5},
56 {2, 2, 11, 9, 4, 5, 4, 5},
57 {7, 6, 7, 6, 8, 10, 1, 1},
58 {7, 6, 7, 6, 0, 0, 11, 9},
59 {0, 8, 1, 11, 6, 8, 6, 11},
60 {10, 0, 9, 1, 10, 7, 9, 7},
61 {10, 4, 9, 4, 10, 2, 9, 3},
62 {5, 8, 5, 11, 2, 8, 3, 11},
63 {4, 9, 0, 0, 7, 9, 2, 2},
64 {1, 1, 8, 5, 3, 3, 8, 6},
65 {11, 5, 0, 0, 11, 6, 2, 2},
66 {1, 1, 4, 10, 3, 3, 7, 10}};
67#endif
68#endif
69
73#if DIM == 1
74 CUDA_CALLABLE_MEMBER const unsigned char HilbertTable[1][1] = {{1}}; //TODO: 1D HilbertTable?
75#elif DIM == 2
76 CUDA_CALLABLE_MEMBER const unsigned char HilbertTable[4][4] =
77 {{0,3,1,2},
78 {0,1,3,2},
79 {2,3,1,0},
80 {2,1,3,0}};
81#else
82#ifndef __CUDACC__
83 const unsigned char HilbertTable[12][8] = {{0, 7, 3, 4, 1, 6, 2, 5},
84 {4, 3, 7, 0, 5, 2, 6, 1},
85 {6, 1, 5, 2, 7, 0, 4, 3},
86 {2, 5, 1, 6, 3, 4, 0, 7},
87 {0, 1, 7, 6, 3, 2, 4, 5},
88 {6, 7, 1, 0, 5, 4, 2, 3},
89 {2, 3, 5, 4, 1, 0, 6, 7},
90 {4, 5, 3, 2, 7, 6, 0, 1},
91 {0, 3, 1, 2, 7, 4, 6, 5},
92 {2, 1, 3, 0, 5, 6, 4, 7},
93 {4, 7, 5, 6, 3, 0, 2, 1},
94 {6, 5, 7, 4, 1, 2, 0, 3}};
95#else
96 __device__ const unsigned char HilbertTable[12][8] = {{0, 7, 3, 4, 1, 6, 2, 5},
97 {4, 3, 7, 0, 5, 2, 6, 1},
98 {6, 1, 5, 2, 7, 0, 4, 3},
99 {2, 5, 1, 6, 3, 4, 0, 7},
100 {0, 1, 7, 6, 3, 2, 4, 5},
101 {6, 7, 1, 0, 5, 4, 2, 3},
102 {2, 3, 5, 4, 1, 0, 6, 7},
103 {4, 5, 3, 2, 7, 6, 0, 1},
104 {0, 3, 1, 2, 7, 4, 6, 5},
105 {2, 1, 3, 0, 5, 6, 4, 7},
106 {4, 7, 5, 6, 3, 0, 2, 1},
107 {6, 5, 7, 4, 1, 2, 0, 3}};
108#endif
109#endif
110
118 CUDA_CALLABLE_MEMBER keyType lebesgue2hilbert(keyType lebesgue, integer maxLevel);
119
120 CUDA_CALLABLE_MEMBER keyType lebesgue2hilbert(keyType lebesgue, int maxLevel, int level);
121
122}
123
140class Tree {
141
142public:
143
144 //TODO: count, start, sorted currently unused (since sort() and computeForces() described by Griebel not used!)
146 integer *count;
148 integer *start;
150 integer *child;
152 integer *sorted;
154 integer *index;
155
157 integer *toDeleteLeaf;
159 integer *toDeleteNode;
160
162 real *minX;
164 real *maxX;
165#if DIM > 1
167 real *minY;
169 real *maxY;
170#if DIM == 3
172 real *minZ;
174 real *maxZ;
175#endif
176#endif
177
181 CUDA_CALLABLE_MEMBER Tree();
182
196 CUDA_CALLABLE_MEMBER Tree(integer *count, integer *start, integer *child, integer *sorted, integer *index,
197 integer *toDeleteLeaf, integer *toDeleteNode,
198 real *minX, real *maxX);
212 CUDA_CALLABLE_MEMBER void set(integer *count, integer *start, integer *child, integer *sorted,
213 integer *index, integer *toDeleteLeaf, integer *toDeleteNode,
214 real *minX, real *maxX);
215
216#if DIM > 1
232 CUDA_CALLABLE_MEMBER Tree(integer *count, integer *start, integer *child, integer *sorted, integer *index,
233 integer *toDeleteLeaf, integer *toDeleteNode,
234 real *minX, real *maxX, real *minY, real *maxY);
250 CUDA_CALLABLE_MEMBER void set(integer *count, integer *start, integer *child, integer *sorted,
251 integer *index, integer *toDeleteLeaf, integer *toDeleteNode,
252 real *minX, real *maxX, real *minY, real *maxY);
253
254#if DIM == 3
272 CUDA_CALLABLE_MEMBER Tree(integer *count, integer *start, integer *child, integer *sorted, integer *index,
273 integer *toDeleteLeaf, integer *toDeleteNode,
274 real *minX, real *maxX, real *minY, real *maxY, real *minZ, real *maxZ);
292 CUDA_CALLABLE_MEMBER void set(integer *count, integer *start, integer *child, integer *sorted,
293 integer *index, integer *toDeleteLeaf, integer *toDeleteNode,
294 real *minX, real *maxX, real *minY, real *maxY,
295 real *minZ, real *maxZ);
296#endif
297#endif
298
305 CUDA_CALLABLE_MEMBER void reset(integer index, integer n);
306
321 CUDA_CALLABLE_MEMBER keyType getParticleKey(Particles *particles, integer index, integer maxLevel,
322 Curve::Type curveType = Curve::lebesgue);
323
336 CUDA_CALLABLE_MEMBER integer getTreeLevel(Particles *particles, integer index, integer maxLevel,
337 Curve::Type curveType = Curve::lebesgue);
338
344 CUDA_CALLABLE_MEMBER integer sumParticles();
345
349 CUDA_CALLABLE_MEMBER ~Tree();
350
351
352};
353
354namespace TreeNS {
355
356 namespace Kernel {
357
374 __global__ void set(Tree *tree, integer *count, integer *start, integer *child, integer *sorted,
375 integer *index, integer *toDeleteLeaf, integer *toDeleteNode,
376 real *minX, real *maxX);
377
387 __global__ void info(Tree *tree, Particles *particles, integer n, integer m);
388
389 __global__ void info(Tree *tree, Particles *particles);
390
391 __global__ void testTree(Tree *tree, Particles *particles, integer n, integer m);
392
393 namespace Launch {
397 void set(Tree *tree, integer *count, integer *start, integer *child, integer *sorted, integer *index,
398 integer *toDeleteLeaf, integer *toDeleteNode,
399 real *minX, real *maxX);
400
404 real info(Tree *tree, Particles *particles, integer n, integer m);
405
406 real info(Tree *tree, Particles *particles);
407
408 real testTree(Tree *tree, Particles *particles, integer n, integer m);
409 }
410
411#if DIM > 1
428 __global__ void set(Tree *tree, integer *count, integer *start, integer *child, integer *sorted,
429 integer *index, integer *toDeleteLeaf, integer *toDeleteNode,
430 real *minX, real *maxX, real *minY, real *maxY);
431
432 namespace Launch {
449 void set(Tree *tree, integer *count, integer *start, integer *child, integer *sorted,
450 integer *index, integer *toDeleteLeaf, integer *toDeleteNode,
451 real *minX, real *maxX, real *minY, real *maxY);
452 }
453
454#if DIM == 3
455
474 __global__ void set(Tree *tree, integer *count, integer *start, integer *child, integer *sorted,
475 integer *index, integer *toDeleteLeaf, integer *toDeleteNode,
476 real *minX, real *maxX, real *minY, real *maxY, real *minZ,
477 real *maxZ);
478
479 namespace Launch {
498 void set(Tree *tree, integer *count, integer *start, integer *child, integer *sorted,
499 integer *index, integer *toDeleteLeaf, integer *toDeleteNode,
500 real *minX, real *maxX, real *minY, real *maxY, real *minZ,
501 real *maxZ);
502 }
503
504#endif
505#endif
506
514 __global__ void sumParticles(Tree *tree);
515
536 __global__ void buildTree(Tree *tree, Particles *particles, integer n, integer m);
537
538 __global__ void prepareSorting(Tree *tree, Particles *particles, integer n, integer m);
539
552 __global__ void calculateCentersOfMass(Tree *tree, Particles *particles, integer n, integer level);
553
568 __global__ void computeBoundingBox(Tree *tree, Particles *particles, integer *mutex,
569 integer n, integer blockSize);
570
580 __global__ void centerOfMass(Tree *tree, Particles *particles, integer n);
581
591 __global__ void sort(Tree *tree, integer n, integer m);
592
605 __global__ void getParticleKeys(Tree *tree, Particles *particles, keyType *keys, integer maxLevel, integer n,
606 Curve::Type curveType = Curve::lebesgue);
607
617 __global__ void globalCOM(Tree *tree, Particles *particles, real com[DIM]);
618
619 namespace Launch {
620
626 real sumParticles(Tree *tree);
627
633 real getParticleKeys(Tree *tree, Particles *particles, keyType *keys, integer maxLevel, integer n,
634 Curve::Type curveType = Curve::lebesgue, bool time=false);
635
641 real buildTree(Tree *tree, Particles *particles, integer n, integer m, bool time=false);
642
643 real prepareSorting(Tree *tree, Particles *particles, integer n, integer m);
644
645 real calculateCentersOfMass(Tree *tree, Particles *particles, integer n, integer level, bool time=false);
646
652 real computeBoundingBox(Tree *tree, Particles *particles, integer *mutex, integer n,
653 integer blockSize, bool time=false);
654
660 real centerOfMass(Tree *tree, Particles *particles, integer n, bool time=false);
661
667 real sort(Tree *tree, integer n, integer m, bool time=false);
668
674 real globalCOM(Tree *tree, Particles *particles, real com[DIM]);
675
676 }
677 }
678}
679
680#if UNIT_TESTING
681namespace UnitTesting {
682 namespace Kernel {
683
684 __global__ void test_localTree(Tree *tree, Particles *particles, int n, int m);
685
686 namespace Launch {
687 real test_localTree(Tree *tree, Particles *particles, int n, int m);
688 }
689 }
690}
691#endif
692
693#endif //MILUPHPC_TREE_CUH
Particles
Particle(s) class based on SoA (Structur of Arrays).
Definition: particles.cuh:50
Tree
Tree class.
Definition: tree.cuh:140
Tree::reset
CUDA_CALLABLE_MEMBER void reset(integer index, integer n)
Reset (specific) entries.
Definition: tree.cu:133
Tree::index
integer * index
index used for counting nodes
Definition: tree.cuh:154
Tree::minZ
real * minZ
bounding box minimal z
Definition: tree.cuh:172
Tree::child
integer * child
children/child nodes or leaves (beneath)
Definition: tree.cuh:150
Tree::getParticleKey
CUDA_CALLABLE_MEMBER keyType getParticleKey(Particles *particles, integer index, integer maxLevel, Curve::Type curveType=Curve::lebesgue)
Get SFC key of a particle.
Definition: tree.cu:157
Tree::maxX
real * maxX
bounding box maximal x
Definition: tree.cuh:164
Tree::sorted
integer * sorted
sorted (indices) for better cache performance
Definition: tree.cuh:152
Tree::count
integer * count
accumulated nodes/leaves beneath
Definition: tree.cuh:146
Tree::sumParticles
CUDA_CALLABLE_MEMBER integer sumParticles()
Sum particles in tree.
Definition: tree.cu:288
Tree::toDeleteNode
integer * toDeleteNode
buffer for remembering old indices for rebuilding tree
Definition: tree.cuh:159
Tree::~Tree
CUDA_CALLABLE_MEMBER ~Tree()
Destructor.
Definition: tree.cu:298
Tree::Tree
CUDA_CALLABLE_MEMBER Tree()
Default constructor.
Definition: tree.cu:55
Tree::minX
real * minX
bounding box minimal x
Definition: tree.cuh:162
Tree::start
integer * start
TODO: describe start.
Definition: tree.cuh:148
Tree::toDeleteLeaf
integer * toDeleteLeaf
buffer for remembering old indices for rebuilding tree
Definition: tree.cuh:157
Tree::getTreeLevel
CUDA_CALLABLE_MEMBER integer getTreeLevel(Particles *particles, integer index, integer maxLevel, Curve::Type curveType=Curve::lebesgue)
Get tree level for a (specific) particle.
Definition: tree.cu:242
Tree::set
CUDA_CALLABLE_MEMBER void set(integer *count, integer *start, integer *child, integer *sorted, integer *index, integer *toDeleteLeaf, integer *toDeleteNode, real *minX, real *maxX)
Setter, passing pointer to member variables.
Definition: tree.cu:65
Tree::maxY
real * maxY
bounding box maximal y
Definition: tree.cuh:169
Tree::minY
real * minY
bounding box minimal y
Definition: tree.cuh:167
Tree::maxZ
real * maxZ
bounding box maximal z
Definition: tree.cuh:174
CUDA_CALLABLE_MEMBER
#define CUDA_CALLABLE_MEMBER
Definition: cuda_utilities.cuh:30
Kernel
Definition: device_rhs.cuh:7
KeyNS
Key (keyType) related functions and kernels.
Definition: subdomain.cuh:35
KeyNS::HilbertTable
const unsigned char HilbertTable[12][8]
Table needed to convert from Lebesgue to Hilbert keys.
Definition: tree.cuh:83
KeyNS::DIM
DIM
Table needed to convert from Lebesgue to Hilbert keys.
Definition: tree.cuh:38
KeyNS::lebesgue2hilbert
CUDA_CALLABLE_MEMBER keyType lebesgue2hilbert(keyType lebesgue, integer maxLevel)
Convert a Lebesgue key to a Hilbert key.
Definition: tree.cu:4
ProfilerIds::Time::tree
const char *const tree
Definition: h5profiler.h:57
TreeNS::Kernel::Launch::calculateCentersOfMass
real calculateCentersOfMass(Tree *tree, Particles *particles, integer n, integer level, bool time=false)
Definition: tree.cu:1944
TreeNS::Kernel::Launch::prepareSorting
real prepareSorting(Tree *tree, Particles *particles, integer n, integer m)
Definition: tree.cu:1939
TreeNS::Kernel::Launch::sort
real sort(Tree *tree, integer n, integer m, bool time=false)
Wrapper for TreeNS::Kernel::sort()
Definition: tree.cu:1962
TreeNS::Kernel::Launch::centerOfMass
real centerOfMass(Tree *tree, Particles *particles, integer n, bool time=false)
Wrapper for TreeNS::Kernel::centerOfMass()
Definition: tree.cu:1957
TreeNS::Kernel::Launch::sumParticles
real sumParticles(Tree *tree)
Wrapper for TreeNS::Kernel::sumParticles()
Definition: tree.cu:1929
TreeNS::Kernel::Launch::globalCOM
real globalCOM(Tree *tree, Particles *particles, real com[DIM])
Wrapper for TreeNS::Kernel::globalCOM()
Definition: tree.cu:1974
TreeNS::Kernel::Launch::testTree
real testTree(Tree *tree, Particles *particles, integer n, integer m)
Definition: tree.cu:1886
TreeNS::Kernel::Launch::computeBoundingBox
real computeBoundingBox(Tree *tree, Particles *particles, integer *mutex, integer n, integer blockSize, bool time=false)
Wrapper for TreeNS::Kernel::computeBoundingBox()
Definition: tree.cu:1949
TreeNS::Kernel::Launch::set
void set(Tree *tree, integer *count, integer *start, integer *child, integer *sorted, integer *index, integer *toDeleteLeaf, integer *toDeleteNode, real *minX, real *maxX)
Wrapper for TreeNS::Kernel::set()
Definition: tree.cu:1869
TreeNS::Kernel::Launch::getParticleKeys
real getParticleKeys(Tree *tree, Particles *particles, keyType *keys, integer maxLevel, integer n, Curve::Type curveType=Curve::lebesgue, bool time=false)
Wrapper for TreeNS::Kernel::getParticleKeys()
Definition: tree.cu:1967
TreeNS::Kernel::Launch::info
real info(Tree *tree, Particles *particles, integer n, integer m)
Wrapper for TreeNS::Kernel::Launch::info()
Definition: tree.cu:1876
TreeNS::Kernel::Launch::buildTree
real buildTree(Tree *tree, Particles *particles, integer n, integer m, bool time=false)
Wrapper for TreeNS::Kernel::buildTree()
Definition: tree.cu:1934
TreeNS::Kernel::prepareSorting
__global__ void prepareSorting(Tree *tree, Particles *particles, integer n, integer m)
Definition: tree.cu:1404
TreeNS::Kernel::set
__global__ void set(Tree *tree, integer *count, integer *start, integer *child, integer *sorted, integer *index, integer *toDeleteLeaf, integer *toDeleteNode, real *minX, real *maxX)
Kernel call to setter.
Definition: tree.cu:1679
TreeNS::Kernel::testTree
__global__ void testTree(Tree *tree, Particles *particles, integer n, integer m)
Definition: tree.cu:1815
TreeNS::Kernel::info
__global__ void info(Tree *tree, Particles *particles, integer n, integer m)
Info Kernel for tree class (for debugging purposes)
Definition: tree.cu:1684
TreeNS::Kernel::getParticleKeys
__global__ void getParticleKeys(Tree *tree, Particles *particles, keyType *keys, integer maxLevel, integer n, Curve::Type curveType=Curve::lebesgue)
Kernel to get all particle's keys.
Definition: tree.cu:1618
TreeNS::Kernel::computeBoundingBox
__global__ void computeBoundingBox(Tree *tree, Particles *particles, integer *mutex, integer n, integer blockSize)
Kernel to compute the bounding box/borders of the tree or rather the particles within the tree.
Definition: tree.cu:302
TreeNS::Kernel::calculateCentersOfMass
__global__ void calculateCentersOfMass(Tree *tree, Particles *particles, integer n, integer level)
Compute center of masses (level wise).
Definition: tree.cu:1416
TreeNS::Kernel::sumParticles
__global__ void sumParticles(Tree *tree)
Kernel call to sum particles within tree.
Definition: tree.cu:440
TreeNS::Kernel::sort
__global__ void sort(Tree *tree, integer n, integer m)
Kernel to sort tree/child indices to optimize cache efficiency.
Definition: tree.cu:1527
TreeNS::Kernel::buildTree
__global__ void buildTree(Tree *tree, Particles *particles, integer n, integer m)
Kernel to construct the tree using the particles within particles
Definition: tree.cu:768
TreeNS::Kernel::globalCOM
__global__ void globalCOM(Tree *tree, Particles *particles, real com[DIM])
Compute center of mass for all particles.
Definition: tree.cu:1645
TreeNS::Kernel::centerOfMass
__global__ void centerOfMass(Tree *tree, Particles *particles, integer n)
Definition: tree.cu:1501
TreeNS
Definition: tree.cuh:354
real
double real
Definition: parameter.h:15
keyType
unsigned long keyType
Definition: parameter.h:18
integer
int integer
Definition: parameter.h:17
Curve::Type
Type
Definition: parameter.h:206
Curve::lebesgue
@ lebesgue
Definition: parameter.h:207

milupHPC - include/subdomain_key_tree/tree.cuh Source File
Generated on Wed Aug 31 2022 12:16:52 by Doxygen 1.9.3