Public Data Access Overview / Discussion Forum

Hello,

I have a slightly related follow-up question from here on how to get to the initial particle grid of snapshot 0.
My goal for eg. Illustris-3-Dark (with 455^3 DM particles) is to store the DM particle displacement field vectors of a given snapshot in an array of shape (3,455,455,455). For this I would have to assign each snapshot 0 particle a position in a 455^3 mesh. Eg. particle with ID 1 will be placed at position (0,0,0), particle ID 2 at position (0,0,1) and so on. Ideally this assignment will happen in a way that this mesh position somewhat resembles this particle's true coordinates in snapshot 0. The snapshot 0 particle IDs don't really seem to be ordered in any kind of trivial way though (they are grouped but not simply left to right, or top to bottom).
What is the best solution for my problem here?

Thank you for the help in advance!

Best,
David

Hi David,

Note that snapshot 0 is z=20 (for TNG). There is also the actual IC file (z=127). Which is appropriate here?

I'm looking for the IC files for the illustris-dark sims, I think. If the particles there are placed on a grid, that should make it easy to obtain the displacement fields for any other given snapshot! I suppose the snapshot 0 files are not the IC files then?

Thank you for your quick support!

I just found the IC files on the website, I mistook snapshot 0 for the IC files :(
In what order are subsequent particle IDs placed on the grid? They don't seem to be placed on a regular grid one after the other

You are right, the ICs are not a grid, but instead a "glass" configuration. If you create a uniform 455^3 grid and deposit these particles into grid cells, I expect some cells will have zero, one, and two particles, so this mapping may not be so easy.

The actual ordering should be a 3D tiling of a small glass configuration, i.e. a particle distribution which is then replicated N times in each linear dimension. For TNG100-3, this is a 91^3 glass replicated 5x5x5. For TNG100-1 this is a 182^3 glass replicated 10x10x10.

So is there any way for me to put each particle onto its original point on a pre-IC cube lattice? Not really, since each glass config is a Poisson distribution and thus in each glass config each ID does not correlate to a specific point on a given cube lattice at all, right? :( So storing the displacement field vectors in an array with any kind of relation between IC particle position and array-position is not possible?

I see a few possibilities: (i) you can compute the displacement vector for every single DM particle between its position in the ICs and its position in any snapshot, e.g. z=0. (ii) on a coarser grid, e.g. 256^3 for Illustris-3-Dark or up to 1024^3 for Illustris-1-Dark, you could do a similar calculation as you originally described, where each grid cell would contain several DM particles, and thus the displacement vectors would be average/coarsened.

Yes I think option (ii) might be the only viable option for me since neighboring array entries storing the displacement field (of one, or more particles averaged) should also relate to particles being in close proximity in the simulation. Why did you use glass configs before perturbing the positions instead of just a simple cube lattice?

In any case thank you for your help!!

Hi David,

If you google "glass cosmological initial conditions" you will find many papers on the topic. I am not an expert on the topic, but from my perspective they are better for us in terms of the initial force accuracy and subsequent (especially early, linear) evolution, given the algorithms used in the AREPO gravity solver. I believe this depends on method, e.g. Zeldovich vs 2LPT, as well as starting redshift, i.e. how initially non-linear.

Makes sense, I suppose solvers are either more or less compatible with certain ICs depending on their exact architecture :)

Do you have the Illustris-Dark glass files lying around somewhere? I think using them would allow me to use a finer grid together with option (ii), maybe even a 455^3 grid for the 455^3 illustris-3-Dark particles

Hi David,

I put this here:

www.tng-project.org/files/GlassTile_n91.dat

to maybe make some testing with. You can refer to N-GenIC for details of its structure/usage.

Thank you for all the help! I think that should do it :)

Do you by any chance know how I can read the file? I suppose it is (just like the dummy file in NGenIC) a f77-unformatted binary file, but every try I made so far in Python and also by modifying the C code reading the glass file in NGenIC, produced only nonsense (both for the headers but also the pos-data) :((
If you don't know, no worries! At some point after trying some more stuff I'll find the correct way, no need for you to also spend time on it!

So from what I managed to deduce from reading through the relevant parts of the NGenIC code I still wasn't able to read your fortran sim output glass tile file properly. The header seems to contain 256 bytes, and the particle positions seem to be stored as a 4 byte float, together with the position and particle IDs. This yields a almost uniform particle distribution inside a cube with sidelength 1 (in simulation units I suppose). But the area around the origin seems to have a density around 10^3 OOMs higher than the rest of the cube, which doesn't make sense, plus my particle ID values are also nonsensical.
But only by offsetting the 256 byte header and then reading in the pos and velocity as 4 byte float arrays and the IDs as 8 byte ints, was I able to get to 91^3 particles, so something I am doing has to be somewhat correct.

If you also don't know how the file is formatted exactly maybe can you tell me which code you used to produce the glass tile? I suppose there I'll find the correct datatypes and header offsets.

Thank you again for the help!

Best,
David

Hi David,

The structure I see is:

Ntot = N * N * N

npart=lonarr(6)
massarr=dblarr(6)
time=0.0D
redshift=0.0D
flag_sfr=0L
flag_feedback=0L
npartall=lonarr(6)
flag_cooling= 0L
num_files= 1L
BoxSize = 0.0D

bytesleft=120
la=intarr(bytesleft/2)

BoxSize = 100.0D
npart(1) = Ntot
npartall(1) = Ntot

pos= fltarr(3, Ntot)

for i=0L, N-1 do begin
  for j=0L, N-1 do begin
    for k=0L, N-1 do begin
      pos(0, (i*N+j)*N+k) = (i+0.0)/N * BoxSize
      pos(1, (i*N+j)*N+k) = (j+0.0)/N * BoxSize
      pos(2, (i*N+j)*N+k) = (k+0.0)/N * BoxSize
    endfor
  endfor
endfor

openw,1,fout,/f77_unformatted
writeu,1, npart,massarr,time,redshift,flag_sfr,flag_feedback,npartall,flag_cooling,num_files,BoxSize,la
writeu,1, pos
close,1

Hey Dylan,

yeah, I also tried offsetting the header by summing up the bytes specified in that code, but it gave me nonsense. My next approach now was to simply run NGen-IC without the perturbations, so it would simply place each glass tile 5 times along each dimension. That worked (as in it outputs a uniform field of 455^3 particles) but equal particle IDs between my Ngen-IC output and the Illustris IC catalog don't share even remotely the same positions, such that the correlation between the two flattened pos catalogs is only 0.3 :( My NGen-IC outout Do you still have the ics.param file used to launch NGen-IC? maybe some of my params are not the same which leads to the error. I'm appending my file to this message.

Thanks again for all your help!!

Nmesh            128       % This is the size of the FFT grid used to 
                           % compute the displacement field. One
                           % should have Nmesh >= Nsample.

Nsample          128       % sets the maximum k that the code uses,
                           % i.e. this effectively determines the
                           % Nyquist frequency that the code assumes,
                           % k_Nyquist = 2*PI/Box * Nsample/2
                           % Normally, one chooses Nsample such that
                           % Ntot =  Nsample^3, where Ntot is the
                           % total number of particles


Box              75000.0  % Periodic box size of simulation

FileBase         ics                 % Base-filename of output files
OutputDir        ./ICs/              % Directory for output

GlassFile        GlassTile_n91.dat  % File with unperturbed glass or
                                  % Cartesian grid

TileFac          5                % Number of times the glass file is
                                  % tiled in each dimension (must be
                                  % an integer)


Omega            0.2726       % Total matter density  (at z=0)
OmegaLambda      0.7274       % Cosmological constant (at z=0)
OmegaBaryon      0.0       % Baryon density        (at z=0)
HubbleParam      0.704       % Hubble paramater (may be used for power spec parameterization)

Redshift         127        % Starting redshift

Sigma8           0.9       % power spectrum normalization



SphereMode       1         % if "1" only modes with |k| < k_Nyquist are
                           % used (i.e. a sphere in k-space), otherwise modes with
                           % |k_x|,|k_y|,|k_z| < k_Nyquist are used
                           % (i.e. a cube in k-space)


WhichSpectrum    0         % "1" selects Eisenstein & Hu spectrum,
                   % "2" selects a tabulated power spectrum in
                           % the file 'FileWithInputSpectrum'
                           % otherwise, Efstathiou parametrization is used


FileWithInputSpectrum   input_spectrum.txt  % filename of tabulated input
                                            % spectrum (if used)
InputSpectrum_UnitLength_in_cm  3.085678e24 % defines length unit of tabulated
                                            % input spectrum in cm/h. 
                                            % Note: This can be chosen different from UnitLength_in_cm

ReNormalizeInputSpectrum   1                % if set to zero, the
                                            % tabulated spectrum is
                                            % assumed to be normalized
                                            % already in its amplitude to
                                            % the starting redshift,
                                            % otherwise this is recomputed
                                            % based on the specified sigma8


ShapeGamma       0.21      % only needed for Efstathiou power spectrum 
PrimordialIndex  1.0       % may be used to tilt the primordial index, 
                           % primordial spectrum is k^PrimordialIndex


Seed             123456    %  seed for IC-generator


NumFilesWrittenInParallel 2  % limits the number of files that are
                             % written in parallel when outputting


UnitLength_in_cm          3.085678e21   % defines length unit of output (in cm/h) 
UnitMass_in_g             1.989e43      % defines mass unit of output (in g/cm)
UnitVelocity_in_cm_per_s  1e5           % defines velocity unit of output (in cm/sec)

Hi David,

We have Nmesh 910, Nsample 455, Sigma8 0.809, Seed 5463. Note that we also use a tabulated TF, not Eisenstein & Hu, but perhaps you can get close with the above changes.

I tried with the additional params with almost no change (still only 0.3 correlation between the two catalogs). I think those settings are also mainly important for the perturbations themselves right? In any case its really strange that the particle IDs don't correspond to one another. I suppose I could simply look at each Illustris_ic particle and try to place it in a 455^3 grid according to its posistion. Then I would only need to keep track of which ID is in which grid-cell to calculate the displacement field at a given snapshot

I now also tried simply creating the ICs using your params and glass tile with the unmodified NGen-IC code and I still don't get particles in a similar order as in your Illustris-3 ic catalog. I don't think that a tabulated TF will make a difference there right? Maybe the tile you sent me was from a different simulation?

It is definitely the correct glass file for Illustris-3. Do you have CORRECT_CIC turned on in the N-GenIC makefile? It should be. Other than that, I'm not sure where the remaining differences could be. I think you can see how the positions in the glass should be tiled, and assigned IDs. If you look at the actual Illustris-3 initial conditions file I hope that you can find the correspondence?