

Here are a number of things that I (AFN) think need to be done with 
vine as it evolves into the future. Some of these are large scale 
structural changes to the code, some are lower level detail things.
Some are in between.



Things reasonable for version 1.1/1.2 etc:

Improvements to SPH

  -Add entropy conserving formulation, following some Price/Monaghan
   work. (See also Springel)

  -Add multimaterial capability in robust form (currently present but
   essentially inactive). 

       -Multi-material EOS stuff
       -Strength (e.g. for moon/earth or asteroid collisions etc.
                  see various Benz/Asphaug/Canup/Cameron work)
       -Reaction Chemistry (e.g. nuclear reactions for supernovae,
                            chemical reactions for molecular clouds, 
                            disks, etc)

  -Add radiative transport. Critical step for doing accretion disks
   better, such as in Jovian planet formation stuff.

  -Would be nice to tweak the density and sphaccel loops to be 
   vectorizable. Currently they don't seem to be because of the
   existence of references to arrays like vel(:,jpart). Interesting
   fact is that the intel compiler vectorizes a number of loops
   in vine that other compilers don't. Why?

  -Is Moving Least Squares Particle Hydrodynamics (MLSPH) a reasonable
   thing to try?

  -Almost certainly a good idea to set default ranges for neighbor counts
   more restrictively than now. Better behavior all around, both for hydro
   and gravity.

Improvements to the integrator

  -Add KDK version of leapfrog. Better for stability. Also better for
   parallel scalability, since particles line up on time step increments.
   That means about half the number of cycles through the update loops,
   which scale far less well than the force calculations.

Improvements to the derivative calculations 

  -Add ability to do the gravity and SPH calculations in parallel
   when only a few nodes are active. Will this improve scalability
   by given more processors something to do?

Improvements to the forcelaw module

  -Add Price/Monaghan energy conservative gravitational potential
   formulation for varying smoothing lengths.

  -Go to higher order multipoles, or at least have that as a compile
   time option. Big question is, will this actually be effective
   in practice?

  -When few particles/clumps are active, send group node lists to
   more than one processor rather than sequentially, to improve
   scalability on otherwise unused processors?

Improvements to the tree traversal process

  -Would be nice if I could get activate_tree to work in O(nactive)
   rather than O(npart) operations. It is annoying and painful when
   there are only a few active particles.

  -Would be nice to extract out the MAC criteria from their current
   locations deeply embedded into the traversals, into (e.g.) the
   forcelaw.f90 module, where they actually belong. That makes the
   tree traversal that much more flexible in terms of being able to
   set some other opening criterion for some other force law etc.

  -The tree traversals in vine are vectorizable, according to
   Makino 1990. Are they still vectorizable on modern simd machines? 

Improvements to boundary treatments

  -Fill out the currently blank 'reflneigh' routine, so that we
   can do reflecting boundaries.


Improvements to the tree build process.

   -The parallel sort algorithm in treesort.f90 is not very scalable 
    in vine it appears (or I didn't implement it correctly?). Some
    time ago, I found another algorithm that I think will work a lot
    better. Specifically, the description in:

       Tsigas, P., Zhang, Y. "A Simple, Fast Parallel Implementation
       of Quicksort and its Performance Evaluation on SUN 
       Enterprise 10000" 2003, Proceedings of the 11th Euromicro
       Conference on Parallel Distributed and Network based 
       Processing, p372-381 IEEE press

   -Also could do all ndim sorts at once. This would allow parallelizing
    the makegrid/initlists routines. Both this and the sort improvement
    above would improves build scalability.

   -Separate out sph and nbody particles into their own trees,
    joined only at the root. Generalizes easily to many different
    kinds of particles too. Improves force calculation scalability
    for mixed calculations.

   -Do better in the getnearest subroutine, which is most of the
    build cost. Currently have very tiny loops here. Instead of
    dwall(iw,id) it would be better to have dwall(iw) where iw
    goes from 1..2*ndim. Makes the loops longer. Other things 
    might also improve the speed, but fall in the category of `art'
    and I can't really describe them very well.

   -Could in principle also add oct tree option. Would this be better
    for scalability of the build? What does it do to the traversal 
    speed?

   -The tree revision process currently uses a separate copyin routine.
    Some evidence now suggests this may not be the way to go. Combine?

   -Performance, especially parallel performance, of gethclumps and friends
    is much less than desirable. 

   -Allow partial tree revisions. Is this a good idea?


Things reasonable for a version 2.0

  -Add distributed memory parallelism over the top of the current shared
   memory parallelism. This will allow larger simulations than are commonly
   possible now, since large shared memory machines are comparatively
   rare and clusters much more common. Since clusters are now getting
   lots of chips/cores per node, it will also be a nice fit in the sense
   of needing only one MPI rank per box rather than many, while doing
   shared memory parallelism inside the node.

   Questions: What is the best way to domain decompose? Simple stuff
              like bisecting boxes and etc, as the distributed sph/nbody
              folk already do? That doesn't fit very well with the
              binary tree though...
              
              Or something more complex? Many years ago now, I saw a
              paper that did "Voronoi tesselation" to determine where
              to put particles (Yahagi et al. ApJS 124, 1). This method
              combines distance from some center and work together to create 
              a metric of whether or not to keep a given particle on the 
              current processor or move it to another. Will this method be 
              subject to wierd pathologies where you get reasonable behavior 
              for some tiles, but strange behavior for others? For example,
              would it end up with one processor that had almost no particles
              on it, and the ones it did have all lived on the furthest
              outskirts of the system, as flyaways from the main system?

              Who knows?



