Possibly remove MeshBlockTree

[lroberts36]

As I am working on including geometric multigrid in Parthenon (#911), I am finding more and more that I would like to move away from relying on MeshBlockTree. I could possibly remove MeshBlockTree entirely in that PR, but I wanted to get some feedback before I went down that (possibly drastic) path.

The reasons that I can see for doing this are:

1. The code would be substantially simplified - If we move to storing the LogicalLocations of the mesh in standard library containers (i.e. either a std::map<LogicalLocation, ...> or a std::unordered_map<LogicalLocation, ...>), then all of the MeshBlockTree code and most of the code in bvals_base.cpp is no longer necessary. Generally, I have found the code in bvals to be some of the hardest to follow within Parthenon (see Refactor bvals code #903). With the recent changes to LogicalLocation, there are methods to find lists of possible neighbor blocks for a given block, determine the neighboring block pairs offsets, etc. which are required for building a NeighborBlock object, and other mesh location related utilities. As a result, it is very straightforward to find all of the neighbors of a given block by just finding a list of possible neighbors, and then seeing which of those neighbors are in the mesh.
2. It allows easy generalization to including multi-grid coarsened meshes in a Mesh - Storing the levels of the mesh in a vector of std::maps makes it very easy to find parent blocks on the coarser mesh and daughter blocks on a finer mesh, as well as find neighbor blocks on the current mesh. This could of course be done if MeshBlockTree was extended, but this seems somewhat more complex to me.
3. It may make it easier to move to having a distributed mesh structure - Currently, the structure of the mesh (i.e. the MeshBlockTree) is duplicated on all ranks. For very large numbers of ranks, this can result in the mesh structure taking up a significant fraction of each nodes memory. As a result, we may eventually want to move to having the mesh structure not be entirely duplicated on each rank. It seems easier to me to do this if the local mesh structure is stored in a std::map than if it is stored in a MeshBlockTree (although I admittedly haven't spent a huge amount of time thinking about how to do this).

Possible drawbacks:

1. The code is working fine using MeshBlockTree currently.
2. I am not certain of the performance implications of moving away from a MeshBlockTree (although it seems likely that we can make a std::unordered_map perform faster than the current octree structure).
3. The swarms code relies on inheriting from BoundaryBase and BoundaryCommunication, so there will be some refactoring that has to be done there.

[Yurlungur]

IMO if there's good reason to refactor this code, this seems perfectly reasonable to me. A tree datastructure doesn't have to be represented with pointer chasing, and your space filling curve/morton number code already is the biggest lift for moving us off of pointers.

I might worry about using actual std::unordered_map, at least with standard hashing functionality, though. I think we have an essentially perfect hash in the morton number, so my preference would be to implement our own hashing (which std::unordered_map lets one do).

[lroberts36]

There is already a hashing overload built for LogicalLocation, since this type is not default hashable (see the bottom of logical_location.hpp). That being said, for really deep AMR structures we don't have a perfect hash with the Morton number since it takes up more bytes than a uint64_t which is what the hash function must return. This is only an issue for more than 21 levels of refinement though. There is a TODO comment in the code musing about how to deal with this though.

[Yurlungur]

If we're that deep, I think we have other problems, so not disturbing IMO

[pgrete]

Interesting idea!

Here are my 50 cents on the items raised

• I'm generally in favor of refactoring if it simplifies/cleans up code. If it also enables/simplifies other pieces or new methods then I'm also happy to sacrifice a couple percent in performance.
• Regarding the duplication of the MeshBlockTree, is this actually a concern right now? This was a concern in the original Enzo, which motivated the forrest-of-octree approach in Enzo-E. So far I considered our/Athena++ current approach as already on the lightweight side wrt to mesh management.
• 21+ levels are not out of the ordinary especially when it comes to zoom-in cosmological simulations. So we should definitely keep this in mind (though it'd also required adaptive timestepping to be practically usable).
• Given that the swarm code (to my knowledge) has not been tested a large scale (or optimized for performance), it might a good opportunity to do this in tandem -- though @brryan should probably weigh in here

Some questions:

• In my mind this approach would also simplify the implementation of adaptive timestepping. Do you concur with that?
• How big of a lift would a clean replacement be? (Or in other words, is it feasible to just do it and compare to the current approach in terms of performance and memory footprint?)

[lroberts36]

In my mind this approach would also simplify the implementation of adaptive timestepping. Do you concur with that?

I haven't thought about this aspect at all really. Assuming that implementing this would rely on having separate meshblock lists similar to GMG, I think it would be helpful.

How big of a lift would a clean replacement be? (Or in other words, is it feasible to just do it and compare to the current approach in terms of performance and memory footprint?)

I think the tree object in Mesh only gets used in ~20 places in the code, so that is easy to replace if we move away from the SearchAndSetNeighbors code. The lift seems like porting over the swarm stuff, but maybe it wouldn't be so bad for somebody who understood that part of the code.

[lroberts36]

Regarding the duplication of the MeshBlockTree, is this actually a concern right now? This was a concern in the original Enzo, which motivated the forrest-of-octree approach in Enzo-E. So far I considered our/Athena++ current approach as already on the lightweight side wrt to mesh management.

I have not carefully though about it, but when weakly scaling a simulation there certainly has to be a number of ranks at which the memory footprint of the mesh on a rank overwhelms the memory footprint of meshblocks on a rank. I think @Yurlungur mentioned to me that he actually ran into this issue at some point.

[Yurlungur]

I have not carefully though about it, but when weakly scaling a simulation there certainly has to be a number of ranks at which the memory footprint of the mesh on a rank overwhelms the memory footprint of meshblocks on a rank. I think @Yurlungur mentioned to me that he actually ran into this issue at some point.

I did... When doing a weak scaling test with AMR on Trinity. This was with Athena++