Binarized octree generation for Cartesian adaptivemesh refinement around immersed geometriesJaber J. Hasbestan1, Inanc Senocak2Department of Mechanical Engineering and Materials Science,University of Pittsburgh, Pittsburgh, PA 15261, USAAbstractWe revisit the generation of balanced octrees for adaptive mesh refinement(AMR) of Cartesian domains with immersed complex geometries. In a recentshort note [Hasbestan and Senocak, J. Comput. Phys. vol. 351:473-477 (2017)],we showed that the data-locality of the Z-order curve in hashed linear octreegeneration methods may not be perfect because of potential collisions in thehash table. Building on that observation, we propose a binarized octree gener-ation method that complies with the Z-order curve exactly. Similar to a hashedlinear octree generation method, we use Morton encoding to index the nodes ofan octree, but use a red-black tree in place of the hash table. Red-black treeis a special kind of a binary tree, which we use for insertion and deletion ofelements during mesh adaptation. By strictly working with the bitwise repre-sentation of the octree, we remove computer hardware limitations on the depthof adaptation on a single processor. Additionally, we introduce a geometry en-coding technique for rapidly tagging the solid geometry for refinement. Ourresults for several geometries with different levels of adaptations show that thebinarized octree generation outperforms the linear octree generation in terms ofruntime performance at the expense of only a slight increase in memory usage.We provide the current AMR capability as open-source software.Keywords: AMR, Hash table, Octree, Red-black tree, Z-order curve1. IntroductionMesh adaptation is a powerful technique in the solution of multiphysics prob-lems governed by coupled partial differential equations (PDEs). By refining thespatial resolution where needed, researchers can realize substantial savings incomputer memory and solution turnaround time. Unlike the solution algorithm,mesh adaptation is not governed by a set of governing equations. Therefore, dif-ferent approaches and strategies have been proposed for mesh adaptation [1].1Postdoctoral Research Associate, jaber@pitt.edu2Associate Professor, senocak@pitt.eduPreprint submitted to J. Comput. Phys. December 4, 2017arXiv:1712.00408v1 [cs.CE] 1 Dec 2017In this study our focus is on the so-called structured adaptive mesh refinement(AMR) method. Berger and Oliger [2] pioneered the use of patch-based adapta-tion for the solution of hyperbolic PDEs with refinement in both space and time.Berger and Colella [3] further improved the AMR algorithm and used an esti-mate of the local truncation error to guide the refinement decisions. Althoughthese pioneering works were restricted to an explicit finite-difference method,they are one of the early works that introduced efficient use of data-structuresin scientific computing and the strategy developed in these works is often re-ferred to as Berger-Oliger-Colella adaptive refinement. This unique patch-basedstrategy has been adopted in Boxlib [4], PARAMESH [5], AMREx [6], Chombo[7], ENZO [8], SAMRAI [9, 10] and AMRClaw [11] packages. In a patch-basedAMR, a hierarchy of patches are nested within a coarse level mesh that occupiesthe entire domain. A refinement ratio is assigned between successive patches.During refinement, patches are tagged based on a specified criterion for refine-ment. Ghost-cells are used at patch boundaries to compute the solution acrossthe patches.Implementation of an AMR algorithm is major software development effort.It is common to refer to the end product as an AMR infrastructure. In build-ing this software infrastructure, the design goal is to introduce an abstractionfor the meshing part such that developers of the solution algorithm need notknow the details of the computer implementation. For static-mesh problems,the efficiency of the AMR infrastructure may not be critical. But for large-scaleproblems in which the refinement criteria is dynamic, time spent on meshing,while the simulation is evolving, can be prohibitively expensive. Additionally,the memory hierarchy on modern computing hardware favors computer imple-mentations that exploit data-locality and use less memory. It is worth men-tioning that these performance issues related to adaptive meshing have beenthe focus of several studies in the computer graphics field, driven by real-timeapplications [12]. There are a number of open-source packages that provideAMR capabilities for meshing and solving PDEs. Dubey et al. [13] survey arepresentative set of block-structured AMR infrastructures.An alternative strategy to the patch-based AMR approach is a tree-basedstrategy. Octree is a hierarchical data structure that has been the subject ofmany studies [14, 15]. The tree structure has also been adopted in the solutionof partial differential equations (PDEs) with adaptive resolutions. In its classicalrepresentation, an octree keeps track of parent-child relations with pointers. Inthe three-dimensional case, eight pointers are needed to store the neighbors.Sundar et al. [16] developed parallel algorithms to generate balanced linearoctrees. A linear octree is a representation of an octree in which a linear arrayis used instead of a tree data structure to store the leaf nodes [17, 18]. It iscommon to use space-filling curves to represent linear octrees. Sundar et al.used Morton encoding as one of the many ways to generate a space-filling curvebecause of its ease of computer implementation. Their method addressed theso-called ripple-effect associated with balancing the octree.Burstedde et al. [19] pursued the linear octree representation approach de-scribed in Sundar et al. [16] and developed p4est, which is an AMR packagefor multiple connected octrees. Authors refer to these connected octrees asforest-of-octrees. An encoding scheme was introduced to enable inter-octreeconnectivity. Octrees and octants were numbered using a space-filling Z-ordercurve, which was effective in parallel load balancing. p4est has been shownscale across many processors. In one exercise, a mesh with 5.13× 1011 elementswas created using 220,320 processes. Isaac et al. [20] introduced recursive algo-rithms for this forest-of-octrees approach to traverse a linear octree and furtherimprove the parallel scalability of p4est aimed at improving the 2:1 balancealgorithm.Another approach is the building-cube method (BCM) [21, 22], which isa combination of explicit tree and unstructured grid methodology. In BCM,each cube has its own equally-spaced Cartesian grid. All cubes have the samenumber of Cartesian mesh. Therefore, the resolution is determined by the levelof adaptation of the cube.Having a balanced tree is important in numerical solution of PDEs as itenables one to gradually traverse from low-resolution to high-resolution withoutsignificantly degrading numerical accuracy, especially for simulation of multi-scale problems such as turbulent flows. A balanced-tree is also needed to enablea generic implementation of numerical schemes for the solution of PDEs in acode because coarse to fine blocks obeys the same rule. The common approachto enforce a balanced tree is to confine the maximum level of the element ateach face to 2d−1 where d is the dimension of the computational space (i.e. 2or 3). The difficulty of the 2:1 balancing algorithm emanates from the fact thatthe domain of influence of any element under consideration for refinement mightextend beyond the immediate neighboring elements. This issue has been coinedas the ripple effect in the literature [16]. Different methods have been proposedto improve the efficiency of enforcing the 2:1 balance condition.Tu and Ohallaron [23] demonstrate 2:1 balanced refinement of octrees. Theirsolution consists of two steps; balance by parts and prioritize the ripple effect.The key idea in their approach is to bulk load most of the data on memoryand enforce the balance locally. They load a small region to a pointer-based(sub)octree in memory and resolve the 2:1 constraint and the ripple effect with-out worrying about octants outside of the region. This is due to the fact thatmost impact caused by a tiny octant tagged to be refined is localized in a smallregion. Their balance-by-parts algorithm is described as follows. Imagine a mov-ing window inside the 3D domain, at any given moment, the octants that residein this window are retrieved from disk and cached in a (sub)octree. Then, theyadjust the data structure to enforce 2-to-1 constraint in memory, the contenton disk is updated accordingly.Ripple effect is even more challenging in distributed systems due to thecomplicated communication pattern. Sundar et al. [16] proposed a bottom-upapproach to enforce the balance condition. This approach eliminates the needfor synchronization in parallel implementation. Isaac et al. [20] improved thebalance algorithm in the p4est package. In their approach, they first reducethe input set to a compressed equivalent and iteratively add a sparse set ofoctants to that list to enforce the balance.Aside from solving the PDEs, octrees are extensively used in computergraphics. Unlike in the solution of PDEs, balanced trees are not required incomputer graphics applications, and the desired resolution can be much lowerthan resolutions needed for turbulent flow simulations. Driven by interactivityrequirements of computer games and ray tracing applications, the implementa-tion of octree on central processing units (CPUs) or graphics processing units(GPUs) can be highly optimized. Morton codes are used in building octrees,mapping multi-dimensional data to one dimension and thereby enhancing data-locality. Use of a hash table with Morton encoding is common in generatingthe octrees. Meagher [14] utilized octree data structure for modeling the three-dimensional complex geometries including geometries with holes and disjointedparts and sculptured surfaces. Frisken and Perry [24] demonstrate efficientmethods for tree structure traversals.Lewiner et al. [12] developed a fast pointerless method to generate dual oc-trees. They proposed strategies for dual generation of static or dynamic pointer-less octrees, where memory usage is lowered by several factors compared to theusual recursive generation methodologies. Note that pointerless representationis another terminology for hashed linear octrees. However, their study is notrelated to solving PDEs and hence, the generated mesh is not necessarily bal-anced. Lewiner et al. provide a good review on different approaches to octreegeneration in the same work.Schwarz and Seidel [25] describe octree implementation on GPUs. Baertet al. [26] uses bit-interleaving in constructing the Morton code to map everyvoxel center to a one-dimensional space to generate out-of-core octrees for aCartesian mesh, the method presented there is not directly applicable to AMR,as they use integer indices in the mapping. The bit-interleaving introduced inthat study uses fixed point arithmetic to store the center of the cube. Bit-interleaving is also introduced by Tero Karras [27]. The problem with theseapproaches for solving PDEs is that one needs to always scale the geometry tounit length and extracting the connectivity is not trivial.Different versions of space-filling curves have been used for partitioning theoctrees [28, 29, 30, 31]. Z-order curve, also referred to as Morton code [32],is often preferred because of its ease of computer implementation [16]. Animportant feature of using Morton code is that it ensures that the elementswith closer Morton codes are stored close to each other in the computer memory.This feature makes the search for a given key efficient, since one only needs tolook for the elements close to the current key.The method used to search for a key is also critical for performance. In arecent short note [33], we have shown that a red-black tree helps preserve the Z-order curve exactly, whereas the Z-order can be broken because of collisions in ahashed linear octree. A red-black tree also supports deep-level mesh adaptationfor extreme-scale problems because the method is based on bitwise comparisonsonly.Building on those useful features of the red-black tree data structure, wepresent a fast, low-storage AMR algorithm that complies with the Z-order curvestrictly. Our presentation includes the algorithms for 2:1 balance, geometryencoding for immersion of complex geometries, and the tricks and trades ofcomputer implementation that exploits bitwise representation heavily. In ourapproach, we do not store the neighborhood connectivity explicitly. To thebest of our knowledge, our open-source software, which we call rebl-AMR,is the first Cartesian adaptive mesh generation package that make use of thefundamental red-black tree data structure at its core. We present numerousexamples to demonstrate the performance of red-black tree version and compareits computational performance to a hashed linear octree. For completeness,we compare the performance results from both methods against a graph-basedunstructured mesh approach as well. The present rebl-AMR [34] package isavailable as open-source software.2. Octree representationsAn octree is constructed by isotropic refinement where a rectangular volumeis subdivided in half in x, y and z directions. Refinement of each elementgenerates eight more elements in three dimensions. We refer to these elements assiblings. Complicated geometries can be represented using octrees by recursivelysubdividing the elements. Different types of data structures have been used tostore the tree data. The selection of the data structures is a crucial step inoctree generation algorithms because it will impact the execution time as wellas memory usage. For dynamic adaptive mesh generation problems where thegeometry deforms or moves in time, it is desirable to have constant-in-timeinsertion and removal of mesh elements.One of the early approaches in octree generation is the use of parent-childrelation in creating a so-called explicit octree, where use of pointers to childrenis standard. This approach is the classical way of storing the tree and uses eightpointers for octree representation [35]. A sample four-level tree is presentedin Fig. 4. The number of pointers needed to represent a large adaptivelyrefined mesh would be unnecessarily large as well and would also complicatethe computer implementation. Pointerless implicit trees have been proposed asan alternative to the explicit tree representation [17, 18].We present three methods to generate octrees. The first method is thehashed linear octree generation. Although this method is very efficient, there isno guarantee to preserve the Z-order curve once it is used in conjunction with ahash table because of collision, which we have demonstrated recently [33]. Thesecond method is unique to the present work, which we call it the binarizedoctree generation. Similar to a hashed linear octree generation method, werepresent octree nodes on a Z-order curve and use a red-black tree for elementsearch, insertion and deletion operations. A red-black is a special kind of binary-tree and bitwise representation is used to implement the Z-order curve. Becausewe use a binary tree data structure to store the leaves of the octree instead ofan array, we refer to this approach as binarized octree generation. Compared tothe classical approach to generate an octree, both the hashed linear octree andthe binarized octree use much less memory.In the third approach, we generate an octree as an unstructured mesh. Theconnectivity is generated for every element using a graph-based approach. Thismethod updates the connectivity directly without the need to actually store thetree information explicitly. The unstructured mesh approach is provided forcompleteness. We describe that approach in Hasbestan and Senocak [36] and itis not the main focus of the current work.Z-order curve, also referred to as the Morton code, is common to both thelinear and binarized octree generation methods. Therefore, we first explain theconstruction of Morton codes and how we manipulate them for neighborhoodconnectivity information.2.1. Morton code generationMorton encoding is a mapping from a multi-dimensional space to one di-mension [32]. When generating a Morton code, first, a bit code is constructedfor every node. This node is then converted to an integer, if needed. The nodes,once laid out, follow a Z-order curve, which enhances data-locality.Figure 1 presents the Morton ordering for a mesh in three dimensions, whichis refined one level only. Note that we colored the four lower octants on theirbottom surface while coloring the top four octants on their top surfaces to avoidvisual clutter. The bit representation on each surface is valid for the individualoctant. The Morton code for each octant is constructed by specifying its relativelocation index in x, y, z directions in the Cartesian coordinates, respectively.000 100001 101011 111110010XYZFigure 1: Morton code construction at a single level based on element location in x, y, zdirectionsLet us now demonstrate refinement. At each level, the Morton code corre-sponds to three bits in three-dimensional space and two bits in two-dimensionalspace. The Morton code after refinement is generated by concatenating everybits at every level, starting from the left with the top level node in the treeand placing each node’s bit representation to the right as we traverse to thebottom of the tree. A two-dimensional mesh with two levels of adaptation ispresented in Fig. 2 to illustrate Morton code generation with refinement. Here,we considered a two-dimensional case for ease of presentation. The final code isgenerated by concatenating two bits per level of adaptation. The length of thecodes are extended to 4 bits by zero padding to be consistent.xy000001001000110010101001 1011Figure 2: Morton code representation of a two-dimensional mesh with two levels of adaptation.An octree for a three-dimensional mesh with three levels of adaptation isshown in Fig. 3. Here, the procedure is illustrated with coordinated colors fornode 23. The Morton code for node number 23 is 001, 000, 101 as shown in Fig.4. Note that the change of color in Fig. 4 coordinates with Fig. 3 to highlightplacement of bits from left to right.Morton code enables us to keep the elements that have closer keys allocatedclose to each other in the hash table. It is important to mention that the cacheperformance also improves due to this locality. There are several methods forMorton encoding [27, 12].The procedure described above limits the number of levels for adaptationto 10 for natively 32-bit systems and to 21 for natively 64-bit systems, Ourstatement is valid if the AMR algorithm at hand has to convert the Mortoncode to an integer number because it would produce an overflow.We should note that in the hashed linear octree approach, any hash functionthat does not use the entire bits and only selects some portion of the Mortoncode would be able to attain adaptation levels deeper than the aforementionedvalues at the expense of performance due to collisions in the hash table. As wehave shown in Hasbestan and Senocak [33] Z-order curve is no longer retainedonce collisions occur. Lewiner et al. [12] discuss several approaches pertinent tothe use of hash table with Morton codes.As we discuss later in this paper, we can overcome this limitation on depthof adaptation without a performance hit by using the red-black tree data struc-ture in the AMR algorithm. C++ bitset class enables defining the numberof required bits beyond what is natively supported by the hardware, and thered-black tree operates on bit comparison and therefore there is no need for000(0)000(1)100(2)010(3)001(4)000(9)000(17)100(18)010(19)001(20)110(21)011(22)101(23)111(24)100(10)010(11)001(12)101(13)011(14)110(15)111(16)101(5)011(6)110(7)111(8)Figure 3: Morton code construction procedure illustrated for the node 23, The numbers inthe parentheses represent the node index as an integer number001,000,101Figure 4: Morton code constructed for node 23 after three levels of refinementconverting the Morton code to an integer.2.2. Hashed linear octree generationData-locality is key to superior memory performance. As we mentionedearlier, one option to enhance data-locality is to store data using space-fillingcurves. The use of Z-order curve [32] is popular because of its easy computerimplementation. Use of hash tables is also common in pointerless linear octreegeneration. A hash table potentially provide O(1) access to its elements owingto the use of array as the underlying data-structure. However, O(1) access isonly guaranteed as long as there are no collisions in the hash table. In practice,this requirement is difficult to fulfill and collisions do happen when the hashfunction generates the same index for different keys. Typically, a remainderoperator is used as the hash function to fit the keys inside a given array.The state-of-the-art in linear octree generation is to use a hash table withMorton codes. This choice combines the efficiency of a hash table with data-locality offered by the Z-order curve to further improve the performance of theoctree generation. Morton codes are used as keys to the hash table. [26, 12].The works by Sundar et al. [16] and Burstedde et al. [19] fall under this category.There are several ways to construct the Morton code. In our approach,Morton codes are generated recursively from the octree hierarchy and carry allthe information about the octree with them for each node of the octree. Weexplain the details in section 3.In this part of our implementation, we make use of hash maps, which is awell-known associative container data structure in the field of computer science.A hash map uses a hash function to produce an index for each key. An efficientimplementation of a hash map unordered map is readily available as partof the C++ standard library. C++ unordered map class, is an associativecontainer. It associates every key to a mapped value. The mapped values arestored in buckets in such a way that the elements with key values close to eachother are also stored in nearby buckets, which is good for performance. Whensolving a PDE, each key is mapped to a pointer to the solution vector, whichis dependent on the numerical formulation and the discretization method. Forexample, for a cell-centered finite-volume scheme, this pointer will point to anarray of unknowns in each cell, whereas in the case of a finite-element methodit will point to the unknowns at the nodes of each element.The hash map data structure is schematically illustrated in Fig. 5. (We usethe term hash table to refer to the fundamental data structure, and hash mapto refer to the C++ class.) The only restriction of a hash map is that the keyhas to be unique. The uniqueness condition is readily satisfied for the presentapplication since every node has a unique Morton code. The only disadvantageof the hash map associative container is that rehashing might occur if the finalsize of the container is not well predicted. We note that rehashing is an expensiveoperation computationally. If the size of the array in hash map is big enoughdata-locality can be preserved and perfect hashing would be possible. However,this is hard to achieve in practice for dynamic problems.Figure 5: Illustration of hash map data structure, each key is mapped to an index by hashfunction and each index corresponds to a pointer held in a bucketTheoretically, a hash table can have O(1) time for insertion, removal andsearch operations overits elements. However, this is only possible when there areno collisions, which is known as perfect hashing. Let us explore the conditionsfor perfect hashing in the realm of linear octree representation. First of all, thehash function maps each key to an integer index. It is a common practice to usethe remainder operator in the hash function to fit a given number of elements toan array, because a hash table uses array as its underlying data structure. Notethat the operation here is limited by the maximum integer number representableon the computer architecture at hand. Consider the mesh given in Fig. 6(a).In order to have a perfect hash, we first convert the Morton code to integer toproduce an index and then allocate enough memory. For this simple case thesize of the array should be at least 13To illustrate the collision issue with the hash table approach, we assume thatthe initial size of the hash table is 8. Now let us try to fit the 7 elements shownin Fig. 6(a) in the 8 slots by using the remainder operator in the hash function.For elements with hash values of 4 (4 mod 8 = 4) and 12 (12 mod 8 = 4) thecalculated index is the same. This implies a collision will occur in the fifth slot.In other words, both elements are assigned to the same slot in the hash table.To conform to the Z-order curve and have zero collisions at the same time,the size of the array should be at least 13, which is a substantial increase inallocated memory. We also need to allocate extra space in the array for possiblerefinements that might occur in the next step. Obviously, this remedy to avoidcollision in the hash table can be prohibitively expensive in terms of memorywhen a mesh with a deep-level of adaptation is desired. Subsequently, one hasto relax the no-collision condition and try to use the remainder operator to fitthe mesh elements in the array. Once the no-collision condition is relaxed, O(1)access time is no longer guaranteed. Another draw back of using a hash tableis that rehashing might occur if one does not predict the final octree size with agood approximation, and rehashing is computationally expensive. Any collisionwill also break the Z-order curve compliance, leading to performance loss. Tohave a Z-order curve with a regularly increasing index (i.e., 0, 1, 2, ...) a sortingprocess is inevitable, which is not readily available in a hashed linear octreegeneration method. To address this issue we propose to use the red-black treein place of the hash table.xy000001001000110010101001 1011(a)xy04812109 11(b)Figure 6: Illustration of mesh element indexing that arise in the hash table data structure,(a) Morton code representation of a mesh with two levels of adaptation, (b) Hash valuescorresponding to the Morton code from (a)2.3. Binarized octree generationSelf-balancing binary search trees, such as a red-black tree, can be favoredover a hash table in certain applications, because the worst-case performanceof a red-black tree for insertion, removal and search operations is guaranteed tobe O(logN) . A red-black tree also keeps its elements sorted as opposed to theout-of-order storage of elements in a hash table. Red-black trees are widely usedin different areas of the computer science field [37]. Its applications ranges fromcomputational geometry to process scheduling in Linux kernels. It guaranteesthat the path from the root to the farthest node is no more than twice as longas the path from the root to the nearest node [35, 38]. In the present study,we propose to adopt the red-black tree in place of the hash table in the designof the AMR algorithm. Our goal is to create an AMR algorithm that does notlose its Z-order curve compliance at deep levels of adaptation.A red-black tree is a self-balancing binary search tree with an extra bit fordefining the color of a node. The color can be either red or black for a node.The color bit is needed to keep the tree balanced. A red-black tree obeys thefollowing rules:• The root node is always black• A red node can not have a red child• Every path from a given node to any of its descendant nodes contains thesame number of black nodesThe memory usage of a red-black tree is much less than an explicit octree.For a large mesh, data-locality due to strict compliance with the Z-order curveoffsets the O(logN) element access cost of a red-black tree, as we have demon-strated in [33].An associated container that uses a red-black tree data-structure is readilyavailable in the C++ standard template library (STL) as Set/Map classes.Set class is a special case of the Map class in which the associated value is thekey itself. The red-black tree needs a comparison function to insert or removeelements. We use the Morton code representation of the octree nodes as thecomparison criteria for insertion of the new elements. Note that to comparetwo keys represented by bits, conversion of the key to an integer number is notrequired anymore with a red-black tree. We make use of this feature of thered-black tree in conjuction with the bitset class in C++ STL to generatemeshes with deep levels of adaptation. We use the fast bitwise operation (xor)to perform this comparison.By adopting a red-black in the AMR algorithm, we accomplish the following:1. Remove the limitation on maximum level of adaptation. With a red-blacktree there is no need to convert the Morton code to an integer index.Red-black tree only requires the comparison of the bits.2. Preserve the order of the Z-curve exactly and ensure data-locality in thecomputer implementation for improved performance.2.4. Graph-based unstructured mesh approachThe Cartesian adaptive mesh generation can also be implemented using agraph-based data structure by treating the mesh as an unstructured grid topol-ogy. In this method, the connectivity for each element has to be stored. Fora large mesh with a deep level of adaptation this approach requires a lot ofmemory because of the need to store connectivity. We have presented an un-structured approach in sufficient detail in [36] and released it as an open-sourcesoftware [39]. Here, we will make use of it to provide a performance perspectiveto the other two approaches that we presented above.3. AMR infrastructure methodsWhen numerically solving a PDE on an adaptively refined mesh, the AMRlibrary needs to provide certain mesh information, such as connectivity to neigh-boring elements, coordinates of the elements’ centroid and vertices. In thissection, we explain the methods to compute such information nimbly by manip-ulating the Morton code representation of an octree. The methods presented inthis section are common to both the hashed linear and binarized octree gener-ation methods.3.1. Adaptation levelWe need to know the current adaptation level of the element under consid-eration as well as the neighboring elements when constructing a 2:1 balancedoctree. The Morton code representation of the octree helps a great deal to ex-tract the level of adaptation without extra storage, because we only need tocheck the non-zero triple digit.As an example, let us again consider the node 23 in Fig. 4. The Mortoncode representation of node 23 is given in Fig. 3. On a 16-bit system, Mortoncode would be as shown in Fig. 7. Note that the Morton code is padded tothe right with zeros to create a 16-bit system representation. Starting from therightmost digit and scanning to the left, the last non-zero element occurs atlevel 3, therefore, the level of adaptation for this element is 3. Note that in thepresent study, element level 0 is the root element and we do not include it inthe bit representation.001,000,101,000,000,0Figure 7: Morton code for node 23 in Fig. 3 on a prototypical 16-bit system.For the special case of an element represented with all zeros as Mortoncode, we need to employ a trick. This special case of all zero bit elementbecomes important while enforcing the 2:1 balance. Because all digits are zeros,we cannot extract the adaptation level, but we can manipulate the siblings toextract this information. Note that the maximum bitset instance size is presetbefore the refinement as this is a template class where the size of the object is thetemplate parameter. One can start by flipping any one of the digits belongingto the maximum size. Note that by flipping the bit we are now operating onthe sibling now. After this we check for the existence of such a key. If it existsthen the element under consideration and the sibling belong to the same level.3.2. Edge lengthEdge length is required to calculate volume and surface integrals. Once theadaptation level is known following the above procedure, it is straightforwardto calculate the edge length resulting from an isotropic refinement in three-dimensional space, which is simply a generalization of one-dimensional bisectionprocedure.Given the edge length at the top level l0, the length of the segment at eachadaptation level i can be calculated as follows:li =l02i, i = 1, 2, ..., nLevel. (1)3.3. Centroid coordinatesIn the solution of fluid flow equations, whether a collocated or staggeredarrangement of the variables is pursued, the centroid of the element needs tobe known. Given the initial element’s centroid and length of the edges of theinitial element in x, y and z directions, it is possible to calculate the centroid ofthe current element as follows:xi = x0i − ((−1)bit[3j+i]12j+2)li, (2)j = 0, 1, ..., nLevel − 1,i = 1, 2, 3.where x0i implies the xi component of the centroid of the initial element. Notethat summation is implied on the j index. To illustrate with an example, thecoordinates of the centroid for the Morton code given in Fig. 4 are calculatedas follows:xc = x0 + l0(−14− 18+116),yc = y0 + l0(−14− 18− 116), (3)zc = z0 + l0(+14− 18+116).3.4. Element verticesWith an element’s centroid and edge length at hand, the vertices of eachelement is obtained by adding and subtracting the half of the edge length to thecentroid at each direction as follows:x = xc ± 0.5× l0,y = yc ± 0.5× l0, (4)z = zc ± 0.5× l0.3.5. Connectivity to neighborsConnectivity to surrounding elements is essential information when solvingPDEs, whether it is for calculating the flux of a quantity for enforcing conserva-tion principles or constructing the Jacobian (or stiffness) matrices for implicitschemes. The way connectivity information is stored and accessed in computermemory affects the overall performance tremendously. In this section, we dis-cuss the construction of connectivity that is unique to our implementation. Weconsider a two-dimensional example for ease of illustration, but our implemen-tation is for the three-dimensional case. Before proceeding further, we defineour terminology.3.5.1. SiblingsSiblings of each element, can be found by flipping the corresponding bit inx, y, z directions, respectively. In this operation, only one bit is flipped. Flippingthe bit in this context implies the modification of the bit to true if it is false,and to false when it is originally set to true. An example is shown in Fig. 8. Forthe element highlighted in red, the Morton code in two-dimensional space, usingthe procedure explained in previous sections, can be constructed as (10,01,01).The siblings for this element can be identified by flipping the corresponding bitin x and y directions. Therefore, the Morton code for the element located in theright hand side is (10,01,11) and for the element underneath is (10,01,00). Thesiblings are depicted in orange and yellow in Fig. 8, respectively. The elementhighlighted in cyan is a non-neighbor sibling.Figure 8: Illustration of neighbor identification from Morton code. Morton codes for differentcolors are: (00,11,11) for green; (10,01,01) for red; (01,10,00) for yellow; (10,01,11) for orange;(11,00) for magenta; (10,00,10) for cyan.3.5.2. Off-branch neighborsThe method we presented in section 2.1 for storing an octree as Mortoncodes enables the identification of off-branch neighbors by only manipulatingthe bits. Here, the word off-branch refers to the neighbors that are not siblingsto a given element. Finding non-local siblings in each direction, in this casethe green and purple element in Fig. 8, consists of two steps. In the first step,we need to detect the level at which the bit was flipped at the correspondingdirection. The technique to accomplish this task is presented in Algorithm 1.Algorithm 1 Detection of bit-flip leveldefine T : set of octree Morton codesgiven M ∈ Tdefine k = 1, 2, 3 : for (x, y, or z) directions, respectivelydefine d : dimension of the computational spacedefine f : flip-leveldefine Ml : level of MInitialize : bolean=M(d*l+k)for i = Ml − 1 : 0 doif bolean!=M(d*i+k) thenf = ibreak;end ifend forReturn fLet us describe it with an example. Consider the element highlighted in redin Fig. 8. Morton code for this element can easily be calculated as (10,01,01).Our goal is to construct the Morton code depicted in green without an explicitsearch. This requires displacement in the x−direction (i.e., k = 1 in 1). Thereare two steps involved in this task:• The location of the red element in Fig. 8 at the finest level is 01, so thebit in x−direction is 0 (or false). Moving up one level at the tree towardsparents, the location of this element is found as 01, which corresponds tofalse as well. Proceeding one level further, the bit corresponding to thatlocation becomes 1 (true). So this level is the bit-flip level. Therefore,moving from left to right, the bits corresponding to the x−direction is(1, 0, 0). So the level at which flip has occurred goes all the way up to thefirst level.• The second step is to use the information embedded in the bit representa-tion such that the code for the off-branch neighbor is directly extracted.To obtain the Morton code for off-branch neighbor (i.e., the green elementin Fig. 8) in x−direction, we flip all the bits in the x−direction up to thebit-flip level. This operation yields the Morton code as (00,10,11).Note that by following the above procedure, there is no need for an ex-plicit search to find the neighboring elements, because the Morton code for eachelement carries with itself those information in a compact fashion.Next, we explain two special cases in which the off-branch neighbor has adifferent level of adaptation. Because we work with a balanced tree, such specialcases can easily be resolved.000000 01 10 1100 01 10 110100 01 10 111000 0100 01 10 1110 111100 01 10 11Figure 9: Representation of the mesh given in Fig. 8 as a tree. Matching colors are used inboth figures.3.5.3. Case IThe off-branch neighbor has an adaptation level lower than the elementunder consideration. This is the case for the element highlighted in red in Fig.8. The goal is to construct the Morton code for the element highlighted inmagenta. Again all that needs to be done is to eliminate the last two bits sothat the level is reduced by one. This will lead to Morton code of (11,00) forthe magenta element in Fig. 8. Note that the length of each code is consistentwith the level of refinement.3.5.4. Case IIThe off-branch neighbor has an adaptation level higher than the elementunder consideration. To illustrate this case, consider the element highlightedin magenta. The goal here is to go in the reverse direction and construct theMorton code of the off-branch neighbors in the y-direction (i.e., the elementshighlighted in red and orange). This adds another step to what has alreadybeen done earlier,• The first step involves the prior two steps given in section 3.5.2. TheMorton code for the magenta element is given as (11,00), the flip hasoccurred at level 2, therefore flipping the bits of this code will result in(10,01),• Since the level of this element is higher than the element under consider-ation, it is required to repeat the last two bits, this will give (10,01,01),this element is one of the neighbors. To construct the second neighborall we have to do is find the sibling of this code in x-direction, which is(10,01,11).In three-dimensional space, two more neighbors are added to the list. Theseneighbors are constructed by alternating the bit in x− and z− directions.3.6. Detecting boundary elementsBoundary elements are those elements that only include siblings as part oftheir neighborhood connectivity. Boundary elements do not possess any off-branch neighbors. We make use of bitwise operations to identify boundaryelements for a given Morton code. The following theorem can be used to detectboundary elements from their Morton code representation.Theorem 3.1. An element is a boundary element at a given direction if andonly if its Morton code is free of bit-flips in that direction.Elements with Morton codes (10,11), (10,10) do not have any off-branchneighbors in the x-direction. Hence, they are boundary elements. If we traversein the y-direction, elements with codes (10,00) and (10,10) are also boundaryelements. The advantage of our approach is that no bits are used for taggingthe boundary elements since this criteria can efficiently be extracted from theMorton code at hand.4. Balancing the octreeDuring the process to generate an octree with dynamic adaptation, the 2:1balance may break. Such imbalance propagates, affecting immediate neighborsand neighbors of neighbors. This issue has been coined as the ripple effect inthe literature of balanced octree generation [16]. Our computer implementationled to the algorithm described in 2. It starts with an initial given queue, andextends that queue to maintain the balance for the entire tree.There are several ways to implement the NB /∈ Q section of the code. Werefer the reader to Algorithm 2 for definitions of NB and Q. The NB /∈ Qoperation can be prohibitively expensive if not handled in a proper way. In theprocess of balancing, if any of the neighbors of the element under considerationviolates the condition, before adding this element to the end of queue, one needsto check if that element is already in the refinement list. When using array datastructure or the Vector class, the cost of this operation can be calculated togrow as O(N ∗ N), where N is the number of elements tagged for refinement.The cost is prohibitively large for a mesh on the order of several hundred millionelements. To alleviate this problem and reduce the O(N ∗N) cost to m ∗O(N),where m is the number of iterations done by in the while loop due to the rippleeffect, we use a hash map for this list. We note here that for both the hashedlinear octree and the binarized octree generation methods, a hash table datastructure is used to store the list of the elements tagged for refinement. WeAlgorithm 2 2:1 balanced refinement-list generationdefine P = 1define Q : List of elements tagged to be refineddefine NB : Neighbordefine NBl : Level of neighbordefine Cl : Level of C elementdefine istart = 0; iend = Q.size()while P > 0 dodefine P = 0for i = istart : iend doC = Q(i)if NBl < Cl & NB /∈ Q thenappend NB to Qset P = 1end ifend foristart = iend; iend = Q.size()end whileReturn Qachieve significant reduction in time spent for the balance algorithm with thismodification. Timings for the Stanford bunny geometry using the Vector andMap in search are presented in Table 1. Once again, this example illustrates thevalue of selecting the right data structure in scientific computing.Table 1: Comparison of total runtime for two different data structures and their effect on thesearch algorithms in the 2:1 balance, the Stanford bunny geometry is used in this analysislevel 9 10 11 12 13mesh size 118,357 437,557 1,199,290 2,830,815 5,867,800time [s](Hash Map) 0.8 2.57 4.71 10.47 26.0time [s](Vector) 1.35 5.40 19.02 63.0 185.375. Tagging elements for refinementThe mesh generation and adaptation process requires the identification ofthe elements which intersects the immersed geometry. The brute-force approachto obtain the list of tagged elements for refinement would be through a search.However, a brute-force approach can become expensive for a large mesh. Witha brute-force search, the computational cost grows with m × n, where n is thenumber of points of the geometry and m number of mesh elements at the currentlevel. and ensure data-locality in the computer implementationOne way to alleviate this bottleneck and reduce the cost of search is tovoxelize the geometry and perform the search on the specific voxel instead ofthe entire geometry. This method improves search proportional to number ofvoxels. Alternatively, by taking advantage of the Morton code representation,it is possible to completely remove this search stage using geometry encoding.Removing the search here implies dropping the cost from O(m × n) to O(1) ifhash table used as an underlying container for the list of tagged elements.Coordinates of the centroid of the STL triangles are used to determine if anelement resides in a given voxel. Normals of the surfaces can also be assignedin a similar fashion.5.1. Voxelization of immersed geometryVoxelization of the immersed geometry comes handy to restrict certain nu-merical operations around the geometry. One could consider voxelization asa refined bounding box wrapping around the geometry. The map values arepointers to the geometry file portion that is segmented by the given element.We store the solid geometry in STereoLithography (STL) file format.The level of refinement in voxelization is much lower than the solution mesh.It is specified independent of the solution according to the solid geometry athand. The voxel class inherits from the base tree class as shown in Listing 1.We note that voxel need not be a 2:1 balanced octree.Figure 10: Voxelized representation of the Stanford bunny geometry.Voxelized view of the Stanford bunny is presented in Fig. 10. One doesnot need this many voxels to accelerate the search. A higher refinement level isadopted here for illustration purposes.Listing 1: Base tree class and inherited voxel classtemplate < s i z e t N, typename value>class Tree{template < s i z e t N1 , typename value>friend class Hdf5Xmf ; //wri te out in hdf5 formatprotected :r e a l r oo t l eng th [ 3 ] ; // leng th of the f i r s t generat ionr e a l roo tcoords [ 3 ] ; // centeroid of the f i r s t generat ionb i t s e t rootkey = 0 ; // root va lue i s always se t as 0uint npx ; // d i s c r i t i z a t i o n in x d i r ec t i onuint npy ; // d i s c r e t i z a t i o n in y d i r ec t i onuint npz ; // d i s c r e t i z a t i o n in y d i r ec t i on// s e l e c t e i t h e r unordered map or mapprivate :unordered map, va lue ∗ , hash over r ide>mesh ;map, va lue ∗ , compare overr ide>mesh ;unordered map, int> r e f i n e l i s t ;unordered map, int> d e r e f i n e l i s t ;public :Tree ( r e a l ∗ length , r e a l ∗ coords , u int nx , u int ny , u int nz ) ;}template < s i z e t N, typename value>class Voxel : public Tree{private :u int maxlevel ; //maximum l e v e l of refinementuint numMax; //number of elementsunordered map, va lue∗>mesh ;public :Voxel ( r e a l ∗ length , r e a l ∗ coords )// more member funct ions}The voxelization mesh has a similar structure as the solution mesh exceptthat it need not be balanced. Ancestor of each element can be found usingthe bitwise operations. Therefore, the box that includes the correspondinggeometric data can be identified quickly. In flow computations around complexgeometry, face normals of the solid geometry are also needed for the calculationof fluxes in finite-volume/discontinuous Galerkin methods or implementationof reconstruction schemes in the immersed boundary method. Voxelization isvery useful in retrieving such geometric information. If one is only interested inmesh generation without any solution, the more efficient method of geometricencoding using bisection search is recommended, which we discuss next.5.2. Geometry encodingEach point of the geometry may be encoded in a similar fashion as themesh elements. Given the coordinates of each point and the maximum level ofadaptation, it is possible to generate a Morton code for each point at the highestlevel of adaptation. In other words, the finest level element that withholds thepoint is detected and saved in the class tree. This is basically equivalent tousing the well-known bisection search in three-dimensional space and is onlyperformed once regardless of the maximum level of adaptation. This way thesearch stage for tagging elements is eliminated altogether. Geometry encoding isillustrated for a one-dimension case in Algorithm 3. The three-dimensional casecan be encoded by using bit interleaving in each directions. In this algorithm,we accumulate true or false for each time we bisect a line and generate a codefor a given point.Algorithm 3 Geometry encoding in one-dimensiondefine maxLevel : maximum desired level of refinement which is a presetvaluedefine Xmin : minimum x-coordinate of the computational domaindefine Xmax : maximum x-coordinate of the computational domaindefine xc : center coordinatedefine xp : x-coordinate of a given pointsdefine b : instance of std::bitsetfor i = 0 : maxLevel doxc = 0.5× (Xmax + Xmin)if xp < xc thenb(i) = falseXmax = xcelseb(i) = TrueXmin = xcend ifend forThe expense to encode the geometry is O(N), where N is the number oftriangular elements used to represent the solid geometry that we store as anSTL file. The clear advantages of geometry encoding are that it does not growas the solution mesh grows in size, and it is also cheaper in terms of time andmemory than generating another coarse tree for voxelization. Figure 11 showsthe improvement in runtime due to geometry encoding and compare the runtimeagainst the voxelization approach.6. ResultsWe present examples of 2:1 balanced octree generation for several STL ge-ometries. The test cases include: Stanford bunny, which is a well-known testcase from the computer graphics community; a sports car; a wind turbine; com-plex terrain of Mount Hood, OR. AMR for these geometries are presented inFig. 12, demonstrating the versatility of the current meshing capability.0 1 2 3 4 5 601020304050Million ElementsTotalruntime[s]Geometric EncodingVoxel level 6Figure 11: Effect of geometry encoding on overall execution time as a function of number ofCartesian mesh elements. Stanford bunny geometry with 88,754 surface triangles is used forthis analysisTable 2: Total mesh size for different adaptation levelslevel 10 12 14 16Bunny (no. of triangles = 88,754)No. of elements [millions] 0.43 2.83 10.52 23.96Sports car (no. of triangles = 844,804)No. of elements [millions] 0.11 1.84 16.18 74.096.1. Performance analysisIn this section, runtime and memory consumption of the hashed linear octreeand binarized octree generation methods are investigated. We also demonstratedeep-level adaptation for the wind turbine geometry using the binarized octreemethod. Deep-level implies refinement levels that are greater than 21. All thecalculations are carried out on a single compute-node of a Linux cluster with 64GB of DDR3 1600 MHz memory and an Intel Xeon E5-2670 2.60 GHz processor.The number of elements at different levels for the Stanford bunny and thesports car are presented in Table 2. These numbers are independent from thedata structure utilized in the AMR algorithm.Here, using the two Stanford bunny and a sports car geometries, the run-time for binarized octree and hashed linear octree are presented for levels 10,12, 14 and 16 in Tables 3 and 4, respectively. The percentage of slow-down due(a) Stanford bunny (b) Sports car(c) Wind turbine (d) Mt. Hood, ORFigure 12: Illustration of the 2:1 balanced Cartesian adaptive mesh generation for complexgeometries from different fields.to the use of a hash table is also presented for further clarification. For bothof the geometries, the Z-compliant red-black tree is superior in terms of totalruntime.Table 3: Comparison of total run-time. The Stanford bunny geometry is used for this analysislevel 10 12 14 16Bunny (No. of triangles = 88,754)Z-curve compliant Binarized octree [seconds] 3.3 21.14 69.75 145.87Z-curve enhanced Hashed linear octree [seconds] 3.7 26.85 87.85 180.01Memory consumption for the binarized octree and the hashed linear octreeare presented for levels 12, 10, 14 and 16 in Tables 5 and 6 for the same geome-tries discussed in the previous section. Valgrind is used for memory profiling[40].According to the results, hashed linear octree generation method use lessmemory than the binarized octree method, but the difference in memory usageTable 4: Comparison of total run-time. The sports car geometry is used for this analysislevel 10 12 14 16Car (No. of triangles = 844,804)Z-curve compliant binarized octree [seconds] 1.58 19.68 164.56 696.54Z-curve enhanced hashed linear octree [seconds] 1.64 24.73 260.88 753.58is not substantial enough to cause any concern.Table 5: Comparison of memory consumption. The Stanford bunny geometry is used for thisanalysis.level 10 12 14 16Bunny (Ntriangles = 88,754)Z-curve compliant binarized octree [MB] 61.4 281.8 792.9 1500Z-curve enhanced hashed linear octree [MB] 59.6 261.7 675.5 1300extra memory usage [%] 3.0 7.6 17.4 15.4Table 6: Comparison of memory consumption. The sports car geometry is used for thisanalysis.level 10 12 14 16Car (No. of triangles = 844,804)Z-curve compliant binarized octree [MB] 106.3 341.1 1800 5800Z-curve enhanced hashed linear octree [MB] 106.3 305.7 2000 5100extra memory usage [%] 0.0 11.6 -11.1 13.7Unstructured mesh approaches use much more memory than their tree-basedAMR counterparts. The total run-time using puAMR [36] on a single processoris compared against the performance of the binarized octree generation methodin Table 7. In this comparison, we are not able to go beyond 16 levels of adap-tion, because the unstructured approach runs out of computer memory. Theother important factor to mention here is that in the unstructured approach,the search time for tagging elements for refinement is much longer than the bi-narized octree generation that use geometric encoding for the search algorithm.The slow performance of the unstructured approach is expected and we includeit in this study to highlight the superior performance of tree-based approachesfor AMR.Table 7: Comparison of runtime for binarized octree against the graph-based unstructuredmesh approach. Stanford bunny geometry is used for this analysis.level 7 8 9 10 11mesh size 8163 31242 118357 437,557 1,199,290Z-curve compliantBinarized octree[seconds] 0.16 0.21 0.438 2.57 4.71Graph-based[seconds] 3.16 4.71 9.02 32.69 615.12Table 8 demonstrates the deep-level adaption that is inherently feasible withthe binarized octree generation approach. We use the wind turbine geometrywith 70,406 surface triangles. The mesh around this geometry is presented inFig. 12(c).Table 8: Number of mesh elements as a function of adaptation level for a binarized octreegeneration method. Wind turbine geometry with 70,406 surface triangles is used. Note thathashed linear octree method limits the maximum level of adaption to 21 on 64-bit systems.level 20 21 22 23mesh size 1,246,330 2,925,973 6,007,583 10,714,0127. ConclusionsThere are various algorithms for adaptive mesh refinement. Performance ofthese different approaches can vary greatly. In the present work, we proposedthe binarized octree generation method by combining Morton encoding witha red-black tree, which is special kind of self balancing binary tree. In ourcomputer implementation, we benefited from the C++ STL bitset Set/Mapclasses. The C++ STL bitset enables us to use an arbitrary length for the bitrepresentation regardless of the hardware architecture. Additionally, the red-black tree operates on comparison of bits only. Because of these two features, weare able generate meshes with arbitrary levels of adaptation on a single proces-sor. On the contrary, the hashed linear octree generation method is restrictedto 21 levels for 64-bit systems, beyond which an overflow is inevitable becauseof conversion to an integer number in the hash function. However, one can relaxthe Z-order curve compliance and generate arbitrary orders of adaptation at theexpense of data-locality and loss of runtime performance.Another advantage of using the red-black tree in octree generation is that itkeeps its elements sorted by default and preserves the Z-order curve exactly. Inthe hashed linear octree generation method, the order can get broken because ofcollisions in the hash table. An important implication of this strict compliancewith the Z-order curve is that the generated mesh is ready to be partitioned forparallel computations without any need for sorting its elements.We exploited the bitwise representation of an octree and developed tech-niques to extract neighborhood connectivity from the Morton codes withoutexplicitly storing the connectivity information. Additionally, we developed ageometry encoding technique to rapidly tag the elements for refinement.Our results show that the binarized octree technique is faster than the hashedlinear octrree technique at the expense of a slight increase in memory usage. Thecurrent work is released as open-source code which includes both the binarizedand the hashed linear octree generation methods [34].AcknowledgmentsThis material is based upon work supported by the National Science Foun-dation under Grant No. 1440638 and Army Research Office under Grant No.W911NF-17-1-0564.References[1] I. Babushka, J. Chandra, J. E. Flaherty, Adaptive computational methodsfor partial differential equations, vol. 16, SIAM, 1983.[2] M. J. Berger, J. Oliger, Adaptive mesh refinement for hyperbolic partialdifferential equations, J. Comput. Phys. 53 (3) (1984) 484–512.[3] M. Berger, P. Colella, Local adaptive mesh refinement for shock hydrody-namics, J. Comput. Phys. 82 (1) (1989) 64 – 84.[4] BoxLib, https://ccse.lbl.gov/BoxLib/, 2017.[5] P. MacNeice, K. M. Olson, C. Mobarry, R. De Fainchtein, C. Packer,PARAMESH: A parallel adaptive mesh refinement community toolkit,Comput. Phys. Commun. 126 (3) (2000) 330–354.[6] BoxLib, https://ccse.lbl.gov/AMReX/, 2017.[7] M. Adams, Chombo software package for adaptive mesh refinement appli-cations design document .[8] G. L. Bryan, M. L. Norman, B. W. O’Shea, T. Abel, J. H. Wise, M. J.Turk, D. R. Reynolds, D. C. Collins, P. Wang, S. W. Skillman, et al.,Enzo: An adaptive mesh refinement code for astrophysics, Astrophys. J.211 (2) (2014) 19.[9] SAMRAI, https://computation.llnl.gov/projects/samrai,2017.[10] R. D. Hornung, S. R. Kohn, Managing application complexity in the SAM-RAI object-oriented framework, Concurr. Comput. 14 (5).[11] M. J. Berger, R. J. LeVeque, Adaptive Mesh Refinement using Wave-Propagation Algorithms for Hyperbolic Systems, SIAM J. Numer. Anal.35 (1998) 2298–2316.[12] T. Lewiner, V. Mello, A. Peixoto, S. Pesco, H. Lopes, Fast generation ofpointerless octree duals, in: Comput. Graph. Forum, vol. 29, 1661–1669,2010.[13] A. Dubey, A. Almgren, J. Bell, M. Berzins, S. Brandt, G. Bryan, P. Colella,D. Graves, M. Lijewski, F. Lﬄer, B. OShea, E. Schnetter, B. V. Straalen,K. Weide, A survey of high level frameworks in block-structured adaptivemesh refinement packages, J. Paral. Distrib. Computing 74 (12) (2014) 3217– 3227.[14] D. Meagher, Geometric modeling using octree encoding, Comput. Gr. Im-age Process. 19 (2) (1982) 129–147.[15] H. Samet, The quadtree and related hierarchical data structures, ACMComput. Surv. 16 (2) (1984) 187–260.[16] H. Sundar, R. S. Sampath, G. Biros, Bottom-up construction and 2: 1balance refinement of linear octrees in parallel, SIAM J. Sci. Comput. 30 (5)(2008) 2675–2708.[17] I. Gargantini, An Effective Way to Represent Quadtrees, Commun. ACM25 (12) (1982) 905–910.[18] I. Gargantini, Linear octtrees for fast processing of three-dimensional ob-jects, Comput. Gr. Image Process. 20 (4) (1982) 365–374.[19] C. Burstedde, L. C. Wilcox, O. Ghattas, p4est: Scalable algorithms for par-allel adaptive mesh refinement on forests of octrees, SIAM J. Sci. Comput.33 (3) (2011) 1103–1133.[20] T. Isaac, C. Burstedde, O. Ghattas, Low-cost parallel algorithms for 2: 1octree balance, in: Parallel & Distributed Processing Symposium (IPDPS),2012 IEEE 26th International, IEEE, 426–437, 2012.[21] N. Ishikawa, D. Sasaki, K. Nakahashi, Large-scale distributed computationusing building-cube method, in: 49th AIAA Aerospace Sciences Meeting,754, 2011.[22] T. Ishida, S. Takahashi, K. Nakahashi, Efficient and robust cartesian meshgeneration for building-cube method, Journal of Computational Scienceand Technology 2 (4) (2008) 435–446.[23] T. Tu, D. R. Ohallaron, Balance refinement of massive linear octrees, Tech.Rep., Technical Report CMU-CS-04-129, Carnegie Mellon School of Com-puter Science, 2004.[24] S. F. Frisken, R. N. Perry, Simple and efficient traversal methods forquadtrees and octrees, J. Graph. Tools 7 (3) (2002) 1–11.[25] M. Schwarz, H.-P. Seidel, Fast parallel surface and solid voxelization onGPUs, in: ACM Trans. Graph., vol. 29, ACM, 179, 2010.[26] J. Baert, A. Lagae, P. Dutre´, Out-of-core construction of sparse voxel oc-trees, in: Comput. Graph. Forum, vol. 33, Wiley Online Library, 220–227,2014.[27] Tero Karras, Thinking Parallel III, https://devblogs.nvidia.com/parallelforall/thinking-parallel-part-iii-tree-construction-gpu/,2012.[28] P. M. Campbell, K. D. Devine, J. E. Flaherty, L. G. Gervasio, J. D. Teresco,Dynamic octree load balancing using space-filling curves, Williams CollegeDepartment of Computer Science, Tech. Rep. CS-03-01 .[29] J. E. Flaherty, R. M. Loy, M. S. Shephard, B. K. Szymanski, J. D. Teresco,L. H. Ziantz, Adaptive local refinement with octree load balancing for theparallel solution of three-dimensional conservation laws, J. Paral. Distrib.Comput. 47 (2) (1997) 139–152.[30] P. Brunet, I. Navazo, Solid representation and operation using extendedoctrees, ACM Trans. Graph. 9 (2) (1990) 170–197.[31] H. T. Vo, C. T. Silva, L. F. Scheidegger, V. Pascucci, Simple and efficientmesh layout with space-filling curves, J. Graph. Tools 16 (1) (2012) 25–39.[32] G. M. Morton, A computer oriented geodetic data base and a new techniquein file sequencing, International Business Machines Company, 1966.[33] J. J. Hasbestan, I. Senocak, A short note on the use of the redblack treein Cartesian adaptive mesh refinement algorithms, J. Comput. Phys. 351(2017) 473 – 477.[34] J. J. Hasbestan, I. Senocak, rebl-AMR: Red-black tree based adaptive meshrefinement for immersed-boundary simulations, https://github.com/GEM3D/Hash_vs_RB, version 0.1, DOI = 10.5281/zenodo.836102, 2017.[35] A. Drozdek, Data structures and algorithms in C++, Cengage Learning,2012.[36] J. J. Hasbestan, I. Senocak, A parallel adaptive mesh refinement softwarefor complex geometry flow simulations, in: 23rd AIAA ComputationalFluid Dynamics Conference, 3301, 2017.[37] L. J. Guibas, R. Sedgewick, A dichromatic framework for balanced trees,in: Foundations of Computer Science, 1978., 19th Annual Symposium on,IEEE, 8–21, 1978.[38] B. Schneier, Red-black trees, Dr. Dobb’s Journal 17 (4) (1992) 42–46.[39] J. J. Hasbestan and I. Senocak , puAMR: Parallel adaptive mesh refinement(AMR) with unstructured connectivity, https://github.com/GEM3D/AMR/puAMR, version 1.0.0, DOI = 10.5281/zenodo.571216, 2017.[40] N. Nethercote, J. Seward, How to shadow every byte of memory used bya program, in: Proceedings of the 3rd International Conference on VirtualExecution Environments, ACM, 65–74, 2007.

GEM3D/Hash_vs_RB: rbAMR

Abstract

Extracted data

GEM3D/Hash_vs_RB: rbAMR

Abstract

Extracted data

Related items

Related items