BradGreig / Hybrid21CM

1 stars 3 forks source link

CoreLightConeModule build_model_data always the same lightcone #32

Open NGillet opened 5 years ago

NGillet commented 5 years ago

Hi,

I am trying to do a likelihood_NeutralFraction MCMC : 1) CoreLightConeModule

core_LC = mcmc.CoreLightConeModule( 
    redshift = 5.5,                 
    max_redshift = 12.0,             
    user_params = dict(       
        HII_DIM = 50,         
        BOX_LEN = 125.0
    ),
    z_step_factor=1.04,   
    regenerate=False,
    nchunks=4,
    do_spin_temp=True,
    change_seed_every_iter=False,
) 

2) likelihood

likelihood_NeutralFraction = mcmc.LikelihoodNeutralFraction()  

3) MCMC:

chain = mcmc.run_mcmc(
    core_LC, 
    likelihood_NeutralFraction, #[likelihood_planck, likelihood_NeutralFraction], #, likelihood_greig],
    datadir='data',          # Directory for all outputs
    model_name=model_name,   # Filename of main chain output
    params=dict(             # Parameter dict as described above.
        F_STAR10   = [ -1.5, -2.5, -0.5, 1.0  ],
        F_ESC10    = [  -1.0, -2.0,  0.0, 1.0  ],
               ),

    walkersRatio=4,          # The number of walkers will be walkersRatio*nparams
    burninIterations=0,      # Number of iterations to save as burnin. Recommended to leave as zero.
    sampleIterations=10,  # Number of iterations to sample, per walker.
    threadCount=8,          # Number of processes to use in MCMC (best as a factor of walkersRatio)
    continue_sampling=False, # Whether to contine sampling from previous run *up to* sampleIterations.
) 

The problem is that the return lnL is always 0. I looked at the core and likelihood module and found that the actual modeled neutral fraction is identical between walkers and iterations:

In CoreLightConeModule I print the params (they do change as expected), but after the computation of the lightcone lightcone.global_xHI is always the same !

    def build_model_data(self, ctx):
        # Update parameters
        astro_params, cosmo_params = self._update_params(ctx.getParams())

        print( 'IN BUILD astro_params :', astro_params.F_ESC10, astro_params.ALPHA_ESC )

        print( self.redshift[0], self.max_redshift, self.z_step_factor, self.initial_conditions_seed )

        # Call C-code
        lightcone = p21.run_lightcone(

            redshift=self.redshift[0],
            max_redshift=self.max_redshift,

            astro_params=astro_params, 
            cosmo_params=cosmo_params, 

            user_params=self.user_params,
            flag_options=self.flag_options,

            z_step_factor=self.z_step_factor,
            regenerate=self.regenerate,
            random_seed=self.initial_conditions_seed,
            write=self.io_options['cache_ionize'],
            direc=self.io_options['cache_dir'],
        )

        print( 'IN BUILD xHI :', lightcone.global_xHI )

        ctx.add('lightcone', lightcone)
>>> 
IN BUILD astro_params : -0.422045384281 -0.5
IN BUILD astro_params : -1.15392864038 -0.5
IN BUILD astro_params : -1.29628662729 -0.5
IN BUILD astro_params : -1.9993342796 -0.5
5.5 12.0 1.04 396696065585
5.5 12.0 1.04 396696065585
5.5 12.0 1.04 396696065585
5.5 12.0 1.04 396696065585
IN BUILD astro_params : -1.98190265686 -0.5
IN BUILD astro_params : -0.762196428892 -0.5
5.5 12.0 1.04 396696065585
5.5 12.0 1.04 396696065585
IN BUILD astro_params : -1.06946089891 -0.5
5.5 12.0 1.04 396696065585
IN BUILD astro_params : -0.565129474849 -0.5
5.5 12.0 1.04 396696065585
IN BUILD xHI : [ 0.96435517  0.95072514  0.93239552  0.90946776  0.8801595   0.84389597
  0.80032223  0.74773973  0.68238193  0.60641336  0.51849377  0.42002162
  0.31888756  0.22795989  0.15128654  0.09266903  0.04975107  0.02413185
  0.00956937]IN BUILD xHI : [ 0.96435517  0.95072514  0.93239552  0.90946776  0.8801595   0.84389597
  0.80032223  0.74773973  0.68238193  0.60641336  0.51849377  0.42002162
  0.31888756  0.22795989  0.15128654  0.09266903  0.04975107  0.02413185
  0.00956937]
lnprob :  0
lnprob :  0
IN BUILD xHI : [ 0.96435517  0.95072514  0.93239552  0.90946776  0.8801595   0.84389597
  0.80032223  0.74773973  0.68238193  0.60641336  0.51849377  0.42002162
  0.31888756  0.22795989  0.15128654  0.09266903  0.04975107  0.02413185
  0.00956937]
lnprob :  0
IN BUILD xHI : [ 0.96435517  0.95072514  0.93239552  0.90946776  0.8801595   0.84389597
  0.80032223  0.74773973  0.68238193  0.60641336  0.51849377  0.42002162
  0.31888756  0.22795989  0.15128654  0.09266903  0.04975107  0.02413185
  0.00956937]
lnprob :  0
IN BUILD xHI : [ 0.96435517  0.95072514  0.93239552  0.90946776  0.8801595   0.84389597
  0.80032223  0.74773973  0.68238193  0.60641336  0.51849377  0.42002162
  0.31888756  0.22795989  0.15128654  0.09266903  0.04975107  0.02413185
  0.00956937]

I am stuck here, I do not know why because the actual input params are different and should change the resulting xHI.

I try to switch of regenerate=False or to put the random_seed generator. I get a warning for the first and I am not sure the second work too.

steven-murray commented 5 years ago

Hi @NGillet , thanks for the report. Can you copy/paste this to https://github.com/21cmFAST/21CMMC/issues/new ?

BradGreig commented 5 years ago

@NGillet , the default flag_options struct has USE_MASS_DEPENDENT_ZETA set to False. This should be set to True if you are using F_STAR and F_ESC.