gandalfcode / gandalf

GANDALF (Graphical Astrophysics code for N-body Dynamics And Lagrangian Fluids)
GNU General Public License v2.0
44 stars 12 forks source link

NeighbourList, Dust #127

Closed rbooth200 closed 7 years ago

rbooth200 commented 7 years ago

Here are some more changes based upon the neighbour manager concept. Giovanni suggested earlier that we were paying a heavy price in copying the particles unnecessarily in the dust tests, which he implemented as a temporary fix while we were getting the neighbour manager working. I can now indeed confirm that those copies were coming at a high cost. The changes introduced here reduce the running time of the dust tests from 138s to 99s on my laptop with OpenMP and 4 threads.

What I've done here is to introduce a small NeighbourList template class. This class is a simple container for a single neighbour list, e.g. the direct list. It also hides the iteration over the indices of the neighbour array from the user and allows us to return the list of neighbours as if it were a vector. It's also cheap to copy around a NeighbourList since it only holds a reference to the internal arrays of the NeigbourManager. The implementation is safe since the underlying vectors are hidden so the users can't accidentally kill the list of neighbours or mess up the indexing, but they are free to modify the neighbour particles.

On top of this class I've implemented an iterator, NeighbourIterator, which might come in useful, but is not currently being used.

I think this pull request goes a long way towards implementing our ideal interface for neighbour lists. There are still a few open questions though: 1) How should we deal with the indices of the neighbours original location in memory? It would be simple enough to add a reference to neib_idx to the NeighbourList if we want, so we could provide a uniform index-based way of accessing this information. However, this approach would deviate from a simple standard library container type interface, and it's not obvious how to extend this to use the iterators (or whether we even should). 2) How should we implement the interface for the gravity neighbours, which need three different neighbour lists? The current implementation can't be default constructed, so we can't just pass three of these objects in and have the function fill them. I can see two plausible ways to get around this: 1) We could get around this by holding pointers rather than references. I'm reluctant to do this though, because its likely to introduce lots of edge cases that would need guarding against and that seems like a bad idea. 2) Splitting up the 'generation' of the neighbour lists for a particular particle from the return of them. To do this, either the user would have to pre-generate the reduced lists for a particle by calling a function, or we could hide the generation and force the user to pass in the particle and masks each time they want to generate a list. We'd have to either recalculate the lists multiple times, or try to cache the input and only recalculae when it changes.

Thoughts?

giovanni-rosotti commented 7 years ago

How should we deal with the indices of the neighbours original location in memory? It would be simple enough to add a reference to neib_idx to the NeighbourList if we want, so we could provide a uniform index-based way of accessing this information. However, this approach would deviate from a simple standard library container type interface, and it's not obvious how to extend this to use the iterators (or whether we even should).

I don't know what you mean with "uniform index-based way", but in any case, I see two ways:

  1. I think that the best solution is to costruct another NeighbourList object which returns references to the real particles. Currently this cannot be done because NeighbourList expects vectors while the real particles are in plain arrays, but I don't see why it cannot be done in the future... I think this is the best solution to avoid class proliferation or interface complciations
  2. Extend the class and add begin_real and end_real functions, which return a different kind of iterator. The disadvantage is that when eventually we will switch to c++11 we won't be able to use range-based for. Unless I am mistaken, in c++98 the names begin and end are just conventions; I don't think that the STL or the language actually use them anywhere.

How should we implement the interface for the gravity neighbours, which need three different neighbour lists? The current implementation can't be default constructed, so we can't just pass three of these objects in and have the function fill them. I can see two plausible ways to get around this: We could get around this by holding pointers rather than references. I'm reluctant to do this though, because its likely to introduce lots of edge cases that would need guarding against and that seems like a bad idea.

i don't like it either

Splitting up the 'generation' of the neighbour lists for a particular particle from the return of them. To do this, either the user would have to pre-generate the reduced lists for a particle by calling a function, or we could hide the generation and force the user to pass in the particle and masks each time they want to generate a list. We'd have to either recalculate the lists multiple times, or try to cache the input and only recalculae when it changes.

i'll think more about this... if these are the alternatives, we should hide the generation and cache the results. but maybe there are simpler ways. what if we just go for the simplest solution of returning a struct with three neighbour lists inside? it's basically what we are doing already (but currently we return only three integers). i don't know if the compiler can do return by value optimization in this case, but these objects are so cheap to copy that we shouldn't care...

rbooth200 commented 7 years ago

So after more thinking, problem 2 does have an obvious solution within the current framework. I see no reason why we can't just return a struct containing the three different neighbour lists. If desired I can list a few other good arguments as to why returning a struct of neighbour lists is the best solution.

rbooth200 commented 7 years ago

I'll implement the struct based solution this evening. We should chat briefly about getting the true particles.

rbooth200 commented 7 years ago

Ok, I believe I have a working solution that fulfills the requirements we discussed.