Implementing a quadtree using arrays

Question

I am trying to write code for simulating n-body problem using the Barnes-Hut tree algorithm. I plan on using CUDA in the future and thus want my quadtree data structure to not be composed of heap objects.

From the paper "An Efficient CUDA Implementation of the Tree-Based Barnes Hut n-Body Algorithm" by Martin Burtscher and Keshav Pingali (sorry couldn't find link) the authors state:

Dynamic data structures such as trees are typically built from heap objects, where each heap object contains multiple fields, e.g., child-pointer and data fields, and is allocated dynamically. Because dynamic allocation of and accesses to heap objects tend to be slow, we use an array-based data structure. Accesses to arrays cannot be coalesced if the array elements are objects with multiple fields, so we use several aligned scalar arrays, one per field, as outlined in Figure 6.6. As a consequence, our code uses array indexes instead of pointers to tree nodes.

I understand the part about aligned scalar arrays (i.e. SOA vs AOS idiom in parallel computing), but unfortunately the authors do not explain how one constructs the quadtree using arrays.

My question is how does one implement a quadtree data structure (with methods for inserting spacial points) using arrays? I know how to implement a quadtree in the traditional way using node structs and pointers for children etc. Can someone provide a reference that details how to do this using arrays. Even information on how to implement a binary tree (or any tree really) using arrays might be useful here.

Keaton

en.wikipedia.org/wiki/Binary_tree#Arrays – Autochthon

Answer 1

Implementing a binary tree using array is really simple, firstly we are going to index the array from 1, i.e the root node is going to be 1, and then the

left child will be : leftChildIndex = 2 * parentIndex;

right child will be at : rightChildIndex = 2 * parentIndex + 1;

Now if one want to find the parent of current node : parentIndex = currIndex/2;

I have written a c++ code that does the preorder traversal of the tree :

#include<iostream>

using namespace std;

int binaryTree[20], lengthOfTree;
int leftChild(int idx){ return 2*idx; }
int rightChild(int idx){ return 2*idx+1; }
int parentIndex(int idx){ return idx/2; }

void traverseTree(int idx){
    if(idx >= lengthOfTree) return;
    cout << binaryTree[idx] << " ";
    traverseTree(leftChild(idx));
    traverseTree(rightChild(idx));
}

int main(){

    lengthOfTree = 15;
    for(int i = 1;i <= lengthOfTree;i++){
        cin >> binaryTree[i];  
    }
    traverseTree(1);
    cout << endl;

return 0;
}

Link to solution on Ideone : http://ideone.com/ZpTJCa

----------------------------------------------------------------------------

Indexing a quad tree might be a little more complex, so what we can do is again index the tree from 1, for each node we can find the level of that node, eg : the level of node 1 is 0, level of node 4 is 1, level of node 11 is 2.

pseudoCode for finding level : O(log n)

int findLevel(int nodeNo){
    int level = 0;
    int currNode = 1;
    while(currNode < nodeNo){
        currNode = currNode + pow(4, level++);
    }
    return level;
}

Similarly the leftmost node of the current level and the rightmost node of the currentlevel can be calculated using the above pseudocode, then to find the 4 children of the current node we can do :

1st child of current node : child1 = (rightmostNode - currentNode) + 4 * (currentNode - leftmostNode);

2nd child of current node : child2 = child1 + 1;

3rd child of current node : child3 = child2 + 1;

4th child of current node : child4 = child3 + 1;

Also you can create a mapping for finding the parent.

Houri

This is a good idea although I think it is only efficient for complete binary trees. I could extend this idea for quad trees, but if many of the child nodes are null then there will be a lot of wasted memory I think? For example I think the array size would have to grow something like 2^d where d is the tree depth. For a quad tree this gets worse (now 4^d)... – Keaton

Yes that is correct a lot of nodes would be wasted while creating a quad tree, and the tree would grow exponentially, that would be very inefficient. – Houri

Answer 2

A quad tree, represented as an array is called a "Linear Quadtree". Using this term you will find some literature.

A paper from Hannan Samet recommends to first implement an application using a tradional quadtree, then check whether a linear quadtree approach will work. Not all applications can use a linear quadtree.

"(with methods for inserting spatial points)"

Such linear approaches often demand a static quad tree, that is one which does not change its content. Same applies for GPU applications, they demand a (huge) static data set, where millions of operations are executed. The access (uploading time) to the GPU is relative slow, so the type of the application should one where the most data will not change (for some time).

Yurt

Answer 3

One solution might be to pack the quadtree into a binary tree, for example by using space-filling curves. I find the z-curve (morton order / z-order) easiest to use. For a z-ordering you need to interleave the bits of your coordinates so that for example two 64 value (x,y) are interleaved into a single 128 bit value. The 128 bit value can then be store in a binary tree or a trie. In C++ there should be opcode for efficient interleaving of bits (I think it's called 'packing'?).

You can also do the interleaving with floating point values, see Section 3.3 in the linked PDF. It shows how floating point values can be quickly converted to a integer-format and back without loss of precision. The example code (Java) is taken from here, I believe in C++, instead of using Double.doubleToRawLongBits, you could simply cast a float to an integer and then apply the same transformation as shown here:

public static long toSortableLong(double value) {
    long r = Double.doubleToRawLongBits(value);
    return (r >= 0) ? r : r ^ 0x7FFFFFFFFFFFFFFFL;
}

public static double toDouble(long value) {
    return Double.longBitsToDouble(value >= 0.0 ? value : value ^ 0x7FFFFFFFFFFFFFFFL);
}

EDIT

The code for interleaving bits could look like this (from here):

unsigned short x;   // Interleave bits of x and y, so that all of the
unsigned short y;   // bits of x are in the even positions and y in the odd;
unsigned int z = 0; // z gets the resulting Morton Number.

for (int i = 0; i < sizeof(x) * CHAR_BIT; i++) // unroll for more speed...
{
  z |= (x & 1U << i) << i | (y & 1U << i) << (i + 1);
}

The above solution is the obvious approach. There are other versions in the link above that may be even faster.

As I understand, modern CPUs also have operations for interleaving (also called 'interlacing'), for example with 'shuffling'. Another algorithm and some more info can be found here

EDIT end

Once you have the binary tree, you can apply another mapping to get an array. Have also a look at Binary Heaps, they are usually implemented as arrays.

Kashgar

An Z order interleaving transforms a 2 dimensional coordinate into a one dimensional, as you have correctly explained. Unfortunately the code does not match that. soemthing like long zOrderIndex(int x, int y) is needed. – Yurt

No, the code describes the transformation that needs to be performed on floating point values before you interleave them. I'll update the answer with the interleaving code. – Kashgar