an octree is a tree data 
 structure 
 that recursively and adaptively subdivides the 3d space into octants , 
 and is often used as a container for cartesian space data . here , we 
 use octrees to store all atomic centers ( or 3d atomic cartesian coordinates ) 
 of a molecular system , organized recursively based on their spatial 
 locality .
octrees have been used in the past for fast multipole - like 
 approximations of force fields in particle simulations and in molecular solvation energy and force calculations . in this work , we show that octrees can 
 be used for space - efficient 
 maintenance of nonbonded neighborhood lists ( i.e. , nblists ) of a molecular 
 system and can be updated very efficiently as the configuration of 
 the system changes .
a special type of octree , as defined below , will 
 be used for the purpose .
octrees are more space - efficient , update - efficient , 
 and cache - efficient in implicit maintenance of nblists compared to 
 explicitly maintaining them . also , unlike explicit nblists
, a single 
 octree can be used for any distance cutoff required for evaluating 
 pairwise interactions among atoms
. an octree   is called 
 ( , )-admissible 
 provided no leaf of   contains more 
 than  points and 
 each internal node has more 
 than / points , 
 where   > 0 is 
 an integer and   
 1 .
figure 1 shows an example of a ( 3,1)-admissible 
 quadtree , which is a 2d variant of octrees
. 
 the quadtree is constructed by recursive subdivision of the box bounding 
 the initial set of points ( i.e. , circle centers ) into nonempty quadrants . 

observe that in the resulting tree each internal node contains more 
 than 3/1 = 3 points , but no leaf contains more than 3  1 = 3 .
given a set p of points ( or atom centers ) in three dimensions , representing 
 a simulated molecular system , consisting of arbitrary number of chemically 
 connected components , a ( , )-admissible 
 octree can be constructed 
 by first finding a box ( or cube ) bounding p and then 
 recursively subdividing the initial cube into smaller nonempty subcubes 
 until each subcube encloses at most  points . 

it can be shown that such a tree can be constructed in   ( n 
 log n ) time , where n = |p| ( see supporting 
 information ) . in 
 this section ,
we show 
 that a ( , )-admissible 
 octree can be modified 
 easily to use space linear in the number of atoms in the molecular 
 system it stores .
the trick is to use a contracted octree , which is obtained from a standard octree   by directly 
 connecting each node of   to its nearest 
 proper ancestor with more 
 than one child .
thus , each internal node of a contracted octree 
 has at least two children .
all algorithms we have described so far 
 can be easily modified to work on contracted octrees without any asymptotic 
 increase in running times .
it can be shown that a contracted ( , )-admissible 
 octree   storing a 
 molecular system containing n atoms uses (n ) space , where   is a positive 
 integer and   
 1 ( as shown in the supporting information ) . given two 
 octrees a and b storing molecules a and b ( possibly a = b and 
 thus a = b ) , respectively , 
 the pairwise interactions between atoms of a and b can be computed as follows .
suppose we are only interested 
 in computing interactions between a pair of atoms provided their centers 
 lie within a given distance cutoff d. we perform 
 a simultaneous recursive traversal of a and b starting from 
 their root nodes .
suppose at some point we are at node u of a and node v of b. if the two 
 cubes corresponding to u and v are 
 separated by a distance ( surface - to - surface ) larger than d , then there can be no interaction . otherwise ,
if both nodes are 
 leaves , we consider each pair p  u , q  v , 
 and compute the interaction between p and q provided the two atom centers lie within distance d. if neither of the two conditions above hold , we move 
 to the children of u and/or v and 
 compute the interactions recursively . supporting 
 information figure 2
it can be shown that all pairwise interactions between atoms within 
 distance cutoff d can be computed in   ( nd ( d + ) ) time , where  is 
 the time to compute one pairwise interaction ( see supporting information ) .
for any fixed value of d the running time of interaction calculation increases with the increase 
 of .
this happens 
 because the larger the value 
 of , the larger 
 the region of unnecessary exploration . 

however , as computing interactions is a recursive function , the smaller 
 the value of , the deeper 
 the recursion and the larger 
 the recursion overhead .
hence , one must chose a value of   that balances 
 between these two types of 
 overheads .
the key 
 novelty of our work 
 is that we allow fast updates of the octree under conformational changes 
 of the molecular system , during the course of simulation .
figure 2 illustrates the basic idea behind such updates 
 using the ( 3,1)-admissible quadtree from 1 .
figure 2a shows the original state of the quadtree with the point 
 ( i.e. , circle center ) to be updated marked with solid black outline . 

the updated position of that point is shown in figure 2b but without updating the structure of the quadtree .
observe 
 that the point has moved out of leaf u and internal 
 node x ( i.e. , parent of u ) but remains 
 inside node w ( i.e. , parent of x ) and all its ancestors .
hence , in order to keep the structure of 
 the quadtree correct , the point must be removed from nodes u and x. figure 2c shows the quadtree after the removals .
observe that since leaf u is now empty it no longer appears in the quadtree and 
 since the number of points inside x is less than 
 4 it has now become a leaf .
also observe that the updated point has 
 now moved to the top - right quadrant of node w , which 
 was originally empty .
so , the new children ( leaf v ) of w are created corresponding to that quadrant 
 containing the updated point .
detailed description 
 of the method together with pseudocode is given in the supporting information .
updating the ( 3,1)-admissible 
 quadtree from 1 : ( a ) the point to 
 be updated is marked with solid black outline .
( b ) the point is moved 
 to a new location , and so it moves out of leaf u and 
 internal node x but remains inside node w and all its ancestors .
( c ) node u becomes empty , 
 and so , it is removed .
node x now contains only 3 
 points , and so , it becomes a leaf .
the updated point moves to the 
 top - right quadrant of node w , which was empty before , 
 and so , a new leaf v containing only the updated 
 point is created .
we prove in the supporting information that a ( , 2)-admissible 
 octree storing a molecular 
 system of n atoms supports updates of atomic positions 
 in   ( log n ) amortized time 
 each , where   is a positive 
 integer .
however , observe 
 that in molecular dynamics simulations atoms often move very short 
 distances at every step , and the number of nodes on the shortest path 
 between the leaves corresponding to the old and the new centers of 
 the atom is very small .
the update time will often be proportional 
 to the length of this path , and thus much better than   ( log n ) .
we reimplemented 
 the bycc method , and incorporated it into our molecular simulation 
 library after a rigorous testing and validation phase that verified 
 compliance with the charmm implementation . in order to avoid the costly 
 update of nblists after each change of atomic position , we update 
 them only when at least one atom has moved by more than 0.5  
 since the last update . the space requirements of explicit nblists 
 and octrees are compared in section 2.1 .
section 2.3 reports the running times and cache performance 
 of a limited - memory bfgs minimizer using 
 the two data structures for maintaining molecular systems under conformational 
 changes . all algorithms 
 were implemented in c , and tested on intel xeon 5680 , using the icc 
 compiler with optimization flag 
we ran 
 our experiments on the following six proteins of varying sizes downloaded 
 from the protein data bank , and prepared using 
 charmm c34b2:1clv ( 4680 atoms ) , 2hmi ( 13  758 atoms ) , 1suv ( 25  118 atoms ) , 2xrp ( 33  602 atoms ) , 3lue ( 46  200 atoms ) , 
 and 2cgt ( 70  847 
 atoms ) .
we have found that ( 60 , 2 ) octrees work best for molecular 
 simulation 
 application . in our experiments
, we compare the cost of maintaining 
 molecules 
 under conformatonal changes using nblists and octrees . for this purpose
, 
 we use these data structures to store molecules that undergo changes 
 in atomic positions during energy optimization using the limited - memory 
 bfgs ( lbfgs ) algorithm .
lbfgs is a quasi - newton 
 optimization method that uses a limited memory variation of the broyden  fletcher  goldfarb  shanno 
 ( bfgs ) algorithm for solving nonlinear optimization problems .
we have used our own 
 c port of jorge nocedal s original fortran implementation of 
 the lbfgs algorithm .
plot a compares the sizes 
 of of nblists and octrees constructed 
 for six protonated proteins containing 460071  000 atoms . 

plot 
 b zooms in on the octree curve in plot a. for each protein included 
 in our experiments , we constructed nblists for each integral distance 
 cutoff ranging from 9 to 20  . figure 3 compares the sizes of these data structures . in figure 3a , we plot the space used by the octrees as well 
 as by the nblists for 12 , 15 , and 18  distance cutoff .
figure 3b zooms in on the octree curve in figure 3a . as evident from the plots both octrees and nblists 
 for any given distance cutoff use space linear in the number of atoms 
 in the molecular system .
however , size of an nblist grows as the distance 
 cutoff increases . also , the size of an octree is much smaller than 
 that of an nblist for the same system .
for the six proteins in the 
 plots , octree size varied from 55 kb ( for 1clv ) to 1 mb ( for 2cgt ) , 
 while nblist size varied from 3.5 mb ( for 1clv ) to around 50 mb ( for 
 2cgt ) for 12  distance cutoff , and from 9 mb to around 130 mb 
 for 18  cutoff .
overall , nblists used from 45 to 65 times more 
 space than octrees for 12  cutoff , and the range increased to 
 120165 for 18  .
it turns out that for the proteins 
 and the d values 
 included in our experiments 3.7nd is 
 a lower bound on the sizes of nblists while 1.7nd is an upper bound . in comparison , octrees use only 
  ( n ) space .
as can be seen , octrees are orders 
 of magnitude more memory efficient than nblists .
figure 4 compares the time required to complete 
 100 steps of the lbfgs minimizer using octrees and nblists as the 
 distance cutoff for computing pairwise interactions grows . for smaller 
 distance cutoffs the overhead of recursion in calculation of interactions 
 slows down octree - based computation .
however , that overhead gradually 
 diminishes compared to the cost of actual computation as the cutoff 
 value increases .
also , nblists constructed for large cutoff values 
 are often too large to fit in the cache , and as a result , computations 
 slow down due to costly cache misses . octrees , on the other hand , 
 incur very few cache misses as they are often small enough to fit 
 into the cache .
figure 4a plots the ratio of the running 
 time of the nblists - based lbfgs minimizer to that of the octree - based 
 minimizer when minimizing van der waals ( vdw ) energy only .
octrees 
 were up to 20% slower than nblists for distance cutoffs 9 11 
  , comparable to nblists for cutoff values 12 and 13  and 
 ran faster than nblists for all systems when the cutoff value reached 
 14  .
supporting information figure 6 reveals one of the major reasons behind the improving performance 
 of octrees . as the distance cutoff increased
, the sizes of nblists 
 also increased , leading to increasing number of misses at various 
 levels of the caches . on the other hand , much smaller sizes of octrees 
 and
their cache - friendly recursive traversal for energy evaluation 
 meant that the number of cache misses incurred by the octree - based 
 minimizer grew at a much slower rate .
time required for 100 steps of lbfgs minimization 
 involving nblists 
 and octrees : ( a ) minimize in ( and hence compute ) vdw energy for various 
 distance cutoffs , where both data structures use the same cutoff value ; 
 ( b ) minimize in vdw plus hbond energy , where nblists use the same 
 cutoff value ( given in the horizontal axis ) for both types of energy 
 while octrees use a fixed 3  cutoff for hbond and the same cutoff 
 as nblists for vdw ( given in the horizontal axis ) ; ( c ) minimize in 
 vdw plus hbond plus coul , where nblists use the same cutoff value 
 ( given in the horizontal axis ) for all types of energy while octrees 
 use a fixed 3  cutoff for hbond , a fixed 9  cutoff for 
 vdw , and the same cutoff as nblists for coul ( given in the horizontal 
 axis ) .
figure 4b is similar to figure 4a but uses hydrogen 
 bonding ( hbond ) energy in addition 
 to vdw energy .
the nblists - based minimizer evaluated the hbond energy 
 using the same nblist constructed for evaluating vdw energy .
however , 
 since the same octree works for all distance cutoffs , the octree - based 
 minimizer always evaluated the hbond energy using a 3  distance 
 cutoff irrespective of what cutoff value was used for vdw energy . 
 because of the lower overhead in evaluating hbond energy , octrees 
 achieved even better speedup with respect to nblists than in figure 4a . for example , for 16  cutoff , octrees were 
 1030% faster than nblists in figure 4a , but in figure 4b , the speedup values improved 
 to 2040% .
figure 4c is for minimizing 
 the sum of vdw , 
 hbond , and coulomb electrostatic ( coul ) energy .
while the nblists - based 
 minimizer used the same distance cutoff for evaluating all three types 
 of energy , the octree - based minimizer always used a 3  cutoff 
 for hbond , and 9  for vdw .
octrees always outperformed nblists 
 and ran 80100% faster ( i.e. , 1.82.0 times faster ) 
 for 16  cutoff ( see figure 4c ) .
figure 5 shows a breakdown of the total 
 time spent in minimizing 2cgt into the time required for updating 
 the data structure and the time needed for computing the energy .
observe that octree 
 update times are almost insensitive to distance cutoffs while time 
 required for updating nblists increases with the cutoff value . even 
 for a 9  cutoff
, octrees can be updated more than 20 times faster 
 than nblists , and for 20  , updating octrees can be more than 
 100 times faster .
though octree - based energy computation is slower 
 than nblists for smaller distance cutoffs , this relative speed improves 
 and octrees start to beat nblists for large cutoffs .
plots show a breakdown 
 of the total time required for running 100 
 steps of lbfgs minimization on 2cgt into update time and energy computation 
 time .
please see 4 for the cutoff values used in parts a , b , and 
 c.
we 
 have developed memory and cache efficient octree implementation 
 of implicit nblists , which has almost 2 orders of magnitude lower 
 memory usage as compared to existing implementations .
in addition 
 to tremendous reduction in memory requirements the code runs approximately 
 1.5 faster than nblists in practical use case scenarios .
the algorithm 
 described in the paper is suitable for effective parallelization and 
 vectorization , as described and demonstrated in the supporting information section 2.2 , and thus can be applied 
 to molecular dynamics ( md ) .
current implementation is directly useful 
 for simulations requiring multiple simultaneous minimizations , such 
 as multistart monte carlo minimization ( mcm ) .
 molecular 
 mechanics and dynamics simulations use distance based 
 cutoff approximations for faster computation of pairwise van der waals 
 and electrostatic energy terms . 
 these approximations traditionally 
 use a precalculated and periodically updated list of interacting atom 
 pairs , known as the  nonbonded neighborhood lists  or 
 nblists , in order to reduce the overhead of finding atom pairs that 
 are within distance cutoff . 
 the size of nblists grows linearly with 
 the number of atoms in the system and superlinearly with the distance 
 cutoff , and as a result , they require significant amount of memory 
 for large molecular systems . 
 the high space usage leads to poor cache 
 performance , which slows computation for large distance cutoffs . 
 also , 
 the high cost of updates means that one can not afford to keep the 
 data structure always synchronized with the configuration of the molecules 
 when efficiency is at stake . 
 we propose a dynamic octree data structure 
 for implicit maintenance of nblists using space linear in the number 
 of atoms but independent of the distance cutoff . 
 the list can be updated 
 very efficiently as the coordinates of atoms change during the simulation . 
 
 unlike explicit nblists , a single octree works for all distance cutoffs . 
 
 in addition , octree is a cache - friendly data structure , and hence , 
 it is less prone to cache miss slowdowns on modern memory hierarchies 
 than nblists . 
 octrees use almost 2 orders of magnitude less memory , 
 which is crucial for simulation of large systems , and while they are 
 comparable in performance to nblists when the distance cutoff is small , 
 they outperform nblists for larger systems and large cutoffs . 
 our 
 tests show that octree implementation is approximately 1.5 times faster 
 in practical use case scenarios as compared to nblists . 
Octrees for Maintaining Molecules
under Conformational
Changes
Results and Discussions
Conclusion