Timothy  0.9
Tissue Modelling Framework
 All Data Structures Files Functions Variables Typedefs Macros
tree.c File Reference

contains tree build functions More...

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include "global.h"
#include "inline.h"
Include dependency graph for tree.c:

Go to the source code of this file.

Functions

void treeBuildInitial (int depth, struct bht_node **iroots, struct bht_node *node)
 
void treeBuild ()
 
void treeFree ()
 

Variables

double tci [8][3]
 

Detailed Description

contains tree build functions

Definition in file tree.c.

Function Documentation

void treeBuild ( )

This function builds an octree for local cells.

Definition at line 181 of file tree.c.

References bht_root, boxmax, boxmin, boxsize, cells, bht_node::center, bht_node::child, bht_node::father, bht_node::isempty, bht_node::leaf, leafs, bht_node::len, lnc, MPIrank, ni, nnodes, nnodes_init, nx, OMPthreads, bht_node::partnum, s, sdim, statOutStep, step, stopRun(), tnc, tree, treeBuildInitial(), cellData::x, cellData::y, and cellData::z.

182 {
183 
184  int i, j;
185  int pow2, idepth, pparts, k, l;
186  int maxTreeSize;
187  struct bht_node *node;
188  struct bht_node **iroots;
189  int nnodes_per_thread = 0;
190 
191  if (lnc == 0)
192  return;
193 #ifdef DEBUG
194  if (MPIrank == 0 && !(step % statOutStep)) {
195  printf(" Tree build...");
196  fflush(stdout);
197  }
198 #endif
199  if (sdim == 3)
200  tnc = 8;
201  if (sdim == 2)
202  tnc = 4;
203 
204  nnodes = 0;
205  nnodes_init = 0;
206 
207  pow2 = 0;
208  k = 1;
209  while (1) {
210  if (k >= OMPthreads || k >= lnc)
211  break;
212  k = k * 2;
213  pow2 += 1;
214  }
215 
216  idepth = 0;
217  l = 1;
218  while (1) {
219  if (l >= k)
220  break;
221  l = l * tnc;
222  idepth += 1;
223  }
224 
225  if (!
226  (iroots =
227  (struct bht_node **) malloc(l * sizeof(struct bht_node *))))
228  stopRun(106, "iroots", __FILE__, __LINE__);
229 
230  pparts = l / k;
231 
232  /* allocate memory for the leafs table */
233  if (!
234  (leafs =
235  (struct bht_node **) malloc(lnc * sizeof(struct bht_node *))))
236  stopRun(106, "leafs", __FILE__, __LINE__);
237 
238  maxTreeSize = lnc * 16;
239 
240  if (!
241  (tree =
242  (struct bht_node *) malloc(maxTreeSize * sizeof(struct bht_node))))
243  stopRun(106, "tree", __FILE__, __LINE__);
244 
245  bht_root = &tree[0];
246 
247  nnodes_init = 1;
248 
249 
250  boxmin[0] = nx;
251  boxmin[1] = nx;
252  boxmin[2] = nx;
253  boxmax[0] = -nx;
254  boxmax[1] = -nx;
255  boxmax[2] = -nx;
256 
257  for (i = 0; i < lnc; i++) {
258  boxmin[0] = (cells[i].x < boxmin[0] ? cells[i].x : boxmin[0]);
259  boxmax[0] = (cells[i].x > boxmax[0] ? cells[i].x : boxmax[0]);
260  boxmin[1] = (cells[i].y < boxmin[1] ? cells[i].y : boxmin[1]);
261  boxmax[1] = (cells[i].y > boxmax[1] ? cells[i].y : boxmax[1]);
262  boxmin[2] = (cells[i].z < boxmin[2] ? cells[i].z : boxmin[2]);
263  boxmax[2] = (cells[i].z > boxmax[2] ? cells[i].z : boxmax[2]);
264  }
265 
266  boxsize = boxmax[0] - boxmin[0];
267  boxsize =
268  (boxmax[1] - boxmin[1] > boxsize ? boxmax[1] - boxmin[1] : boxsize);
269  boxsize =
270  (boxmax[2] - boxmin[2] > boxsize ? boxmax[2] - boxmin[2] : boxsize);
271  boxsize += 0.1;
272 
273  /* set up the root node */
274  bht_root->len = boxsize;
275  bht_root->center[0] = boxmin[0] + (boxmax[0] - boxmin[0]) / 2.0;
276  bht_root->center[1] = boxmin[1] + (boxmax[1] - boxmin[1]) / 2.0;
277  if (sdim == 3)
278  bht_root->center[2] = boxmin[2] + (boxmax[2] - boxmin[2]) / 2.0;
279  if (sdim == 2)
280  bht_root->center[2] = 0.0;
281 
282  bht_root->leaf = 1;
283  bht_root->isempty = 1;
284  bht_root->father = NULL;
285 
286  for (i = 0; i < tnc; i++)
287  bht_root->child[i] = NULL;
288 
289  ni = 0;
290 
291  treeBuildInitial(idepth, iroots, bht_root);
292 
293  ni = nnodes_init;
294 
295  /* stop if there are no particles in this process */
296  if (lnc == 0)
297  return;
298 
299 #pragma omp parallel num_threads(k) private(node,i,j) firstprivate(pparts) shared(ni)
300  {
301  int OMPid;
302 
303  OMPid = omp_get_thread_num();
304 
305  node = bht_root;
306 
307  /* loop over cells - each cell is inserted into its own octree location */
308  for (i = 0; i < lnc; i++) {
309  int mycell = 0;
310  int c;
311  int s;
312  struct bht_node *no;
313 
314  no = node;
315 
316  for (j = 0; j < pparts; j++) {
317  mycell = cellInNode(i, iroots[OMPid * pparts + j]);
318  if (mycell == 1) {
319  no = iroots[OMPid * pparts + j];
320  break;
321  }
322  }
323 
324  if (mycell == 0)
325  continue;
326 
327  /* this loop iterates until given cell p is inserted into its location */
328  while (1) {
329  /* if this node is not a leaf we continue searching */
330  if (no->leaf == 0) {
331  c = whichChildNode(i, no);
332  no = no->child[c];
333  /* if this node is a leaf and it is not empty we need to move cells */
334  } else if (no->leaf == 1 && no->isempty == 0) {
335  for (s = 0; s < tnc; s++)
336  no->child[s] = NULL;
337  no->leaf = 0;
338  c = whichChildNode(no->partnum, no);
339  no->child[c]->partnum = no->partnum;
340  no->child[c]->isempty = 0;
341  leafs[no->partnum] = no->child[c];
342  c = whichChildNode(i, no);
343  no = no->child[c];
344  /* if this node is a leaf and it is empty we insert cell into the node */
345  } else {
346  no->partnum = i;
347  no->isempty = 0;
348  leafs[i] = no;
349  break;
350  }
351  }
352 
353  }
354 
355  } /* end of OpenMP parallel region */
356 
357  free(iroots);
358 #ifdef DEBUG
359  if (MPIrank == 0 && !(step % statOutStep)) {
360  printf("done\n");
361  fflush(stdout);
362  }
363 #endif
364 }
int statOutStep
Definition: global.h:214
int MPIrank
Definition: global.h:134
#define lnc
Definition: global.h:102
float s
Definition: global.h:186
unsigned char leaf
Definition: global.h:284
int nnodes
Definition: global.h:292
int ni
Definition: global.h:297
unsigned char isempty
Definition: global.h:285
int sdim
Definition: global.h:160
int OMPthreads
Definition: global.h:137
int tnc
Definition: global.h:295
int step
Definition: global.h:173
struct bht_node * child[8]
Definition: global.h:283
double x
Definition: global.h:69
double center[3]
Definition: global.h:280
int partnum
Definition: global.h:281
struct bht_node * bht_root
Definition: global.h:289
struct bht_node ** leafs
Definition: global.h:290
struct cellData * cells
Definition: global.h:82
struct bht_node * father
Definition: global.h:282
double z
Definition: global.h:71
double boxmax[3]
Definition: global.h:244
double y
Definition: global.h:70
struct bht_node * tree
Definition: global.h:288
double boxsize
Definition: global.h:245
int nx
Definition: global.h:162
void stopRun(int ierr, char *name, char *file, int line)
Definition: utils.c:72
double len
Definition: global.h:279
void treeBuildInitial(int depth, struct bht_node **iroots, struct bht_node *node)
Definition: tree.c:145
int nnodes_init
Definition: global.h:293
double boxmin[3]
Definition: global.h:244

Here is the call graph for this function:

Here is the caller graph for this function:

void treeBuildInitial ( int  depth,
struct bht_node **  iroots,
struct bht_node node 
)

This function builds an initial tree which is needed for OpenMP parallelization

Definition at line 145 of file tree.c.

References bht_node::center, bht_node::child, bht_node::father, bht_node::isempty, bht_node::leaf, bht_node::len, ni, nnodes_init, bht_node::partnum, sdim, tci, tnc, and tree.

147 {
148  int i, j, k;
149  if (depth != 0) {
150  double shift;
151  shift = node->len * 0.25;
152  for (j = 0; j < tnc; j++) {
153  node->child[j] = &tree[nnodes_init];
154  node->child[j]->len = node->len * 0.5;
155  node->child[j]->center[0] = node->center[0] + tci[j][0] * shift;
156  node->child[j]->center[1] = node->center[1] + tci[j][1] * shift;
157  if (sdim == 3)
158  node->child[j]->center[2] = node->center[2] + tci[j][2] * shift;
159  if (sdim == 2)
160  node->child[j]->center[2] = 0.0;
161  node->child[j]->father = node;
162  node->child[j]->isempty = 1;
163  node->child[j]->leaf = 1;
164  node->leaf = 0;
165  for (k = 0; k < tnc; k++)
166  node->child[j]->child[k] = NULL;
167  nnodes_init++;
168  treeBuildInitial(depth - 1, iroots, node->child[j]);
169  }
170  } else {
171  node->leaf = 1;
172  node->partnum = -1; /* indicates an empty node */
173  iroots[ni] = node;
174  ni++;
175  }
176 }
double tci[8][3]
Definition: tree.c:35
unsigned char leaf
Definition: global.h:284
int ni
Definition: global.h:297
unsigned char isempty
Definition: global.h:285
int sdim
Definition: global.h:160
int tnc
Definition: global.h:295
struct bht_node * child[8]
Definition: global.h:283
double center[3]
Definition: global.h:280
int partnum
Definition: global.h:281
struct bht_node * father
Definition: global.h:282
struct bht_node * tree
Definition: global.h:288
double len
Definition: global.h:279
void treeBuildInitial(int depth, struct bht_node **iroots, struct bht_node *node)
Definition: tree.c:145
int nnodes_init
Definition: global.h:293

Here is the caller graph for this function:

void treeFree ( )

This function dealocates arrays used in tree building

Definition at line 369 of file tree.c.

References leafs, and tree.

370 {
371  free(tree);
372  free(leafs);
373 }
struct bht_node ** leafs
Definition: global.h:290
struct bht_node * tree
Definition: global.h:288

Here is the caller graph for this function:

Variable Documentation

double tci[8][3]
Initial value:
= { {-1, -1, -1}, {-1, 1, -1}, {1, -1, -1}, {1, 1, -1},
{-1, -1, 1}, {-1, 1, 1}, {1, -1, 1}, {1, 1, 1},
}

Definition at line 35 of file tree.c.