Sorting by inital positions

I have added a feature to Walter's code, namely each particle carries about a memory of where it started it's life.   I have done several runs, but all of the following plots refer to initial conditions that were a power law density profile with $\rho \propto x^{-1/2}$.   I have defined 5 logarithmic bins in the particle initial positions.  Below is the phase space diagram, with different colors representing particles that started in different initial positions, the reddish brown ones began life the bin farthest out, and the dark blue began life in the innermost bin.  Each plot is zoomed in from it's predecessor (you can click on the plots to view them larger):

Just for illustration, here are the same types of plots, but for the initial conditions (note the aspect ratio is no longer square):

To get a better idea of how the different colors are distributed around phase space, I did a transformation on these data to fit all the information onto a single plot.  Below I have defined a distance in phase space which is $r=\sqrt{x^2 + v^2}$.  Then instead of plotting $(r,\theta)$, I plot $(r^{0.1}, \theta)$ so that we can see the inner parts. 

Here are the initial conditions in this view:

These plots demonstrate that the outermost bins retain their identity, but the inner bins get mixed in the ordering of their positions.  All of the velocities increase, but the inner bins tend to have only high velocities and small positions.  The phase space orbits for the outer particles and inner particles seem qualitatively different.  This result is robust to the phase of the winding:  I checked that it holds true for the last several time dumps, and also for different runs of the same initial conditions (different runs have different numbers of particles).

Another useful view is to look at scatter plots of $|x_i|$ versus $|x_f|$ and $|v_i|$ versus $|v_f|$.  Here they are, the black lines show $x_f=x_i$ and the analogous for the velocity.

You can see that the positions on average are collapsing inward, as expected, but the innermost particles do not collapse inward as much as some of their outer counterparts, in fact the very innermost ones seem to be thrown out of the central region.  The velocities are all increasing, with the innermost particles seeing the largest increase in velocity.  This basic structure is set up very early in the run, here is the first time dump after the start of the run, at 20 ATUs.

I was puzzled about that tail in the outer particles, I thought maybe looking at these plots in a slightly different way might shed some light.  Here are $x_i$ versus $x_f$ and $v_i$ versus $v_f$ plots, where again I just zoom in rather than doing the polar trick.  The axes have the same ranges as in the beginning of this post (except v_initial, which I scale down in the same proportions):

So you can see the velocity tail in a different form, in the first of these, I am still thinking aobut why it's there.  The innermost bins are doing something qualitatively different than the outer three, in my opinion.

There are still many more plots to make, I will start checking some of the other power laws.  Also, if I want to do warm runs, I need to repeat the modification to Walter's code, adding one more field (and since it's the last available bit, that is it, can't add any more data to carry around without doing something hard).    I will try and make some more plots this weekend, have good one.

Cold runs, recentered

Here are the power law cold runs, done with the updated recentering algorithm.   The initial conditions for these runs is $\rho \propto x^a$ with $a=-0.5, -0.7, -0.8, -0.3, +0.5$  They all have 10k particles.

Also, here is the run with $a=-0.2$, and 80k particles.

The cold runs with slope <-0.5 appear to exhibit the same flattening in the center as the warm runs.  Qualitatively, the warm and cold appear indistinguishable.  However, runs with initial slope >-0.5 do not appear to retain their original scaling, and they do not go to the attractor behavior, even when they are cold.   The run with initial slope >0 appears to do something like the attractor behavior, oddly enough. 

Recentering fixed

I have modified my code that generates the density plots.  Whereas before, the scheme was

1) Bin up so that there are N particles per bin
2) Compute the density
3) Find the position of the density maximum
4) shift all the particle positions so that the density maximum is at x=0
5) Bin up again with N particles per bin
6) Compute the density and plot

This worked relatively well for the runs where the particles began in an equally spaced configuration, but once I started the initial conditions with a power law, suddenly a much higher fraction of the particles wind up in the middle, but the offset from 0 was a similar magnitude as the equally spaced case.  Therefore, for the power law runs, when I got to step 5, sometimes with the new binning, the density maximum was again off-shifted from zero by some number of bins.  I have solved this issue by iterating the process until the shift is at most only one bin, then I stop.   Here are the plots for the warm runs, I think they are significantly improved, and I am also relieved to see the appearance of the flattening in the center of the -0.8 run.  I am somewhat surprised that it did not show up in the center of the -0.5 run.

I think it is advisable to remake several of the key plots using the new algorithm, even for the runs that were of constant initial density.  I will do so, and post them in the next post.

Code doing strange things.

I was investigating the part of my analysis that is responsible for recentering on the middle of the collapsed object.  I have not figured out whether it is working the way I expect it to, or not, but in the process I discovered that the code is behaving strangely in some of the runs.  For example, here is a normal sequence, in the warm sim with initial density profile $x^{-0.8}$, dumping out the particles every 20 ATUs (time units).  This run had 10004 particles in it.

Now here is a run that crashed (in the sense that all the positions were NAN, not in the hanging up on timestep 0, which is the other failure mode).  The only difference between this run and the previous one is that this one had 10007 particles, rather than 10004.  And of course the thermal part of the velocity has a random component.  What happens is very hard to explain:

Note the x axis.  The very next time dump was all NAN.  10010, 10013, 10016, and 10019 particles all looked identical to the blue run above, and ran in a well behaved way until the end.    If you zoom into one of those regions in Dumpout 3 of the crashing run, it looks like this:

So some of the coherence in the structure is still being maintained...   Thoughts anyone? 

Warm IC's on power law density profile.

Okay, I am finding this officially weird.  I initialized the particles with a power law density profile of $x^{-0.5}$ and $x^{-0.8}$.  These were cases which, in the cold runs, maintained their scaling, even though particles collapsed inward.   I added a 5% random velocity perturbation, making these runs slightly warm.  I expected a flattening to happen in the middle, but no such thing.  See here, initial and final states for the the 0.5 case:

And now here for the 0.8 case:

The spikes that are occurring are now starting to get worrysome... it means that either there are two density maxima occurring in this system, or for some reason my code that computes the offsets is not working correctly for these runs.  Since these are averages of 6 runs, it might be only one run which has not been centered properly, but since we really want to be looking at the central regions, it's important to make sure.   One thing to notice: the particles don't seem to be collapsing inward the way they have in other runs, which could again indicate a problem with recentering.  I will look into that next, but I wanted to share these preliminary results.

density profile tests

I have done a few more tests on the initial density profiles.  This work is taxing and slow to proceed because the runs are often very slow, and many of them will not run at all, but when I start them it is hard to tell the one case from the other.  The following experiment had an initial density profile of $x^{-0.2}$ and 80k particles.  I used that many particles, so that the initial density profile would extend several orders of magnitude to the interior.  I did this to test the hypothesis that if the initial scaling is already in place, the system would not go to the $x^{-1/2}$ attractor.  However, the hypothesis was incorrect, as you can see in the plots below (initial and then final density profiles of the particles):

   Each of these took 4.5 days.  In fact I thought the runs had hung (but as an experiment, didn't kill it).  I found that 3 of them had been completed when I got back from the break (I think I will kill the other three, but I'll wait a day, because I am curious.  Run 4 seems to have permanently hung on the 0th step, which is the main failure mode, but this time I know how long I have to wait to find out if it is simply running very very slowly). 

   It's somewhat hard to be quantitative about this, but the final configuration looks to me to scale roughly like the initial scaling in the very interior, but roll to an $x^{-1/2}$ scaling over a few decades until the final rolloff.    This is qualitatively different from the cases with initial power law steeper than $x^{-1/2}$.  For example, before the break I was able to get $x^{-0.7}$ to run.  Here are the initial and final configurations:

I just wanted to make sure that what we had seen in the $x^{-0.8}$ scaling (from my last post) was not a fluke.   I have not successfully been able to get anything with initial slope steeper than $x^{-0.8}$ to run. 

My next steps, I think, will be to look at a warm run with initial scaling $x^{-0.8}$ and/or $x^{-0.5}$, and as a related question, find out what happens if I take one of the old cold runs and add a little velocity perturbation AFTER it has reached its final state.  Also, I have been trying to look at the analytic stuff Scott and Walter sent me, to see if I can see any qualitative difference between slopes steeper and shallower than $x^{-1/2}$. 


Comments anyone?

Different initial density profiles

I found it very tricky to get these working.  They are not all working yet.  But here are a few that did, with 10k particles, 6 runs of each.  The initial velocities are set so that there is a coherent collapse going on, the same as for the equally spaced particles, namely:

V=shot.V(n) = -V*sin(shot.X(n));

These runs take a long time somehow.  Anyway, what happened is sort of interesting, I think.  Here is a plot of the initial conditions, these are all power law densities $\rho = A x^a$.  The slopes I picked to test out were a= -0.5, -0.8, -0.3, -2.0, +0.5  . 

I have not plotted the -2.0 because I have not succeeded in getting that set of initial conditions to run at all yet, I have not diagnosed what the problem could be.    So, the final density profiles are here:

 So, my interpretation of this is as follows: -0.8 seems to stay put at its initial slope, (as obviously does -0.5), whereas  -0.3 and +0.5 seem to drift to the -1/2 attractor.  However it seems unlikely that there is anything special about the particular value -0.8.  The big difference I see is that since I fix the initial spacial interval (which thus determines the normalization A), for a fixed number of particles the -0.8 initial configuration extends much farther into the center than the -0.3 and +0.5.  But all the runs collapse inward so that at the final configuration, they all have lots of particles interior to $10^{-3}$.  My guess is that if there are already particles in the interior, and they follow a different power law slope, then apparently they stay there, but if there are no particles in there to begin with, when they get there they do the $r^-1/2$ behavior.   I am pretty sure that adding a bunch more particles to the -0.3 run would make the initial configuration extend farther inward, and then we could test this.  I might make it -0.2 too, so it's easier to distinguish from -0.5 at the end. 

Incidentally, I also did a -0.8 run where the initial conditions were expanding rather than contracting (but the initial density configuration looks identical to what you see in the first figure).  This is what you get:

 So the stability of the inner slope does not seem to care about what the initial velocities are doing, as long as there is some long range coherence (and not thermal behavior between the particles).   Ignore the spike, my code was having trouble recentering one of the 6 runs.   I might try running some warm initial conditions on these non-uniformly-spaced runs, to see if it still flattens itself out in the central region.