Could you use the Barnes-Hut algorithm iteratively-- with multiple center quadrants? I was wondering if you could use Barnes-Hut simulation beyond what it was originally intended to be. For many Barnes-Hut algorithms, the forces are only considered for a single quadrant, the centroid, or the stellar body. Then, the algorithm branches out from there, affecting areas of influence and quadrants recursively. For instance:

Seems like the above Barnes-Hut algorithm was based on the central body from the animation.
My Question:
Would performing Barnes-Hut iteratively across all bodies, treating each body in-turn as the centroid, result in an accurate representation of an n-body problem where the sum force of gravity of all bodies is considered? Or am I misunderstanding exactly what the Barnes-Hut algorithm is?

If I'm misunderstanding the algorithm, can somebody re-explain exactly how this algorithm works? For anyone who understands programming to some degree, could anyone look at this project and tell me if I'm missing something huge here? It's a Java GitHub implementation of the Barnes-Hut algorithm, but I've iterated it across all bodies (which may be incredibly stupid). Also-- yes, I know that's not how time works. Note: Credit due to original professor, as noted on GitHub.
Also, for those who aren't tech-savvy, can you look at this GIF and see anything inherently wrong? Red is less mass, white is more mass; yellow are two or more collided masses. Once the third yellow dot (combined mass) appears, things get interesting. I can't tell if interesting good, or interesting... bad.

 A: It's reasonably well known that at large enough distances, the force between a cluster of objects and another object is essentially equivalent to force between the center-of-mass of the cluster and the other object:
$$
F_i\sim\sum_j\frac{m_im_j}{\vert r_i-r_j\vert^2}\simeq \frac{m_im_\text{tot}}{\vert r_i-r_\text{com}\vert^2}
$$
where $m_\text{tot}=\sum_jm_j$ and $r_\text{com}$ is the center of mass of the $j$ particles.
The Barnes-Hut algorithm applies this to $n$-body simulations by breaking the domain into reasonably sized clusters and thereby reducing the number of effective bodies to iterate over, rather than summing over $n$ particles $n-1$ times. From what I can see/tell, the $n$-body simulations using the BH algorithm are fairly accurate, though there are better algorithms in terms of speed (e.g., fast multipole method or other hybrid schemes).

Would performing Barnes-Hut iteratively across all bodies, treating each body in-turn as the centroid, result in an accurate representation of an n-body problem where the sum force of gravity of all bodies is considered? Or am I misunderstanding exactly what the Barnes-Hut algorithm is?

Your proposal is essentially wanting to revert back to the summation of $n$ particles $n-1$ times, which would be equivalent to not at all using BH. It'll be accurate, though probably not much more accurate than BH, and it most certainly would be really slow. 
So by questioning whether one can iteratively apply BH to all bodies, I suspect that you are misunderstanding the BH algorithm. The whole point of it is to reduce the workload required in finding the forces on the particles at each time step. You don't really want to increase the number of clusters, you want that minimized because doing an $\mathcal{O}\left(n^2\right)$ update each time step is really, really slow. 
