Parameters of PGD and optimal weighting/cic compensation

[dlanzieri]

Hi @modichirag , I compared the matter power spectrum computed with and without PDG. I need know what do you think about these results!
Initially I used kl, ks, mu, and alpha0 parameters obtained by running this Notebook with the following setting for the N-Body simulation :

nc = 128         
field_size = 5.   
nsteps=40
B=2

This is how the power spectrum looks like for all the snapshots:

Because of the alpha parameter we used for the PGD is scale factor dependent according the following relation:
alpha=alpha0*a**mu,
I tried to see how the amplitude of the power spectrum at small scales changes as a function of alpha0 parameter:

I see that the power spectra that better approximate the analytical one are the power spectra obtained with alpha0=0.008, alpha0=0.009, alpha0=0.01, for bigger and smaller values they start to lose power.
This is how the power spectrum looks like for all the snapshots using alpha0=0.01:

Now, I think that probably my scale translation is not right. The kl,ks obtained with the Vanessa's notebook should be in hMpc^(-1), but I'm not sure.
So, still using the the alpha0 obtained from the notebook (alpha0=0.00269), I tried different scale translation:

kl_array=[kl*0.7*0.5/(nc*B/ box_size),kl*0.5/(nc*B/ box_size),kl*2/(nc*B/ box_size),kl*2*0.7/(nc*B/ box_size) ]
ks_array=[ks*0.7*0.5/(nc*B/ box_size),ks*0.5/(nc*B/ box_size),ks*2/(nc*B/ box_size),ks*2*0.7/(nc*B/ box_size)]

Here the peaks of the PGD kernel :

and the power spectrum of only one snapshot:

From this plot looks like the right scaling factor for kl and ks are obtained multiplying them for h and 0.5 (this depends how we defined k) .
Assuming this is the better rescaling, I used it to see again how the power spectrum at small scales changes as a function of alpha0

And, alpha0=0.01 still seems like the better choice.

[modichirag]

Wow, this happened quite quickly, great job :)
I do not have a very good intuition of PGD so I am not sure how much useful comments I would have regarding the scaling factors. It might be good to check with people who have actually implemented this on Wednesday.
But a few remarks I had -

• it seems to do something, but maybe it will be easier to see the transfer function (ratio) plots, so say ratio of with, and without PGD with the theory spectra. It will also allow us to compare with figures in Biwei's paper, like Fig. 5
  https://www.osti.gov/pages/servlets/purl/1543937
• another good figure will be to show B=2 at the same time since this will tell which is better approach, increasing B, PGD or combination of the two.
• this is a 200 Mpc/h box now, right? I guess the plan is to use 200 Mpc/h box for analysis, and not a 50 Mpc/h box? The reason I am asking is because I do not know if PGD will help with the large scale modes that a small box looses.
• In the figures varying alpha and scaling, all the spectra have higher power than the no PGD case, right? (can it be used as a sanity check, that we never reduce power?)
• Regarding scaling, if we fit the parameter values which is what we will be doing ideally at the end of the day, the scaling issue should get fixed naturally, right? But this exercise is useful to set reasonable prior ranges and also see physically what is happening.
• Lastly, does it take us to the scales we need to push? Is the target k=2 and within 1-5% of the correct power?

[dlanzieri]

Thank you @modichirag !
So, I will do the plot you suggest :
1)See the ratio of matter power spectra compared to the jax power spectrum
2)See the effect of different B and PGD/NO PGD at the same time
3)Check that the power spectrum with PGD with different alpha0 is always higher that the one computed without PGD.

The Box I used is 128 Mpc/h because I had to adapt my setting to the different combinations found in the notebook of Vanessa. If I wanted to use a bigger boxsize, I should use a box size of 256 Mpc/h and only 10 steps for the simulation.
But, if we want to find best fit parameters ourselves, (as I understood from your message) we don't have that problem and in principle we could use the setting that we want.

[EiffL]

We now have a PGD implementation thanks to #90 but we still have an open issue of how to find the best weight function and apply CiC compensation, so I'm leaving this issue open. @modichirag if you want to have a look at that, go for it ;-)

[EiffL]

So @modichirag with Denise we looked at the shape of your custom weight function, and it looked very weird, it actually amplified the loss on small scales as far as I can tell...

[EiffL]

Ok, I think I fixed the decic function:

from flowpm.utils import r2c3d, c2r3d
from flowpm.kernels import fftk
nc = N

def decic(field, order=-2):
    kvec = fftk([nc, nc, nc], symmetric=False)
    kwts = [np.sinc(kvec[i]/(2*np.pi)) for i in range(3)]
    wts = (kwts[0]*kwts[1]*kwts[2])**order
    fieldc = r2c3d(field, norm=nc**3)
    fieldc = fieldc * tf.cast(wts, tf.complex64)
    field = c2r3d(fieldc, norm=nc**3)
    return field

Notebook here

[EiffL]

So, as discussed here LSSTDESC/DifferentiableHOS#10 (comment) it seems like offline PGD is failing, we are asking for too much displacement at once. But on the other hand online PGD seems to be doing quite well.

I am adding the decic function and will be adding back some online PGD code in this branch u/EiffL/pgd_fitting.

[EiffL]

So, last thing we discussed with @modichirag is that we may want to add a tool to fit the PGD params with a high res simulation.
Also, currently, the fitting strategy for online PGD, which starts by fitting at a=0.1 and then refines the params for lower redshifts, seems to get stuck in a sub-optimal fit. Here is an illustration of that:

• Fit at a=0.37:
  

• Fit at a=0.51:
  

We may want to fine tune a little bit the training script so that it can find a better fit at a=0.51.

[EiffL]

This might be good enough for our purposes for now, but @modichirag let me know if you want to dig a little deeper before we move on.

[EiffL]

So.... we are using this implementation of PGD for @dlanzieri project, and we are starting to wonder whether cosmology dependence of PGD is going to be a problem when computing gradients.

@modichirag @VMBoehm @Maxelee have you checked at any point how strongly correlated to Omega_m and sigma8 PGD params are?

[EiffL]

so here is a plot, apologies for it being very ugly, but computes the gradients of the matter power spectrum in a snapshot at z=0.
-- lines are grads computed on a sim that was using PGD, +- are without PGD, solid lines are the theory.

The top lines are derivatives wrt s8, bottom are for omega_c

It looks like the jacobian with respect to sigma_8 is larger than what it should be, and I think it might agree with what you are seeing @dlanzieri right?

[dlanzieri]

Yes, this what I found: