1
/*============================================================================
2
* Search octrees and quadtrees of boxes.
3
*============================================================================*/
6
This file is part of Code_Saturne, a general-purpose CFD tool.
8
Copyright (C) 1998-2011 EDF S.A.
10
This program is free software; you can redistribute it and/or modify it under
11
the terms of the GNU General Public License as published by the Free Software
12
Foundation; either version 2 of the License, or (at your option) any later
15
This program is distributed in the hope that it will be useful, but WITHOUT
16
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17
FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
20
You should have received a copy of the GNU General Public License along with
21
this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
22
Street, Fifth Floor, Boston, MA 02110-1301, USA.
25
/*----------------------------------------------------------------------------*/
27
#if defined(HAVE_CONFIG_H)
28
#include "cs_config.h"
31
/*----------------------------------------------------------------------------
32
* Standard C library headers
33
*----------------------------------------------------------------------------*/
40
/*----------------------------------------------------------------------------
42
*----------------------------------------------------------------------------*/
45
#include <bft_printf.h>
47
/*----------------------------------------------------------------------------
49
*----------------------------------------------------------------------------*/
51
#include "fvm_config_defs.h"
52
#include "fvm_parall.h"
54
#include "fvm_box_priv.h"
56
/*----------------------------------------------------------------------------
57
* Header for the current file
58
*----------------------------------------------------------------------------*/
60
#include "fvm_box_tree.h"
62
/*----------------------------------------------------------------------------*/
67
} /* Fake brace to force Emacs auto-indentation back to column 0 */
69
#endif /* __cplusplus */
71
/*=============================================================================
72
* Local Macro and Type definitions
73
*============================================================================*/
75
#define FVM_BOX_TREE_MAX_BUILD_LOOPS 50
77
/* Structures for each octant or quadrant */
78
/*----------------------------------------*/
80
/* If the type is BOX_TREE_NODE, the ordering of children is defined as follows,
81
using notation B: bottom, U: up, E: east, W: west, S: south, N: north.
83
octant: 0: BSW, 1: BSE, 2: BNW, 3: BNE, 4: USW, 5: USE, 6: UNW, 7: UNE
84
quadrant: 0: SW, 1: SE, 2: NW, 3: NE
90
_Bool is_leaf; /* True for leaf nodes */
92
fvm_morton_code_t morton_code; /* Level and coordinates in the grid
93
according to Morton encoding */
95
fvm_lnum_t n_boxes; /* Number of associated bounding boxes */
96
fvm_lnum_t start_id; /* Position of the first box_id */
100
/* Structure used to manage statistics */
104
unsigned max_level_reached; /* Max level number reached */
106
fvm_lnum_t n_leaves; /* Number of leaves in the tree */
107
fvm_lnum_t n_boxes; /* Number of boxes to locate in the tree */
108
fvm_lnum_t n_linked_boxes; /* Number of linked boxes in the tree */
109
fvm_lnum_t n_spill_leaves; /* Number of leaves where n_boxes > threshold */
111
fvm_lnum_t min_linked_boxes; /* Minimum number of boxes for a leaf */
112
fvm_lnum_t max_linked_boxes; /* Maximum number of boxes for a leaf */
114
} fvm_box_tree_stats_t;
116
/* Main box tree structure */
117
/*-------------------------*/
119
struct _fvm_box_tree_t {
121
int n_children; /* 8, 4, or 2 (2^dim) */
123
fvm_morton_int_t max_level; /* Max. possible level */
124
fvm_lnum_t threshold; /* Max number of boxes linked to a
125
node if max_level is not reached */
126
float max_box_ratio; /* Max n_linked_boxes / n_boxes value */
128
fvm_box_tree_stats_t stats; /* Statistics related to the structure */
130
fvm_lnum_t n_max_nodes; /* Current max. allocated nodes */
131
fvm_lnum_t n_nodes; /* Number of nodes (including leaves) */
133
_node_t *nodes; /* Array of nodes (root at index 0) */
135
fvm_lnum_t *child_ids; /* Ids of associated children
136
(size: 2^dim * n_max_nodes) */
137
fvm_lnum_t *box_ids; /* List of associated box ids.
138
size = stat.n_linked_boxes */
140
int n_build_loops; /* Number of loops required to build */
142
#if defined(HAVE_MPI)
143
MPI_Comm comm; /* Associated MPI communicator */
147
/*=============================================================================
148
* Static global variables
149
*============================================================================*/
151
/*============================================================================
152
* Private function definitions
153
*============================================================================*/
155
/*----------------------------------------------------------------------------
156
* Get minimum coordinates for a given box.
159
* box_set <-- pointer to box set structure
160
* box_id <-- id of box
163
* pointer to minimum box coordinates
164
*---------------------------------------------------------------------------*/
166
inline static const fvm_coord_t *
167
_box_min(const fvm_box_set_t *boxes,
170
return boxes->extents + box_id*boxes->dim*2;
173
/*----------------------------------------------------------------------------
174
* Get maxmum coordinates for a given box.
177
* box_set <-- pointer to box set structure
178
* box_id <-- id of box
181
* pointer to maximum box coordinates
182
*---------------------------------------------------------------------------*/
184
inline static const fvm_coord_t *
185
_box_max(const fvm_box_set_t *boxes,
188
return boxes->extents + box_id*boxes->dim*2 + boxes->dim;
191
/*----------------------------------------------------------------------------
192
* Test for intersection between two bounding boxes.
195
* extents <-- array of box extents
196
* id_0 <-- id of first box
197
* id_1 <-- id of second box
201
*---------------------------------------------------------------------------*/
204
_boxes_intersect_3d(const fvm_coord_t *extents,
208
const fvm_coord_t *e0 = extents + id_0*6;
209
const fvm_coord_t *e1 = extents + id_1*6;
211
if ( e0[0] > e1[3] || e1[0] > e0[3]
212
|| e0[1] > e1[4] || e1[1] > e0[4]
213
|| e0[2] > e1[5] || e1[2] > e0[5])
220
_boxes_intersect_2d(const fvm_coord_t *extents,
224
const fvm_coord_t *e0 = extents + id_0*4;
225
const fvm_coord_t *e1 = extents + id_1*4;
227
if ( e0[0] > e1[2] || e1[0] > e0[2]
228
|| e0[1] > e1[3] || e1[1] > e0[3])
235
_boxes_intersect_1d(const fvm_coord_t *extents,
239
const fvm_coord_t *e0 = extents + id_0*2;
240
const fvm_coord_t *e1 = extents + id_1*2;
242
if ( e0[0] > e1[1] || e1[0] > e0[1])
248
/*----------------------------------------------------------------------------
249
* Update octree stat structure (min, max, mean, box ratio, ...)
252
* bt <-> pointer on the fvm_box_tree_t structure to deal with
253
* node_id <-- node on which we collect data
254
*----------------------------------------------------------------------------*/
257
_update_tree_stats(fvm_box_tree_t *bt,
262
const _node_t *node = bt->nodes + node_id;
264
if (node->is_leaf == false) {
265
int n_children = bt->n_children;
266
const fvm_lnum_t *_child_ids = bt->child_ids + node_id*bt->n_children;
267
for (i = 0; i < n_children; i++)
268
_update_tree_stats(bt, _child_ids[i]);
271
else { /* leaf node */
273
fvm_box_tree_stats_t s = bt->stats;
276
s.n_linked_boxes += node->n_boxes;
278
if (node->n_boxes > bt->threshold)
279
s.n_spill_leaves += 1;
281
s.min_linked_boxes = FVM_MIN(s.min_linked_boxes, node->n_boxes);
282
s.max_linked_boxes = FVM_MAX(s.max_linked_boxes, node->n_boxes);
283
s.max_level_reached = FVM_MAX(s.max_level_reached, node->morton_code.L);
289
/*----------------------------------------------------------------------------
290
* Define box_tree->stat structure (min, max, mean, box ratio, ...)
293
* bt <-> pointer to the box-tree structure
294
*----------------------------------------------------------------------------*/
297
_get_box_tree_stats(fvm_box_tree_t *bt)
302
/* Initialize statistics */
304
bt->stats.max_level_reached = 0;
306
bt->stats.n_leaves = 0;
307
bt->stats.n_spill_leaves = 0;
308
bt->stats.n_linked_boxes = 0;
310
bt->stats.min_linked_boxes = INT_MAX;
311
bt->stats.max_linked_boxes = 0;
313
/* Recursively update stats, starting from root */
315
if (bt->nodes != NULL)
316
_update_tree_stats(bt, 0);
319
/*----------------------------------------------------------------------------
320
* Get the coordinates in the grid for the current point at this level.
323
* level <-- level on which we want the coordinates
324
* coords <-- coords of the point to translate in the octree grid.
325
* XYZ <-> pointer to the X, Y, Z coordinates in the grid
326
*----------------------------------------------------------------------------*/
329
_get_grid_coords_3d(fvm_morton_int_t level,
330
const double coords[3],
333
fvm_morton_int_t refinement = 1 << level;
335
XYZ[0] = coords[0] * refinement;
336
XYZ[1] = coords[1] * refinement;
337
XYZ[2] = coords[2] * refinement;
341
_get_grid_coords_2d(fvm_morton_int_t level,
342
const double coords[2],
345
fvm_morton_int_t refinement = 1 << level;
347
XY[0] = coords[0] * refinement;
348
XY[1] = coords[1] * refinement;
352
_get_grid_coords_1d(fvm_morton_int_t level,
353
const double coords[1],
356
fvm_morton_int_t refinement = 1 << level;
358
X[0] = coords[0] * refinement;
361
/*----------------------------------------------------------------------------
362
* Return true if a leaf intersects a box, false otherwise.
365
* morton_code <-- Morton code of the leaf
366
* min_box <-- coordinates of min. point of the bounding box
367
* max_box <-- coordinates of max. point of the bounding box
371
*----------------------------------------------------------------------------*/
374
_node_intersect_box_3d(fvm_morton_code_t morton_code,
375
const fvm_coord_t min_box[3],
376
const fvm_coord_t max_box[3])
379
double min_oct[3], max_oct[3];
381
for (i = 0; i < 3; i++) {
382
min_oct[i] = (double)morton_code.X[i];
383
max_oct[i] = (double)(morton_code.X[i] + 1);
386
if ( min_box[0] > max_oct[0] || min_oct[0] > max_box[0]
387
|| min_box[1] > max_oct[1] || min_oct[1] > max_box[1]
388
|| min_box[2] > max_oct[2] || min_oct[2] > max_box[2])
395
_node_intersect_box_2d(fvm_morton_code_t morton_code,
396
const fvm_coord_t min_box[2],
397
const fvm_coord_t max_box[2])
400
double min_oct[2], max_oct[2];
402
for (i = 0; i < 2; i++) {
403
min_oct[i] = (double)morton_code.X[i];
404
max_oct[i] = (double)(morton_code.X[i] + 1);
407
if ( min_box[0] > max_oct[0] || min_oct[0] > max_box[0]
408
|| min_box[1] > max_oct[1] || min_oct[1] > max_box[1])
415
_node_intersect_box_1d(fvm_morton_code_t morton_code,
416
const fvm_coord_t min_box[1],
417
const fvm_coord_t max_box[1])
419
double min_oct, max_oct;
421
min_oct = (double)morton_code.X[0];
422
max_oct = (double)(morton_code.X[0] + 1);
424
if (min_box[0] > max_oct || min_oct > max_box[0])
430
/*----------------------------------------------------------------------------
431
* Split a node into its children and evaluate the box distribution.
434
* bt <-> pointer to the box tree being built
435
* boxes <-- pointer to the associated box set structure
436
* node_id <-- id of the node to split
439
* the number of associations between nodes and their children
440
*----------------------------------------------------------------------------*/
443
_evaluate_splitting_3d(fvm_box_tree_t *bt,
444
const fvm_box_set_t *boxes,
448
fvm_morton_code_t min_code, max_code;
449
fvm_morton_code_t children[8];
451
int n_linked_boxes = 0;
453
const _node_t node = bt->nodes[node_id];
454
const fvm_morton_int_t next_level = node.morton_code.L + 1;
456
assert(boxes->dim == 3);
458
/* Define a Morton code for each child */
460
fvm_morton_get_children(3, node.morton_code, children);
462
/* Loop on boxes associated to the node_id */
464
for (j = 0; j < node.n_boxes; j++) {
466
fvm_coord_t min_grid_coord[3], max_grid_coord[3];
468
fvm_lnum_t box_id = bt->box_ids[node.start_id + j];
469
const fvm_coord_t *box_min = _box_min(boxes, box_id);
470
const fvm_coord_t *box_max = _box_max(boxes, box_id);
472
min_code = fvm_morton_encode(3, next_level, box_min);
473
max_code = fvm_morton_encode(3, next_level, box_max);
475
if ( fvm_morton_compare(3, min_code, max_code)
476
== FVM_MORTON_DIFFERENT_ID) {
478
_get_grid_coords_3d(next_level, box_min, min_grid_coord);
479
_get_grid_coords_3d(next_level, box_max, max_grid_coord);
481
for (i = 0; i < 8; i++) {
482
if (_node_intersect_box_3d(children[i], min_grid_coord, max_grid_coord))
487
else { /* Box is included in the same octant */
489
assert( fvm_morton_compare(3, max_code, min_code)
490
== FVM_MORTON_EQUAL_ID);
492
for (i = 0; i < 8; i++) {
493
if ( fvm_morton_compare(3, min_code, children[i])
494
== FVM_MORTON_EQUAL_ID) {
498
} /* End of loop on children */
500
} /* If min_code and max_code in the same leaf */
502
} /* End of loop on boxes */
504
return n_linked_boxes;
508
_evaluate_splitting_2d(fvm_box_tree_t *bt,
509
const fvm_box_set_t *boxes,
513
fvm_morton_code_t min_code, max_code;
514
fvm_morton_code_t children[4];
516
int n_linked_boxes = 0;
518
const _node_t node = bt->nodes[node_id];
519
const fvm_morton_int_t next_level = node.morton_code.L + 1;
521
assert(boxes->dim == 2);
523
/* Define a Morton code for each child */
525
fvm_morton_get_children(2, node.morton_code, children);
527
/* Loop on boxes associated to the node_id */
529
for (j = 0; j < node.n_boxes; j++) {
531
fvm_coord_t min_grid_coord[2], max_grid_coord[2];
533
fvm_lnum_t box_id = bt->box_ids[node.start_id + j];
534
const fvm_coord_t *box_min = _box_min(boxes, box_id);
535
const fvm_coord_t *box_max = _box_max(boxes, box_id);
537
min_code = fvm_morton_encode(2, next_level, box_min);
538
max_code = fvm_morton_encode(2, next_level, box_max);
540
if ( fvm_morton_compare(2, min_code, max_code)
541
== FVM_MORTON_DIFFERENT_ID) {
543
_get_grid_coords_2d(next_level, box_min, min_grid_coord);
544
_get_grid_coords_2d(next_level, box_max, max_grid_coord);
546
for (i = 0; i < 4; i++) {
547
if (_node_intersect_box_2d(children[i], min_grid_coord, max_grid_coord))
552
else { /* Box is included in the same quadrant */
554
assert( fvm_morton_compare(2, max_code, min_code)
555
== FVM_MORTON_EQUAL_ID);
557
for (i = 0; i < 4; i++) {
558
if ( fvm_morton_compare(2, min_code, children[i])
559
== FVM_MORTON_EQUAL_ID) {
564
} /* End of loop on children */
566
} /* If min_code and max_code in the same leaf */
568
} /* End of loop on boxes */
570
return n_linked_boxes;
574
_evaluate_splitting_1d(fvm_box_tree_t *bt,
575
const fvm_box_set_t *boxes,
579
fvm_morton_code_t min_code, max_code;
580
fvm_morton_code_t children[2];
582
int n_linked_boxes = 0;
584
const _node_t node = bt->nodes[node_id];
585
const fvm_morton_int_t next_level = node.morton_code.L + 1;
587
/* Define a Morton code for each child */
589
fvm_morton_get_children(1, node.morton_code, children);
591
/* Loop on boxes associated to the node_id */
593
for (j = 0; j < node.n_boxes; j++) {
595
fvm_coord_t min_grid_coord[1], max_grid_coord[1];
597
fvm_lnum_t box_id = bt->box_ids[node.start_id + j];
598
const fvm_coord_t *box_min = _box_min(boxes, box_id);
599
const fvm_coord_t *box_max = _box_max(boxes, box_id);
601
min_code = fvm_morton_encode(1, next_level, box_min);
602
max_code = fvm_morton_encode(1, next_level, box_max);
604
if ( fvm_morton_compare(1, min_code, max_code)
605
== FVM_MORTON_DIFFERENT_ID) {
607
_get_grid_coords_1d(next_level, box_min, min_grid_coord);
608
_get_grid_coords_1d(next_level, box_max, max_grid_coord);
610
for (i = 0; i < 2; i++) {
611
if (_node_intersect_box_1d(children[i], min_grid_coord, max_grid_coord))
616
else { /* Box is included in the same quadrant */
618
assert( fvm_morton_compare(1, max_code, min_code)
619
== FVM_MORTON_EQUAL_ID);
621
for (i = 0; i < 2; i++) {
622
if ( fvm_morton_compare(1, min_code, children[i])
623
== FVM_MORTON_EQUAL_ID) {
627
} /* End of loop on children */
629
} /* If min_code and max_code in the same leaf */
631
} /* End of loop on boxes */
633
return n_linked_boxes;
636
/*----------------------------------------------------------------------------
637
* Evaluate the box distribution over the leaves of the box tree to help
638
* determine if we should add a level to the tree structure.
641
* bt <-> pointer to the box tree being built
642
* boxes <-- pointer to the associated box set structure
643
* node_id <-- id of the starting node
644
* build_type <-- layout variant for building the tree structure
645
* next_level_size --> size of box_ids for the next level
646
*----------------------------------------------------------------------------*/
649
_count_next_level(fvm_box_tree_t *bt,
650
const fvm_box_set_t *boxes,
652
fvm_box_tree_sync_t build_type,
653
fvm_lnum_t *next_level_size)
657
_node_t *node = bt->nodes + node_id;
659
if (node->is_leaf == false) {
661
assert(bt->child_ids[bt->n_children*node_id] > 0);
663
for (i = 0; i < bt->n_children; i++)
664
_count_next_level(bt,
666
bt->child_ids[bt->n_children*node_id + i],
672
else { /* if (node->is_leaf == true) */
673
if ( node->n_boxes < bt->threshold
674
&& node_id != 0 /* Root node is always divided */
675
&& build_type == FVM_BOX_TREE_ASYNC_LEVEL)
676
*next_level_size += node->n_boxes;
678
else { /* Split node and evaluate box distribution between its children */
680
*next_level_size += _evaluate_splitting_3d(bt, boxes, node_id);
681
else if (boxes->dim == 2)
682
*next_level_size += _evaluate_splitting_2d(bt, boxes, node_id);
683
else if (boxes->dim == 1)
684
*next_level_size += _evaluate_splitting_1d(bt, boxes, node_id);
689
/*----------------------------------------------------------------------------
690
* Test if we have to continue the building of the box tree.
693
* bt <-> pointer to the box tree being built
694
* boxes <-- pointer to the associated box set structure
695
* build_type <-- layout variant for building the tree structure
696
* next_size --> size of box_ids for the next tree if required
699
* true if we should continue, false otherwise.
700
*----------------------------------------------------------------------------*/
703
_recurse_tree_build(fvm_box_tree_t *bt,
704
const fvm_box_set_t *boxes,
705
fvm_box_tree_sync_t build_type,
706
fvm_lnum_t *next_size)
709
fvm_lnum_t _next_size = 0;
711
_Bool retval = false;
713
#if defined(HAVE_MPI)
716
MPI_Comm comm = boxes->comm;
718
if (comm != MPI_COMM_NULL)
719
MPI_Comm_size(comm, &n_ranks);
723
bt->n_build_loops += 1;
728
/* To avoid infinite loop on tree building */
730
if (bt->n_build_loops > FVM_BOX_TREE_MAX_BUILD_LOOPS)
733
/* A sufficient accuracy has been reached */
735
if (bt->stats.max_level_reached == bt->max_level)
738
/* Algorithm is converged. No need to go further */
740
if ( bt->stats.n_spill_leaves == 0
741
&& bt->stats.max_level_reached > 0)
744
#if defined(HAVE_MPI)
745
if (n_ranks > 1 && build_type == FVM_BOX_TREE_SYNC_LEVEL) {
747
MPI_Allreduce(&state, &global_state, 1, MPI_INT, MPI_MIN, comm);
748
state = global_state; /* Stop if all ranks require it */
756
/* Limit, to avoid excessive memory usage */
758
_count_next_level(bt,
760
0, /* Starts from root */
764
if (bt->stats.n_boxes > 0)
765
box_ratio = (_next_size*1.0)/bt->stats.n_boxes;
769
if (bt->stats.max_level_reached > 0 && box_ratio > bt->max_box_ratio)
774
#if defined(HAVE_MPI)
775
if (n_ranks > 1 && build_type == FVM_BOX_TREE_SYNC_LEVEL) {
777
MPI_Allreduce(&state, &global_state, 1, MPI_INT, MPI_MAX, comm);
778
state = global_state; /* Stop as as soon as any rank requires it */
782
/* If no condition is encoutered, we have to continue */
784
*next_size = _next_size;
792
/*----------------------------------------------------------------------------
793
* Create a box tree by copying from another.
796
* dest <-> pointer to destination box tree
797
* src <-- pointer to source box tree
798
*----------------------------------------------------------------------------*/
801
_copy_tree(fvm_box_tree_t *dest,
802
const fvm_box_tree_t *src)
804
assert(dest != NULL && src != NULL);
806
memcpy(dest, src, sizeof(fvm_box_tree_t));
808
BFT_MALLOC(dest->nodes, dest->n_max_nodes, _node_t);
809
BFT_MALLOC(dest->child_ids, dest->n_max_nodes*dest->n_children, fvm_lnum_t);
810
BFT_MALLOC(dest->box_ids, (dest->stats).n_linked_boxes, fvm_lnum_t);
812
memcpy(dest->nodes, src->nodes, dest->n_nodes * sizeof(_node_t));
813
memcpy(dest->child_ids,
815
dest->n_nodes * src->n_children * sizeof(fvm_lnum_t));
817
memcpy(dest->box_ids,
819
(dest->stats).n_linked_boxes * sizeof(fvm_lnum_t));
822
/*----------------------------------------------------------------------------
823
* Destroy a fvm_box_tree_t structure.
826
* bt <-- pointer to pointer to fvm_box_tree_t structure to destroy
827
*----------------------------------------------------------------------------*/
830
_free_tree_arrays(fvm_box_tree_t *bt)
835
BFT_FREE(bt->child_ids);
836
BFT_FREE(bt->box_ids);
839
/*----------------------------------------------------------------------------
840
* Create a new node from a Morton code and associate it to a tree.
843
* bt <-> pointer to the box tree being built
844
* morton_code <-- Morton identification number related to the node
845
* node_id <-- id of the starting node
846
*----------------------------------------------------------------------------*/
849
_new_node(fvm_box_tree_t *bt,
850
fvm_morton_code_t morton_code,
858
node = bt->nodes + node_id;
860
if (morton_code.L > bt->max_level)
861
bft_error(__FILE__, __LINE__, 0,
862
_("Error adding a new node in box tree (%p).\n"
863
"Max level reached. Current level: %u and Max level: %d\n"),
864
bt, morton_code.L, bt->max_level);
866
node->is_leaf = true;
867
node->morton_code = morton_code;
870
node->start_id = -1; /* invalid value by default */
872
for (i = 0; i < bt->n_children; i++)
873
bt->child_ids[node_id*bt->n_children + i] = -1;
876
/*----------------------------------------------------------------------------
877
* Split a node into its children and define the new box distribution.
880
* bt <-> pointer to the box tree being built
881
* next_bt <-> pointer to the next box tree being built
882
* boxes <-- pointer to the associated box set structure
883
* node_id <-- id of the starting node
884
* shift_ids <-> first free position free in new box_ids
885
*----------------------------------------------------------------------------*/
888
_split_node_3d(fvm_box_tree_t *bt,
889
fvm_box_tree_t *next_bt,
890
const fvm_box_set_t *boxes,
892
fvm_lnum_t *shift_ids)
895
fvm_morton_code_t min_code, max_code;
896
fvm_morton_code_t children[8];
898
fvm_lnum_t n_linked_boxes = 0;
899
fvm_lnum_t _shift_ids = *shift_ids;
900
fvm_lnum_t n_init_nodes = next_bt->n_nodes;
901
_node_t split_node = next_bt->nodes[node_id];
903
const _node_t node = bt->nodes[node_id];
904
const fvm_morton_int_t next_level = node.morton_code.L + 1;
906
assert(bt->n_children == 8);
908
/* Add the leaves to the next_bt structure */
910
if (n_init_nodes + 8 > next_bt->n_max_nodes) {
911
assert(next_bt->n_max_nodes > 0);
912
next_bt->n_max_nodes *= 2;
913
BFT_REALLOC(next_bt->nodes, next_bt->n_max_nodes, _node_t);
914
BFT_REALLOC(next_bt->child_ids, next_bt->n_max_nodes*8, fvm_lnum_t);
917
/* Define a Morton code for each child and create the children nodes */
919
fvm_morton_get_children(3, node.morton_code, children);
921
for (i = 0; i < 8; i++) {
922
const fvm_lnum_t new_id = n_init_nodes + i;
923
next_bt->child_ids[node_id*8 + i] = new_id;
924
_new_node(next_bt, children[i], new_id);
927
split_node.start_id = 0;
928
split_node.n_boxes = 0;
929
split_node.is_leaf = false;
931
next_bt->nodes[node_id] = split_node;
932
next_bt->n_nodes = n_init_nodes + 8;
934
/* Counting loop on boxes associated to the node_id */
936
for (j = 0; j < node.n_boxes; j++) {
938
fvm_coord_t min_grid_coord[3], max_grid_coord[3];
940
fvm_lnum_t box_id = bt->box_ids[node.start_id + j];
941
const fvm_coord_t *box_min = _box_min(boxes, box_id);
942
const fvm_coord_t *box_max = _box_max(boxes, box_id);
944
min_code = fvm_morton_encode(3, next_level, box_min);
945
max_code = fvm_morton_encode(3, next_level, box_max);
947
if ( fvm_morton_compare(3, min_code, max_code)
948
== FVM_MORTON_DIFFERENT_ID) {
950
_get_grid_coords_3d(next_level, box_min, min_grid_coord);
951
_get_grid_coords_3d(next_level, box_max, max_grid_coord);
953
for (i = 0; i < 8; i++) {
954
if (_node_intersect_box_3d(children[i], min_grid_coord, max_grid_coord))
955
(next_bt->nodes[n_init_nodes + i]).n_boxes += 1;
959
else { /* Box is included in the same octant */
961
assert( fvm_morton_compare(3, max_code, min_code)
962
== FVM_MORTON_EQUAL_ID);
964
for (i = 0; i < 8; i++) {
965
if ( fvm_morton_compare(3, min_code, children[i])
966
== FVM_MORTON_EQUAL_ID) {
967
(next_bt->nodes[n_init_nodes + i]).n_boxes += 1;
970
} /* End of loop on children */
972
} /* If min_code and max_code in the same leaf */
974
} /* End of loop on boxes */
978
for (i = 0; i < 8; i++) {
980
(next_bt->nodes[n_init_nodes + i]).start_id
981
= _shift_ids + n_linked_boxes;
982
n_linked_boxes += (next_bt->nodes[n_init_nodes + i]).n_boxes;
986
_shift_ids += n_linked_boxes;
988
for (i = 0; i < 8; i++)
989
(next_bt->nodes[n_init_nodes + i]).n_boxes = 0;
991
/* Second loop on boxes associated to the node_id: fill */
993
for (j = 0; j < node.n_boxes; j++) {
995
fvm_coord_t min_grid_coord[3], max_grid_coord[3];
997
fvm_lnum_t box_id = bt->box_ids[node.start_id + j];
998
const fvm_coord_t *box_min = _box_min(boxes, box_id);
999
const fvm_coord_t *box_max = _box_max(boxes, box_id);
1001
min_code = fvm_morton_encode(3, next_level, box_min);
1002
max_code = fvm_morton_encode(3, next_level, box_max);
1004
if ( fvm_morton_compare(3, min_code, max_code)
1005
== FVM_MORTON_DIFFERENT_ID) {
1007
_get_grid_coords_3d(next_level, box_min, min_grid_coord);
1008
_get_grid_coords_3d(next_level, box_max, max_grid_coord);
1010
for (i = 0; i < 8; i++) {
1012
if (_node_intersect_box_3d(children[i],
1016
const fvm_lnum_t sub_id = n_init_nodes + i;
1017
const fvm_lnum_t shift = (next_bt->nodes[sub_id]).n_boxes
1018
+ (next_bt->nodes[sub_id]).start_id;
1019
next_bt->box_ids[shift] = box_id;
1020
(next_bt->nodes[sub_id]).n_boxes += 1;
1023
} /* End of loop on children*/
1026
else { /* Box is included in the same octant */
1028
assert( fvm_morton_compare(3, max_code, min_code)
1029
== FVM_MORTON_EQUAL_ID);
1031
for (i = 0; i < 8; i++) {
1033
if ( fvm_morton_compare(3, min_code, children[i])
1034
== FVM_MORTON_EQUAL_ID) {
1036
const fvm_lnum_t sub_id = n_init_nodes + i;
1037
const fvm_lnum_t shift = (next_bt->nodes[sub_id]).n_boxes
1038
+ (next_bt->nodes[sub_id]).start_id;
1040
next_bt->box_ids[shift] = box_id;
1041
(next_bt->nodes[sub_id]).n_boxes += 1;
1045
} /* End of loop on children */
1047
} /* If min_code and max_code in the same leaf */
1049
} /* End of loop on boxes */
1051
#if 0 && defined(DEBUG) && !defined(NDEBUG)
1052
for (i = 1; i < 8; i++) {
1053
_node_t n1 = next_bt->nodes[n_init_nodes + i - 1];
1054
_node_t n2 = next_bt->nodes[n_init_nodes + i];
1055
assert(n1.n_boxes == (n2.start_id - n1.start_id));
1058
== ( next_bt->nodes[n_init_nodes + 8 - 1].start_id
1059
+ next_bt->nodes[n_init_nodes + 8 - 1].n_boxes));
1062
/* Return pointers */
1064
*shift_ids = _shift_ids;
1068
_split_node_2d(fvm_box_tree_t *bt,
1069
fvm_box_tree_t *next_bt,
1070
const fvm_box_set_t *boxes,
1072
fvm_lnum_t *shift_ids)
1075
fvm_morton_code_t min_code, max_code;
1076
fvm_morton_code_t children[4];
1078
fvm_lnum_t n_linked_boxes = 0;
1079
fvm_lnum_t _shift_ids = *shift_ids;
1080
fvm_lnum_t n_init_nodes = next_bt->n_nodes;
1081
_node_t split_node = next_bt->nodes[node_id];
1083
const _node_t node = bt->nodes[node_id];
1084
const fvm_morton_int_t next_level = node.morton_code.L + 1;
1086
assert(bt->n_children == 4);
1088
/* Add the leaves to the next_bt structure */
1090
if (n_init_nodes + 4 > next_bt->n_max_nodes) {
1091
assert(next_bt->n_max_nodes > 0);
1092
next_bt->n_max_nodes *= 2;
1093
BFT_REALLOC(next_bt->nodes, next_bt->n_max_nodes, _node_t);
1094
BFT_REALLOC(next_bt->child_ids, next_bt->n_max_nodes*4, fvm_lnum_t);
1097
/* Define a Morton code for each child and create the children nodes */
1099
fvm_morton_get_children(2, node.morton_code, children);
1101
for (i = 0; i < 4; i++) {
1102
const fvm_lnum_t new_id = n_init_nodes + i;
1103
next_bt->child_ids[node_id*4 + i] = new_id;
1104
_new_node(next_bt, children[i], new_id);
1107
split_node.start_id = 0;
1108
split_node.n_boxes = 0;
1109
split_node.is_leaf = false;
1111
next_bt->nodes[node_id] = split_node;
1112
next_bt->n_nodes = n_init_nodes + 4;
1114
/* Counting loop on boxes associated to the node_id */
1116
for (j = 0; j < node.n_boxes; j++) {
1118
fvm_coord_t min_grid_coord[2], max_grid_coord[2];
1120
fvm_lnum_t box_id = bt->box_ids[node.start_id + j];
1121
const fvm_coord_t *box_min = _box_min(boxes, box_id);
1122
const fvm_coord_t *box_max = _box_max(boxes, box_id);
1124
min_code = fvm_morton_encode(2, next_level, box_min);
1125
max_code = fvm_morton_encode(2, next_level, box_max);
1127
if ( fvm_morton_compare(2, min_code, max_code)
1128
== FVM_MORTON_DIFFERENT_ID) {
1130
_get_grid_coords_2d(next_level, box_min, min_grid_coord);
1131
_get_grid_coords_2d(next_level, box_max, max_grid_coord);
1133
for (i = 0; i < 4; i++) {
1134
if (_node_intersect_box_2d(children[i], min_grid_coord, max_grid_coord))
1135
(next_bt->nodes[n_init_nodes + i]).n_boxes += 1;
1139
else { /* Box is included in the same quadrant */
1141
assert( fvm_morton_compare(2, max_code, min_code)
1142
== FVM_MORTON_EQUAL_ID);
1144
for (i = 0; i < 4; i++) {
1145
if ( fvm_morton_compare(2, min_code, children[i])
1146
== FVM_MORTON_EQUAL_ID) {
1147
(next_bt->nodes[n_init_nodes + i]).n_boxes += 1;
1150
} /* End of loop on children */
1152
} /* If min_code and max_code in the same leaf */
1154
} /* End of loop on boxes */
1158
for (i = 0; i < 4; i++) {
1160
(next_bt->nodes[n_init_nodes + i]).start_id
1161
= _shift_ids + n_linked_boxes;
1162
n_linked_boxes += (next_bt->nodes[n_init_nodes + i]).n_boxes;
1166
_shift_ids += n_linked_boxes;
1168
for (i = 0; i < 4; i++)
1169
(next_bt->nodes[n_init_nodes + i]).n_boxes = 0;
1171
/* Second loop on boxes associated to the node_id: fill */
1173
for (j = 0; j < node.n_boxes; j++) {
1175
fvm_coord_t min_grid_coord[2], max_grid_coord[2];
1177
fvm_lnum_t box_id = bt->box_ids[node.start_id + j];
1178
const fvm_coord_t *box_min = _box_min(boxes, box_id);
1179
const fvm_coord_t *box_max = _box_max(boxes, box_id);
1181
min_code = fvm_morton_encode(2, next_level, box_min);
1182
max_code = fvm_morton_encode(2, next_level, box_max);
1184
if ( fvm_morton_compare(2, min_code, max_code)
1185
== FVM_MORTON_DIFFERENT_ID) {
1187
_get_grid_coords_2d(next_level, box_min, min_grid_coord);
1188
_get_grid_coords_2d(next_level, box_max, max_grid_coord);
1190
for (i = 0; i < 4; i++) {
1192
if (_node_intersect_box_2d(children[i],
1196
const fvm_lnum_t sub_id = n_init_nodes + i;
1197
const fvm_lnum_t shift = (next_bt->nodes[sub_id]).n_boxes
1198
+ (next_bt->nodes[sub_id]).start_id;
1199
next_bt->box_ids[shift] = box_id;
1200
(next_bt->nodes[sub_id]).n_boxes += 1;
1203
} /* End of loop on children*/
1206
else { /* Box is included in the same quadrant */
1208
assert( fvm_morton_compare(2, max_code, min_code)
1209
== FVM_MORTON_EQUAL_ID);
1211
for (i = 0; i < 4; i++) {
1213
if ( fvm_morton_compare(2, min_code, children[i])
1214
== FVM_MORTON_EQUAL_ID) {
1216
const fvm_lnum_t sub_id = n_init_nodes + i;
1217
const fvm_lnum_t shift = (next_bt->nodes[sub_id]).n_boxes
1218
+ (next_bt->nodes[sub_id]).start_id;
1220
next_bt->box_ids[shift] = box_id;
1221
(next_bt->nodes[sub_id]).n_boxes += 1;
1225
} /* End of loop on children */
1227
} /* If min_code and max_code in the same leaf */
1229
} /* End of loop on boxes */
1231
#if 0 && defined(DEBUG) && !defined(NDEBUG)
1232
for (i = 1; i < 4; i++) {
1233
_node_t n1 = next_bt->nodes[n_init_nodes + i - 1];
1234
_node_t n2 = next_bt->nodes[n_init_nodes + i];
1235
assert(n1.n_boxes == (n2.start_id - n1.start_id));
1238
== ( next_bt->nodes[n_init_nodes + 4 - 1].start_id
1239
+ next_bt->nodes[n_init_nodes + 4 - 1].n_boxes));
1242
/* Return pointers */
1244
*shift_ids = _shift_ids;
1248
_split_node_1d(fvm_box_tree_t *bt,
1249
fvm_box_tree_t *next_bt,
1250
const fvm_box_set_t *boxes,
1252
fvm_lnum_t *shift_ids)
1255
fvm_morton_code_t min_code, max_code;
1256
fvm_morton_code_t children[2];
1258
fvm_lnum_t n_linked_boxes = 0;
1259
fvm_lnum_t _shift_ids = *shift_ids;
1260
fvm_lnum_t n_init_nodes = next_bt->n_nodes;
1261
_node_t split_node = next_bt->nodes[node_id];
1263
const _node_t node = bt->nodes[node_id];
1264
const fvm_morton_int_t next_level = node.morton_code.L + 1;
1266
assert(bt->n_children == 2);
1268
/* Add the leaves to the next_bt structure */
1270
if (n_init_nodes + 2 > next_bt->n_max_nodes) {
1271
assert(next_bt->n_max_nodes > 0);
1272
next_bt->n_max_nodes *= 2;
1273
BFT_REALLOC(next_bt->nodes, next_bt->n_max_nodes, _node_t);
1274
BFT_REALLOC(next_bt->child_ids, next_bt->n_max_nodes*2, fvm_lnum_t);
1277
/* Define a Morton code for each child and create the children nodes */
1279
fvm_morton_get_children(1, node.morton_code, children);
1281
for (i = 0; i < 2; i++) {
1282
const fvm_lnum_t new_id = n_init_nodes + i;
1283
next_bt->child_ids[node_id*2 + i] = new_id;
1284
_new_node(next_bt, children[i], new_id);
1287
split_node.start_id = 0;
1288
split_node.n_boxes = 0;
1289
split_node.is_leaf = false;
1291
next_bt->nodes[node_id] = split_node;
1292
next_bt->n_nodes = n_init_nodes + 2;
1294
/* Counting loop on boxes associated to the node_id */
1296
for (j = 0; j < node.n_boxes; j++) {
1298
fvm_coord_t min_grid_coord[2], max_grid_coord[2];
1300
fvm_lnum_t box_id = bt->box_ids[node.start_id + j];
1301
const fvm_coord_t *box_min = _box_min(boxes, box_id);
1302
const fvm_coord_t *box_max = _box_max(boxes, box_id);
1304
min_code = fvm_morton_encode(1, next_level, box_min);
1305
max_code = fvm_morton_encode(1, next_level, box_max);
1307
if ( fvm_morton_compare(1, min_code, max_code)
1308
== FVM_MORTON_DIFFERENT_ID) {
1310
_get_grid_coords_1d(next_level, box_min, min_grid_coord);
1311
_get_grid_coords_1d(next_level, box_max, max_grid_coord);
1313
for (i = 0; i < 2; i++) {
1314
if (_node_intersect_box_1d(children[i], min_grid_coord, max_grid_coord))
1315
(next_bt->nodes[n_init_nodes + i]).n_boxes += 1;
1319
else { /* Box is included in the same segment */
1321
assert( fvm_morton_compare(1, max_code, min_code)
1322
== FVM_MORTON_EQUAL_ID);
1324
for (i = 0; i < 2; i++) {
1325
if ( fvm_morton_compare(1, min_code, children[i])
1326
== FVM_MORTON_EQUAL_ID) {
1327
(next_bt->nodes[n_init_nodes + i]).n_boxes += 1;
1330
} /* End of loop on children */
1332
} /* If min_code and max_code in the same leaf */
1334
} /* End of loop on boxes */
1338
for (i = 0; i < 2; i++) {
1340
(next_bt->nodes[n_init_nodes + i]).start_id
1341
= _shift_ids + n_linked_boxes;
1342
n_linked_boxes += (next_bt->nodes[n_init_nodes + i]).n_boxes;
1346
_shift_ids += n_linked_boxes;
1348
for (i = 0; i < 2; i++)
1349
(next_bt->nodes[n_init_nodes + i]).n_boxes = 0;
1351
/* Second loop on boxes associated to the node_id: fill */
1353
for (j = 0; j < node.n_boxes; j++) {
1355
fvm_coord_t min_grid_coord[2], max_grid_coord[2];
1357
fvm_lnum_t box_id = bt->box_ids[node.start_id + j];
1358
const fvm_coord_t *box_min = _box_min(boxes, box_id);
1359
const fvm_coord_t *box_max = _box_max(boxes, box_id);
1361
min_code = fvm_morton_encode(1, next_level, box_min);
1362
max_code = fvm_morton_encode(1, next_level, box_max);
1364
if ( fvm_morton_compare(1, min_code, max_code)
1365
== FVM_MORTON_DIFFERENT_ID) {
1367
_get_grid_coords_1d(next_level, box_min, min_grid_coord);
1368
_get_grid_coords_1d(next_level, box_max, max_grid_coord);
1370
for (i = 0; i < 2; i++) {
1372
if (_node_intersect_box_1d(children[i],
1376
const fvm_lnum_t sub_id = n_init_nodes + i;
1377
const fvm_lnum_t shift = (next_bt->nodes[sub_id]).n_boxes
1378
+ (next_bt->nodes[sub_id]).start_id;
1379
next_bt->box_ids[shift] = box_id;
1380
(next_bt->nodes[sub_id]).n_boxes += 1;
1383
} /* End of loop on children*/
1386
else { /* Box is included in the same segment */
1388
assert( fvm_morton_compare(1, max_code, min_code)
1389
== FVM_MORTON_EQUAL_ID);
1391
for (i = 0; i < 2; i++) {
1393
if ( fvm_morton_compare(1, min_code, children[i])
1394
== FVM_MORTON_EQUAL_ID) {
1396
const fvm_lnum_t sub_id = n_init_nodes + i;
1397
const fvm_lnum_t shift = (next_bt->nodes[sub_id]).n_boxes
1398
+ (next_bt->nodes[sub_id]).start_id;
1400
next_bt->box_ids[shift] = box_id;
1401
(next_bt->nodes[sub_id]).n_boxes += 1;
1405
} /* End of loop on children */
1407
} /* If min_code and max_code in the same leaf */
1409
} /* End of loop on boxes */
1411
#if 0 && defined(DEBUG) && !defined(NDEBUG)
1412
for (i = 1; i < 2; i++) {
1413
_node_t n1 = next_bt->nodes[n_init_nodes + i - 1];
1414
_node_t n2 = next_bt->nodes[n_init_nodes + i];
1415
assert(n1.n_boxes == (n2.start_id - n1.start_id));
1418
== ( next_bt->nodes[n_init_nodes + 2 - 1].start_id
1419
+ next_bt->nodes[n_init_nodes + 2 - 1].n_boxes));
1422
/* Return pointers */
1424
*shift_ids = _shift_ids;
1427
/*----------------------------------------------------------------------------
1428
* Evaluate the box distribution over the leaves of the box tree when adding
1429
* a level to the tree structure.
1432
* bt <-> pointer to the box tree being built
1433
* next_bt <-> pointer to the next box tree being built
1434
* boxes <-- pointer to the associated box set structure
1435
* node_id <-- id of the starting node
1436
* build_type <-- layout variant for building the tree structure
1437
* shift_ids <-> first free position free in new box_ids
1438
*----------------------------------------------------------------------------*/
1441
_build_next_level(fvm_box_tree_t *bt,
1442
fvm_box_tree_t *next_bt,
1443
const fvm_box_set_t *boxes,
1445
fvm_box_tree_sync_t build_type,
1446
fvm_lnum_t *shift_ids)
1450
fvm_lnum_t _shift_ids = *shift_ids;
1451
const _node_t *cur_node = bt->nodes + node_id;
1453
if (cur_node->is_leaf == false) {
1455
assert(bt->child_ids[bt->n_children*node_id] > 0);
1457
for (i = 0; i < bt->n_children; i++)
1458
_build_next_level(bt,
1461
bt->child_ids[bt->n_children*node_id + i],
1466
else { /* if (node->is_leaf == true) */
1468
if ( cur_node->n_boxes < bt->threshold
1469
&& node_id != 0 /* Root node is always divided */
1470
&& build_type == FVM_BOX_TREE_ASYNC_LEVEL) {
1472
/* Copy related box_ids in the new next_ids */
1474
_node_t *next_node = next_bt->nodes + node_id;
1476
next_node->n_boxes = cur_node->n_boxes;
1477
next_node->start_id = _shift_ids;
1479
for (i = 0; i < cur_node->n_boxes; i++)
1480
next_bt->box_ids[_shift_ids++]
1481
= bt->box_ids[cur_node->start_id + i];
1483
else { /* Split node and evaluate box distribution between its children */
1485
if (boxes->dim == 3)
1491
else if (boxes->dim == 2)
1497
else if (boxes->dim == 1)
1507
/* Prepare return values */
1509
*shift_ids = _shift_ids;
1512
#if defined(HAVE_MPI)
1514
/*----------------------------------------------------------------------------
1515
* Loop on all nodes of the box tree to define an array with Morton codes
1516
* and weights (= number of linked boxes) associated to each leaf
1519
* bt <-> pointer to fvm_box_tree_t structure.
1520
* boxes <-- pointer to the associated box set structure
1521
* node_id <-- id of the current node (to traverse)
1522
* n_leaves <-> current number of leaves in the tree with n_boxes > 0
1523
* leaf_codes <-> Morton code associated to each leaf
1524
* weight <-> number of boxes attached to each leaf
1525
*----------------------------------------------------------------------------*/
1528
_build_leaf_weight(const fvm_box_tree_t *bt,
1530
fvm_lnum_t *n_leaves,
1531
fvm_morton_code_t *leaf_codes,
1536
fvm_lnum_t _n_leaves = *n_leaves;
1538
const _node_t *node = bt->nodes + node_id;
1540
if (node->is_leaf == false)
1541
for (i = 0; i < bt->n_children; i++)
1542
_build_leaf_weight(bt,
1543
bt->child_ids[bt->n_children*node_id + i],
1548
else { /* node is a leaf */
1550
if (node->n_boxes > 0) {
1551
leaf_codes[_n_leaves] = node->morton_code;
1552
weight[_n_leaves] = node->n_boxes;
1557
*n_leaves = _n_leaves;
1560
/*----------------------------------------------------------------------------
1561
* Loop on all nodes of the tree to define an index of ranks related
1565
* bt <-> pointer to fvm_box_tree_t structure.
1566
* distrib <-- structure holding box distribution data
1567
* dim <-- box tree layout dimension (1, 2, or 3)
1568
* node_id <-- id of the current node (to traverse)
1569
* size <-- size of index in which we search
1570
* search_index <-- index on which box distribution is made
1571
* id_rank <-- relation between id and rank
1572
*----------------------------------------------------------------------------*/
1575
_build_rank_to_box_index(const fvm_box_tree_t *bt,
1576
fvm_box_distrib_t *distrib,
1580
fvm_morton_code_t search_index[],
1585
const _node_t *node = bt->nodes + node_id;
1587
if (node->is_leaf == false) {
1589
for (i = 0; i < bt->n_children; i++)
1590
_build_rank_to_box_index(bt,
1593
bt->child_ids[bt->n_children*node_id + i],
1600
if (node->n_boxes > 0) {
1602
int id = fvm_morton_binary_search(size,
1605
int rank = id_rank[id];
1607
distrib->index[rank + 1] += node->n_boxes;
1614
/*----------------------------------------------------------------------------
1615
* Loop on all nodes of the tree to define a list of ranks related
1619
* bt <-> pointer to fvm_box_tree_t structure.
1620
* distrib <-- structure holding box distribution data
1621
* dim <-- box tree layout dimension (1, 2, or 3)
1622
* node_id <-- id of the current node (to traverse)
1623
* counter <-> counter array used to build the list
1624
* size <-- size of index in which we search
1625
* search_index <-- index on which box distribution is made
1626
* id_rank <-- relation between id and rank
1627
*----------------------------------------------------------------------------*/
1630
_build_rank_to_box_list(const fvm_box_tree_t *bt,
1631
fvm_box_distrib_t *distrib,
1634
fvm_lnum_t counter[],
1636
fvm_morton_code_t search_index[],
1641
const _node_t *node = bt->nodes + node_id;
1643
if (node->is_leaf == false) {
1645
for (i = 0; i < bt->n_children; i++)
1646
_build_rank_to_box_list(bt,
1649
bt->child_ids[bt->n_children*node_id + i],
1657
if (node->n_boxes > 0) {
1659
int id = fvm_morton_binary_search(size,
1662
int rank = id_rank[id];
1664
for (i = 0; i < node->n_boxes; i++) {
1666
fvm_lnum_t box_id = bt->box_ids[node->start_id + i];
1667
fvm_lnum_t shift = distrib->index[rank] + counter[rank];
1669
distrib->list[shift] = box_id;
1678
#endif /* defined(HAVE_MPI) */
1680
/*----------------------------------------------------------------------------
1681
* Recursively build an index on boxes which intersect.
1684
* bt <-> pointer to fvm_box_tree_t structure.
1685
* boxes <-- pointer to associated boxes structure
1686
* node_id <-- id of the current node (to traverse)
1687
* count <-> intersection count
1688
*----------------------------------------------------------------------------*/
1691
_count_intersections(const fvm_box_tree_t *bt,
1692
const fvm_box_set_t *boxes,
1698
const fvm_coord_t *box_extents = boxes->extents;
1699
const _node_t *node = bt->nodes + node_id;
1701
if (node->is_leaf == false) {
1703
for (i = 0; i < bt->n_children; i++) /* traverse downwards */
1704
_count_intersections(bt,
1706
bt->child_ids[bt->n_children*node_id + i],
1709
else { /* node is a leaf */
1711
if (boxes->dim == 3) {
1713
for (i = 0; i < node->n_boxes - 1; i++) {
1714
for (j = i+1; j < node->n_boxes; j++) {
1715
fvm_lnum_t id0 = bt->box_ids[node->start_id + i];
1716
fvm_lnum_t id1 = bt->box_ids[node->start_id + j];
1717
if (_boxes_intersect_3d(box_extents, id0, id1)) {
1726
else if (boxes->dim == 2) {
1728
for (i = 0; i < node->n_boxes - 1; i++) {
1729
for (j = i+1; j < node->n_boxes; j++) {
1730
fvm_lnum_t id0 = bt->box_ids[node->start_id + i];
1731
fvm_lnum_t id1 = bt->box_ids[node->start_id + j];
1732
if (_boxes_intersect_2d(box_extents, id0, id1)) {
1741
else if (boxes->dim == 1) {
1743
for (i = 0; i < node->n_boxes - 1; i++) {
1744
for (j = i+1; j < node->n_boxes; j++) {
1745
fvm_lnum_t id0 = bt->box_ids[node->start_id + i];
1746
fvm_lnum_t id1 = bt->box_ids[node->start_id + j];
1747
if (_boxes_intersect_1d(box_extents, id0, id1)) {
1760
/*----------------------------------------------------------------------------
1761
* Recursively build a list on bounding boxes which intersect together.
1764
* bt <-> pointer to fvm_box_tree_t structure.
1765
* boxes <-- pointer to associated boxes structure
1766
* node_id <-- id of the current node (to traverse)
1767
* count <-> intersection count (workin array)
1768
* index <-- index on intersections
1769
* box_g_num <-> global number of intersection boxes
1770
*----------------------------------------------------------------------------*/
1773
_get_intersections(const fvm_box_tree_t *bt,
1774
const fvm_box_set_t *boxes,
1777
fvm_lnum_t box_index[],
1778
fvm_gnum_t box_g_num[])
1782
const fvm_coord_t *box_extents = boxes->extents;
1783
const _node_t *node = bt->nodes + node_id;
1785
if (node->is_leaf == false) {
1787
for (i = 0; i < bt->n_children; i++) /* traverse downwards */
1788
_get_intersections(bt,
1790
bt->child_ids[bt->n_children*node_id + i],
1795
else { /* node is a leaf */
1797
if (boxes->dim == 3) {
1799
for (i = 0; i < node->n_boxes - 1; i++) {
1800
for (j = i+1; j < node->n_boxes; j++) {
1801
fvm_lnum_t id0 = bt->box_ids[node->start_id + i];
1802
fvm_lnum_t id1 = bt->box_ids[node->start_id + j];
1804
if (_boxes_intersect_3d(box_extents, id0, id1)) {
1805
fvm_lnum_t shift0 = box_index[id0] + count[id0];
1806
fvm_lnum_t shift1 = box_index[id1] + count[id1];
1807
box_g_num[shift0] = boxes->g_num[id1];
1808
box_g_num[shift1] = boxes->g_num[id0];
1815
else if (boxes->dim == 2) {
1817
for (i = 0; i < node->n_boxes - 1; i++) {
1818
for (j = i+1; j < node->n_boxes; j++) {
1819
fvm_lnum_t id0 = bt->box_ids[node->start_id + i];
1820
fvm_lnum_t id1 = bt->box_ids[node->start_id + j];
1822
if (_boxes_intersect_2d(box_extents, id0, id1)) {
1823
fvm_lnum_t shift0 = box_index[id0] + count[id0];
1824
fvm_lnum_t shift1 = box_index[id1] + count[id1];
1825
box_g_num[shift0] = boxes->g_num[id1];
1826
box_g_num[shift1] = boxes->g_num[id0];
1834
else if (boxes->dim == 1) {
1836
for (i = 0; i < node->n_boxes - 1; i++) {
1837
for (j = i+1; j < node->n_boxes; j++) {
1838
fvm_lnum_t id0 = bt->box_ids[node->start_id + i];
1839
fvm_lnum_t id1 = bt->box_ids[node->start_id + j];
1841
if (_boxes_intersect_1d(box_extents, id0, id1)) {
1842
fvm_lnum_t shift0 = box_index[id0] + count[id0];
1843
fvm_lnum_t shift1 = box_index[id1] + count[id1];
1844
box_g_num[shift0] = boxes->g_num[id1];
1845
box_g_num[shift1] = boxes->g_num[id0];
1853
} /* End if node is a leaf */
1856
/*----------------------------------------------------------------------------
1857
* Recursively define a counter array on the number of bounding boxes
1858
* associated to a leaf.
1860
* This will be used for displaying a histogram.
1863
* bt <-- pointer to fvm_box_tree_t structure.
1864
* node_id <-- id of the current node (to traverse)
1865
* n_steps <-- number of steps in histogram
1866
* step <-- steps of the histogram
1867
* h_min <-- min. value of the histogram
1868
* counter <-> counter (working array)
1869
*----------------------------------------------------------------------------*/
1872
_build_histogram(const fvm_box_tree_t *bt,
1881
const _node_t *node = bt->nodes + node_id;
1883
if (node->is_leaf == false) {
1885
for (i = 0; i < bt->n_children; i++) /* traverse downwards */
1886
_build_histogram(bt,
1887
bt->child_ids[bt->n_children*node_id + i],
1894
for (i = 0, j = 1; j < n_steps; i++, j++)
1895
if (node->n_boxes < h_min + j*step)
1901
/*----------------------------------------------------------------------------
1902
* Dump a box tree node.
1905
* bt <-- pointer to fvm_box_tree_t structure.
1906
* node_id <-- id of the current node (to traverse)
1907
*----------------------------------------------------------------------------*/
1910
_dump_node(const fvm_box_tree_t *bt,
1915
const char *node_type[] = {"node", "leaf"};
1917
const _node_t *node = bt->nodes + node_id;
1918
const fvm_morton_code_t m_code = node->morton_code;
1922
" level: %3hu - anchor: [ %10u %10u %10u ]\n"
1923
" n_boxes: %3hd - start_id: %u\n"
1925
node_id, node_type[(int)(node->is_leaf)],
1926
m_code.L, m_code.X[0], m_code.X[1], m_code.X[2],
1927
node->n_boxes, node->start_id);
1929
for (i = 0; i < node->n_boxes; i++)
1930
bft_printf(" %d\n", (int)(bt->box_ids[node->start_id + i]));
1932
if (node->is_leaf == false) {
1934
const fvm_lnum_t *c_id = bt->child_ids + bt->n_children*node_id;
1936
if (bt->n_children == 8)
1937
bft_printf(" children_id: %d %d %d %d %d %d %d %d\n",
1938
(int)c_id[0], (int)c_id[1], (int)c_id[2], (int)c_id[3],
1939
(int)c_id[4], (int)c_id[5], (int)c_id[6], (int)c_id[7]);
1940
else if (bt->n_children == 4)
1941
bft_printf(" children_id: %d %d %d %d\n",
1942
(int)c_id[0], (int)c_id[1], (int)c_id[2], (int)c_id[3]);
1943
else if (bt->n_children == 2)
1944
bft_printf(" children_id: %d %d\n",
1945
(int)c_id[0], (int)c_id[1]);
1947
for (i = 0; i < bt->n_children; i++)
1948
_dump_node(bt, c_id[i]);
1952
/*============================================================================
1953
* Public function definitions
1954
*============================================================================*/
1956
/*----------------------------------------------------------------------------
1957
* Create a fvm_box_tree_t structure and initialize it.
1960
* max_level <-- max possible level
1961
* threshold <-- max number of boxes linked to an octant if
1962
* max_level is not reached
1963
* max_box_ratio <-- max n_linked_boxes / n_boxes ratio
1966
* pointer to an empty fvm_box_tree_t structure.
1967
*----------------------------------------------------------------------------*/
1970
fvm_box_tree_create(int max_level,
1972
float max_box_ratio)
1974
fvm_box_tree_t *bt = NULL;
1976
BFT_MALLOC(bt, 1, fvm_box_tree_t);
1981
bft_error(__FILE__, __LINE__, 0,
1982
_(" Forbidden max_level value (%d) in the tree structure\n"),
1986
bft_error(__FILE__, __LINE__, 0,
1987
_(" Forbidden threshold value (%d) in the tree structure\n"),
1990
if (max_box_ratio < 1.0)
1991
bft_error(__FILE__, __LINE__, 0,
1992
_(" Forbidden max_box_ratio value (%f) in the tree structure\n"),
1993
(double)max_box_ratio);
1995
/* Create and initialize tree structure according to its type */
1997
bt->max_level = max_level;
1998
bt->threshold = threshold;
1999
bt->max_box_ratio = max_box_ratio;
2001
#if defined(HAVE_MPI)
2002
bt->comm = MPI_COMM_NULL;
2007
bt->stats.max_level_reached = 0;
2009
bt->stats.n_leaves = 0;
2010
bt->stats.n_spill_leaves = 0;
2011
bt->stats.n_linked_boxes = 0;
2013
bt->stats.min_linked_boxes = INT_MAX;
2014
bt->stats.max_linked_boxes = 0;
2016
/* Initialize nodes */
2018
bt->n_max_nodes = 0;
2025
bt->n_build_loops = 0;
2030
/*----------------------------------------------------------------------------
2031
* Destroy a fvm_box_tree_t structure.
2034
* bt <-- pointer to pointer to fvm_box_tree_t structure to destroy
2035
*----------------------------------------------------------------------------*/
2038
fvm_box_tree_destroy(fvm_box_tree_t **bt)
2040
fvm_box_tree_t *_bt = *bt;
2044
BFT_FREE(_bt->nodes);
2045
BFT_FREE(_bt->child_ids);
2046
BFT_FREE(_bt->box_ids);
2053
/*----------------------------------------------------------------------------
2054
* Get the deepest level allowed by the tree structure.
2057
* bt <-- pointer to fvm_box_tree_t structure.
2060
* deepest allowed level of the tree
2061
*----------------------------------------------------------------------------*/
2064
fvm_box_tree_get_max_level(const fvm_box_tree_t *bt)
2066
return bt->max_level;
2069
/*----------------------------------------------------------------------------
2070
* Assign a set of boxes to an empty fvm_box_tree_t structure.
2072
* The box tree structure must have been created using to fvm_tree_create().
2074
* The depth of the tree is adjusted so that a maximum of max_n_elts boxes
2075
* will be assigned to each leaf, unless this would require going beyond
2076
* the tree's maximum level.
2078
* If max_level = -1, the highest level reachable is FVM_TREE_MAX_LEVEL but
2079
* there is no defined target level.
2082
* bt <-> pointer to fvm_box_tree_t structure.
2083
* boxes <-- pointer to the associated box set structure
2084
* build_type <-- layout variant for building the tree structure
2085
*----------------------------------------------------------------------------*/
2088
fvm_box_tree_set_boxes(fvm_box_tree_t *bt,
2089
const fvm_box_set_t *boxes,
2090
fvm_box_tree_sync_t build_type)
2094
fvm_box_tree_t tmp_bt;
2096
fvm_lnum_t next_box_ids_size = 0, shift = 0;
2097
fvm_coord_t anchor[3] = {0., 0., 0.};
2099
/* Initialization */
2103
bt->n_build_loops = 0;
2105
#if defined(HAVE_MPI)
2106
bt->comm = boxes->comm;
2109
/* Preallocate for the two first levels of a tree */
2111
if (boxes->dim == 3) {
2113
bt->n_max_nodes = 73;
2115
else if (boxes->dim == 2) {
2117
bt->n_max_nodes = 21;
2119
else if (boxes->dim == 1) {
2121
bt->n_max_nodes = 7;
2126
BFT_MALLOC(bt->nodes, bt->n_max_nodes, _node_t);
2127
BFT_MALLOC(bt->child_ids,
2128
bt->n_max_nodes*bt->n_children,
2131
/* Define root node */
2133
_new_node(bt, fvm_morton_encode(boxes->dim, 0, anchor), 0);
2135
/* Initialize bt by assigning all boxes to the root leaf */
2137
BFT_MALLOC(bt->box_ids, boxes->n_boxes, fvm_lnum_t);
2139
for (box_id = 0; box_id < boxes->n_boxes; box_id++)
2140
bt->box_ids[box_id] = box_id;
2142
(bt->nodes[0]).is_leaf = true;
2143
(bt->nodes[0]).n_boxes = boxes->n_boxes;
2144
(bt->nodes[0]).start_id = 0;
2146
bt->stats.n_boxes = boxes->n_boxes;
2148
_get_box_tree_stats(bt);
2150
/* Build local tree structure by adding boxes from the root */
2152
while (_recurse_tree_build(bt,
2155
&next_box_ids_size)) {
2157
/* Initialize next_bt: copy of bt */
2159
_copy_tree(&tmp_bt, bt);
2161
/* Optimize memory usage */
2163
bt->n_max_nodes = bt->n_nodes;
2164
BFT_REALLOC(bt->nodes, bt->n_nodes, _node_t);
2165
BFT_REALLOC(bt->child_ids,
2166
bt->n_max_nodes*bt->n_children,
2169
/* Define a box ids list for the next level of the boxtree */
2171
BFT_REALLOC(tmp_bt.box_ids, next_box_ids_size, fvm_lnum_t);
2174
_build_next_level(bt,
2177
0, /* Starts from root */
2181
assert(shift == next_box_ids_size);
2183
/* replace current tree by the tree computed at a higher level */
2185
_free_tree_arrays(bt);
2186
*bt = tmp_bt; /* Overwrite bt members with those of next_bt */
2188
_get_box_tree_stats(bt);
2190
#if 0 && defined(DEBUG) && !defined(NDEBUG)
2191
bft_printf(" - New box tree level -\n");
2192
fvm_box_tree_dump_statistics(bt);
2195
} /* While building should continue */
2198
#if defined(HAVE_MPI)
2200
/*----------------------------------------------------------------------------
2201
* Compute an index based on Morton encoding to ensure a good distribution
2202
* of boxes among the participating ranks.
2205
* bt <-> pointer to fvm_box_tree_t structure.
2206
* boxes <-- pointer to the associated box set structure
2209
* pointer to newly created fvm_box_distrib_t structure.
2210
*----------------------------------------------------------------------------*/
2213
fvm_box_tree_get_distrib(fvm_box_tree_t *bt,
2214
const fvm_box_set_t *boxes)
2218
int reduce_size = 0;
2219
fvm_lnum_t n_leaves = 0;
2220
int *reduce_ids = NULL;
2221
fvm_morton_code_t *leaf_codes = NULL, *reduce_index = NULL;
2222
fvm_lnum_t *weight = NULL, *counter = NULL;
2224
fvm_box_distrib_t *distrib = NULL;
2227
assert(boxes != NULL);
2229
/* Compute basic box distribution */
2231
distrib = fvm_box_distrib_create(boxes->n_boxes,
2233
(bt->stats).max_level_reached,
2236
if (distrib == NULL)
2239
BFT_MALLOC(leaf_codes, bt->stats.n_leaves, fvm_morton_code_t);
2240
BFT_MALLOC(weight, bt->stats.n_leaves, fvm_lnum_t);
2242
/* Build index for boxes */
2244
_build_leaf_weight(bt,
2250
assert(n_leaves <= bt->stats.n_leaves);
2252
BFT_REALLOC(leaf_codes, n_leaves, fvm_morton_code_t);
2253
BFT_REALLOC(weight, n_leaves, fvm_lnum_t);
2255
/* Compute the resulting Morton index */
2257
fvm_box_set_build_morton_index(boxes,
2263
BFT_FREE(leaf_codes);
2266
/* Compact Morton_index to get an array without "0 element" */
2268
for (i = 0; i < distrib->n_ranks; i++)
2269
if (fvm_morton_a_gt_b(distrib->morton_index[i+1],
2270
distrib->morton_index[i]))
2273
BFT_MALLOC(reduce_index, reduce_size + 1, fvm_morton_code_t);
2274
BFT_MALLOC(reduce_ids, reduce_size, int);
2277
reduce_index[0] = distrib->morton_index[0];
2279
for (i = 0; i < distrib->n_ranks; i++) {
2281
if (fvm_morton_a_gt_b(distrib->morton_index[i+1],
2282
distrib->morton_index[i])) {
2284
reduce_index[reduce_size + 1] = distrib->morton_index[i+1];
2285
reduce_ids[reduce_size++] = i;
2291
/* Define a rank -> box indexed list */
2293
_build_rank_to_box_index(bt,
2296
0, /* starts from root */
2301
for (i = 0; i < distrib->n_ranks; i++)
2302
distrib->index[i+1] += distrib->index[i];
2304
BFT_MALLOC(distrib->list,
2305
distrib->index[distrib->n_ranks], int);
2307
BFT_MALLOC(counter, distrib->n_ranks, fvm_lnum_t);
2309
for (i = 0; i < distrib->n_ranks; i++)
2312
_build_rank_to_box_list(bt,
2315
0, /* starts from root */
2324
BFT_FREE(reduce_ids);
2325
BFT_FREE(reduce_index);
2327
/* Define the final index (without redundancies) and realloc list */
2329
fvm_box_distrib_clean(distrib);
2334
#endif /* defined(HAVE_MPI) */
2336
/*----------------------------------------------------------------------------
2337
* Build an indexed list on boxes to list intersections.
2339
* The index and box_g_num arrays are allocated by this function,
2340
* and it is the caller's responsibility to free them.
2342
* Upon return, box_index[i] points to the first position in box_g_num
2343
* relative to boxes intersecting box i of the boxes set, while
2344
* box_g_num contains the global numbers associated with those boxes.
2347
* bt <-- pointer to box tree structure to query
2348
* boxes <-- pointer to a associated box set
2349
* box_index --> pointer to the index array on bounding boxes
2350
* box_g_num --> pointer to the list of intersecting bounding boxes
2351
*----------------------------------------------------------------------------*/
2354
fvm_box_tree_get_intersects(fvm_box_tree_t *bt,
2355
const fvm_box_set_t *boxes,
2356
fvm_lnum_t *box_index[],
2357
fvm_gnum_t *box_g_num[])
2359
fvm_lnum_t i, list_size;
2361
fvm_lnum_t *counter = NULL;
2362
fvm_lnum_t *_index = NULL;
2363
fvm_gnum_t *_g_num = NULL;
2367
BFT_MALLOC(_index, boxes->n_boxes + 1, fvm_lnum_t);
2369
for (i = 0; i < boxes->n_boxes + 1; i++)
2372
_count_intersections(bt,
2374
0, /* start from root */
2377
/* Build index from counts */
2379
for (i = 0; i < boxes->n_boxes; i++)
2380
_index[i+1] += _index[i];
2382
list_size = _index[boxes->n_boxes];
2384
BFT_MALLOC(_g_num, list_size, fvm_gnum_t);
2386
BFT_MALLOC(counter, boxes->n_boxes, fvm_lnum_t);
2388
for (i = 0; i < boxes->n_boxes; i++)
2393
_get_intersections(bt,
2395
0, /* start from root */
2402
/* Return pointers */
2404
*box_index = _index;
2405
*box_g_num = _g_num;
2408
/*----------------------------------------------------------------------------
2409
* Get global box tree statistics.
2411
* All fields returned are optional: if their argument is set to NULL,
2412
* the corresponding information will not be returned.
2414
* For each field not set to NULL, 3 values are always returned:
2415
* the mean on all ranks (rounded to the closest integer), the minimum,
2416
* and the maximum value respectively.
2418
* In serial mode, the mean, minimum, and maximum will be identical for most
2419
* fields, but all 3 values are returned nonetheless.
2421
* Note that the theoretical memory use includes that of the associated
2425
* bt <-- pointer to box tree structure
2426
* depth --> tree depth (max level used)
2427
* n_leaves --> number of leaves in the tree
2428
* n_boxes --> number of boxes in the tree
2429
* n_threshold_leaves --> number of leaves where n_boxes > threshold
2430
* n_leaf_boxes --> number of boxes for a leaf
2431
* mem_used --> theoretical used memory
2432
* mem_allocated --> theoretical allocated memory
2435
* the spatial dimension associated with the box tree layout (3, 2, or 1)
2436
*----------------------------------------------------------------------------*/
2439
fvm_box_tree_get_stats(const fvm_box_tree_t *bt,
2441
fvm_lnum_t n_leaves[3],
2442
fvm_lnum_t n_boxes[3],
2443
fvm_lnum_t n_threshold_leaves[3],
2444
fvm_lnum_t n_leaf_boxes[3],
2446
size_t mem_allocated[3])
2449
size_t mem_per_node;
2450
size_t s_mean[7], s_min[7], s_max[7];
2451
fvm_box_tree_stats_t s;
2460
if (bt->n_children == 4)
2462
else if (bt->n_children == 2)
2465
/* Prepare array of local values; prior to or in the absence of
2466
MPI communication, mean values are set to local values. */
2468
s_mean[0] = s.n_linked_boxes / s.n_leaves;
2469
/* Round to nearest integer, and not floor */
2470
if (s.n_linked_boxes % s.n_leaves >= s.n_leaves/2)
2473
s_min[0] = s.min_linked_boxes;
2474
s_max[0] = s.max_linked_boxes;
2476
s_mean[1] = s.max_level_reached;
2477
s_mean[2] = s.n_leaves;
2478
s_mean[3] = s.n_boxes;
2479
s_mean[4] = s.n_spill_leaves;
2481
/* Estimate theoretical memory usage */
2483
mem_per_node = sizeof(_node_t) + bt->n_children*sizeof(fvm_lnum_t);
2485
s_mean[5] = sizeof(fvm_box_tree_t);
2486
s_mean[5] += bt->n_nodes * mem_per_node;
2487
s_mean[5] += s.n_linked_boxes * sizeof(fvm_lnum_t);
2489
s_mean[5] += sizeof(fvm_box_set_t);
2490
s_mean[5] += s.n_boxes * ( sizeof(fvm_gnum_t)
2491
+ (dim * 2 * sizeof(fvm_coord_t)));
2493
s_mean[6] = s_mean[5] + (bt->n_max_nodes - bt->n_nodes)*mem_per_node;
2495
/* Pre-synchronize for serial cases (i = 0 already handled) */
2497
for (i = 1; i < 7; i++) {
2498
s_min[i] = s_mean[i];
2499
s_max[i] = s_mean[i];
2502
/* In parallel mode, synchronize values */
2504
#if defined(HAVE_MPI)
2506
if (bt->comm != MPI_COMM_NULL) {
2509
fvm_gnum_t s_l_sum[14], s_g_sum[14];
2511
MPI_Comm_size(bt->comm, &n_ranks);
2513
if (n_ranks > 1) { /* Should always be the case, bat play it safe) */
2515
/* Split value to avoid exceeding fvm_gnum_t limits
2516
(especially when it is a 32-bit value) */
2518
s_l_sum[0] = s.n_linked_boxes/n_ranks;
2519
s_l_sum[7] = s.n_linked_boxes%n_ranks;
2520
for (i = 1; i < 7; i++) {
2521
s_l_sum[i] = s_mean[i]/n_ranks;
2522
s_l_sum[i+7] = s_mean[i]%n_ranks;
2525
MPI_Allreduce(s_l_sum, s_g_sum, 14, FVM_MPI_GNUM, MPI_SUM, bt->comm);
2527
s_mean[0] = s.min_linked_boxes;
2528
MPI_Allreduce(s_mean, s_min, 7, FVM_MPI_GNUM, MPI_MIN, bt->comm);
2529
s_mean[0] = s.max_linked_boxes;
2530
MPI_Allreduce(s_mean, s_max, 7, FVM_MPI_GNUM, MPI_MAX, bt->comm);
2532
/* Specific handling for linked boxes, so as to ensure correct
2533
total using large integers even if we do not know the
2534
corresponding MPI type */
2536
uint64_t s_n = s_g_sum[0]*n_ranks + s_g_sum[7]; /* linked boxes */
2537
uint64_t s_d = s_g_sum[2]*n_ranks + s_g_sum[9]; /* leaves */
2538
uint64_t s_m = s_n / s_d;
2539
/* Round to nearest integer, and not floor */
2540
if (s_n % s_d >= s_d/2)
2544
for (i = 1; i < 7; i++) {
2545
s_mean[i] = s_g_sum[i] + s_g_sum[i+7]/n_ranks;
2546
/* Round to nearest integer, and not floor */
2547
if (s_g_sum[i+7]%n_ranks >= (fvm_gnum_t)n_ranks/2)
2555
/* Set values already in stats */
2557
if (depth != NULL) {
2558
depth[0] = s_mean[1];
2559
depth[1] = s_min[1];
2560
depth[2] = s_max[1];
2563
if (n_leaves != NULL) {
2564
n_leaves[0] = s_mean[2];
2565
n_leaves[1] = s_min[2];
2566
n_leaves[2] = s_max[2];
2569
if (n_boxes != NULL) {
2570
n_boxes[0] = s_mean[3];
2571
n_boxes[1] = s_min[3];
2572
n_boxes[2] = s_max[3];
2575
if (n_threshold_leaves != NULL) {
2576
n_threshold_leaves[0] = s_mean[4];
2577
n_threshold_leaves[1] = s_min[4];
2578
n_threshold_leaves[2] = s_max[4];
2581
if (n_leaf_boxes != NULL) {
2582
n_leaf_boxes[0] = s_mean[0];
2583
n_leaf_boxes[1] = s_min[0];
2584
n_leaf_boxes[2] = s_max[0];
2587
if (mem_used != NULL) {
2588
mem_used[0] = s_mean[5];
2589
mem_used[1] = s_min[5];
2590
mem_used[2] = s_max[5];
2593
if (mem_allocated != NULL) {
2594
mem_allocated[0] = s_mean[6];
2595
mem_allocated[1] = s_min[6];
2596
mem_allocated[2] = s_max[6];
2602
/*----------------------------------------------------------------------------
2603
* Display local statistics about a fvm_box_tree_t structure.
2606
* bt <-- pointer to box tree structure
2607
*----------------------------------------------------------------------------*/
2610
fvm_box_tree_dump_statistics(const fvm_box_tree_t *bt)
2613
fvm_box_tree_stats_t s;
2614
unsigned g_max_level_reached;
2615
fvm_gnum_t n_g_leaves, n_g_boxes, n_g_linked_boxes, n_g_spill_leaves;
2616
fvm_lnum_t g_min_linked_boxes, g_max_linked_boxes;
2617
double mean_linked_boxes, box_ratio;
2619
fvm_gnum_t count[5];
2621
int step = 0, delta = 0;
2622
const int n_steps = 5;
2629
g_max_level_reached = s.max_level_reached;
2630
n_g_leaves = s.n_leaves;
2631
n_g_boxes = s.n_boxes;
2632
n_g_linked_boxes = s.n_linked_boxes;
2633
n_g_spill_leaves = s.n_spill_leaves;
2634
g_min_linked_boxes = s.min_linked_boxes;
2635
g_max_linked_boxes = s.max_linked_boxes;
2637
#if defined(HAVE_MPI)
2639
if (bt->comm != MPI_COMM_NULL) {
2641
fvm_gnum_t l_min[1], g_min[1];
2642
fvm_gnum_t l_max[2], g_max[2];
2643
fvm_gnum_t l_sum[3], g_sum[3];
2645
l_sum[0] = n_g_leaves;
2646
l_sum[1] = n_g_spill_leaves;
2647
l_sum[2] = n_g_linked_boxes;
2649
l_min[0] = g_min_linked_boxes;
2650
l_max[0] = s.max_level_reached;
2651
l_max[1] = g_max_linked_boxes;
2653
MPI_Allreduce(l_sum, g_sum, 3, FVM_MPI_GNUM, MPI_SUM, bt->comm);
2654
MPI_Allreduce(l_min, g_min, 1, FVM_MPI_GNUM, MPI_MIN, bt->comm);
2655
MPI_Allreduce(l_max, g_max, 2, FVM_MPI_GNUM, MPI_MAX, bt->comm);
2657
n_g_leaves = l_sum[0];
2658
n_g_spill_leaves = l_sum[1];
2659
n_g_linked_boxes = l_sum[2];
2661
g_min_linked_boxes = g_min[0];
2662
g_max_level_reached = g_max[0];
2663
g_max_linked_boxes = g_max[1];
2668
/* Redefine final statistics */
2670
mean_linked_boxes = (double)n_g_linked_boxes / (double)n_g_leaves;
2671
box_ratio = (double)n_g_linked_boxes / (double)n_g_boxes;
2673
/* Define axis subdivisions */
2675
for (j = 0; j < n_steps; j++)
2678
delta = g_max_linked_boxes - g_min_linked_boxes;
2682
step = delta/n_steps;
2684
_build_histogram(bt,
2685
0, /* start from root */
2691
} /* max - min > 0 */
2693
/* Print statistics and bounding boxes histogram */
2696
"Box tree statistics:\n\n");
2697
bft_printf(" Number of children per leaf: %d\n"
2698
" Max number of bounding boxes for a leaf: %d\n"
2699
" Max value for box ratio (final/init): %f\n"
2700
" Max level allowed: %d\n\n",
2701
bt->n_children, (int)(bt->threshold),
2702
(double)(bt->max_box_ratio), (int)(bt->max_level));
2704
bft_printf(" Max level reached: %5u\n"
2705
" Number of leaves: %10llu\n"
2706
" Leaves with n_boxes > max_n_boxes: %10llu\n"
2707
" Initial number of boxes: %10llu\n"
2708
" Number of linked boxes: %10llu\n"
2709
" Mean number of leaves per box: %10.4g\n\n",
2710
g_max_level_reached, (unsigned long long)n_g_leaves,
2711
(unsigned long long)n_g_spill_leaves, (unsigned long long)n_g_boxes,
2712
(unsigned long long)n_g_linked_boxes, box_ratio);
2714
bft_printf("Number of linked boxes per box tree leaf:\n"
2715
" Mean value: %10.4g\n"
2716
" min. value: %10llu\n"
2717
" max. value: %10llu\n\n",
2719
(unsigned long long)(s.min_linked_boxes),
2720
(unsigned long long)(s.max_linked_boxes));
2722
if (delta > 0) { /* Number of elements in each subdivision */
2724
for (i = 0, j = 1; i < n_steps - 1; i++, j++)
2725
bft_printf(" %3d : [ %10llu; %10llu [ = %10llu\n",
2727
(unsigned long long)(g_min_linked_boxes + i*step),
2728
(unsigned long long)(g_min_linked_boxes + j*step),
2729
(unsigned long long)(count[i]));
2731
bft_printf(" %3d : [ %10llu; %10llu ] = %10llu\n",
2733
(unsigned long long)(g_min_linked_boxes + (n_steps - 1)*step),
2734
(unsigned long long)(g_max_linked_boxes),
2735
(unsigned long long)(count[n_steps - 1]));
2740
/*----------------------------------------------------------------------------
2741
* Dump an fvm_box_tree_t structure.
2744
* bt <-- pointer to box tree structure
2745
*----------------------------------------------------------------------------*/
2748
fvm_box_tree_dump(fvm_box_tree_t *bt)
2750
fvm_box_tree_stats_t s;
2753
bft_printf("\nBox tree: nil\n");
2757
bft_printf("\nBox tree: %p\n\n", bt);
2759
bft_printf(" n_max_nodes: %d\n\n"
2761
(int)(bt->n_max_nodes), (int)(bt->n_nodes));
2765
/* Print statistics and bounding boxes histogram */
2767
bft_printf(" Number of children per leaf: %d\n"
2768
" Max number of bounding boxes for a leaf: %d\n"
2769
" Max value for box ratio (linked/init): %f\n"
2770
" Max level allowed: %d\n\n",
2771
bt->n_children, (int)(bt->threshold),
2772
(double)(bt->max_box_ratio), (int)(bt->max_level));
2774
bft_printf(" Max level reached: %5u\n"
2775
" Number of leaves: %10llu\n"
2776
" Leaves with n_boxes > max_n_boxes: %10llu\n"
2777
" Initial number of boxes: %10llu\n"
2778
" Number of linked boxes: %10llu\n",
2779
s.max_level_reached,
2780
(unsigned long long)(s.n_leaves),
2781
(unsigned long long)(s.n_spill_leaves),
2782
(unsigned long long)(s.n_boxes),
2783
(unsigned long long)(s.n_linked_boxes));
2785
bft_printf("Bounding boxes related to each leaf of the box tree.\n"
2786
" min. value: %10llu\n"
2787
" max. value: %10llu\n\n",
2788
(unsigned long long)(s.min_linked_boxes),
2789
(unsigned long long)(s.max_linked_boxes));
2794
/*----------------------------------------------------------------------------*/
2798
#endif /* __cplusplus */
added
removed
src/fvm/fvm_box_tree.c
