bccp / nbodykit

Analysis kit for large-scale structure datasets, the massively parallel way
http://nbodykit.rtfd.io
GNU General Public License v3.0
111 stars 60 forks source link

Zeldovich power spectrum #578

Open franlane94 opened 5 years ago

franlane94 commented 5 years ago

plot.pdf

Hi,

I've been trying to calculate the Zeldovich power spectrum and I get some really weird behaviour. I'm not sure if it's a problem with the mcfit functions.

Thanks

rainwoodman commented 5 years ago

It looks strange indeed -- those oscillations. Does smell like something related to cancellation of bessels.

Could you post the portion of your script that's relevant, and warning messages if there are any?

Any comments from @eelregit?

eelregit commented 5 years ago

@franlane94 Did you see my replies to the issue you opened in mcfit here?

@rainwoodman How is new life? :)

You may remember that I told you I changed the default options of mcfit, because many users assumed the opposite without reading the docstring. Nick's pyrsd also have/had some small problem due to this.

I changed lowring from True to False in __init__, and extrap from True to False in __call__. We can fix them in nbodykit by explicitly giving all the keyword arguments.

franlane94 commented 5 years ago

Hi,

Changing the keyword arguments to what you suggested fixed the problem I was having the the P2xi function. I think the problem is arising from the classes calculating the Bessel integrals.

Class to calculate J0 integral

class ZeldovichJ0(mcfit.mcfit):

def __init__(self, k):

    self.l = 0
    UK = mcfit.kernels.Mellin_SphericalBesselJ(self.l)
    mcfit.mcfit.__init__(self, k, UK, q=1., lowring=False)

    # set pre and post factors
    self.prefac = k
    self.postfac = 1 / (2*np.pi)**1.5

They are called with extrap=True. I tried switching around the lowring and extrap arguments but I was still getting the strange oscillations.

Thanks!

eelregit commented 5 years ago

Nick has already fixed those optional arguments, so it shouldn't be affected by my changes of defaults.

Did you get better results previously with the same code, or is this the first time?

FFTLog is a fragile algorithm, at least to me. You never know when these oscillations will show up until you try it. They might arise due to some noise, missing range of the input, or Gibbs ringing. Of course it's great when it works.

I usually play with N, lowring, extrap, and q to see if I can get rid of the oscillations. For extrap, one can try True, False, and also (True, False) or (False, True). For q, you can try 1.5, or 0.5. Can you post what your input P(k) looks like?

franlane94 commented 5 years ago

This is the first time I'm using the code. Here is a plot and txt file of my k and P(k) values

plot.pdf

k_and_P_vals.txt

eelregit commented 5 years ago
class ZJ0(mcfit): 
    def __init__(self, k): 
        UK = Mellin_SphericalBesselJ(0) 
        mcfit.__init__(self, k, UK, q=1, lowring=False) 
        self.prefac = k

class ZJ1(mcfit): 
    def __init__(self, k): 
        UK = Mellin_SphericalBesselJ(1) 
        mcfit.__init__(self, k, UK, q=0, lowring=False) 
        self.postfac = 1 / self.y

worked for me:

image

This is strange ... @rainwoodman

rainwoodman commented 5 years ago

we shall check the version of software used here. maybe we need a new nbodykit release or mcfit release?

On Mon, Apr 15, 2019, 11:11 PM Yin Li notifications@github.com wrote:

class ZJ0(mcfit): def init(self, k): UK = Mellin_SphericalBesselJ(0) mcfit.init(self, k, UK, q=1, lowring=False) self.prefac = k class ZJ1(mcfit): def init(self, k): UK = Mellin_SphericalBesselJ(1) mcfit.init(self, k, UK, q=-1, lowring=False) self.postfac = 1 / self.y

worked for me:

[image: image] https://user-images.githubusercontent.com/7311098/56185988-e6f95800-6058-11e9-8e95-9e31ed0fddb6.png

This is strange ... @rainwoodman https://github.com/rainwoodman

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/bccp/nbodykit/issues/578#issuecomment-483524653, or mute the thread https://github.com/notifications/unsubscribe-auth/AAIbTL0OS5TO3hdyooK4AQd9Gqi-Ib-Pks5vhWmDgaJpZM4cr4Xf .

eelregit commented 5 years ago

@franlane94 I tried to run ZeldovichPower, with the latest version of nbodykit and mcfit, and could not reproduce the problem. Could you update them and try again? Please post your script if the problem persists.

rainwoodman commented 5 years ago

@franlane94 The quickest way to check for the version is

import nbodykit;print(nbodykit.__version__)
import mcfit;print(mcfit.__version__)