[BUG] Non-monotonic reionization when node redshifts are irregular

[samgagnon]

Describe the bug:

To Reproduce:
Steps to reproduce the behavior:

1. Run a lightcone with node_redshifts defined as follows
   node_redshifts = p21c.get_logspaced_redshifts(min_redshift = min_redshift, max_redshift = max_redshift, z_step_factor = 1.02)
node_redshifts = np.concatenate((node_redshifts, np.array([4.9, 5.7, 6.0, 6.6, 7.0, 7.3])))
   and the toml template provided in the attached file sgh25.txt (change file extension to .toml)
2. Read the average neutral fraction as a function of node redshift using lc.global_quantities['neutral_fraction'][::-1]
3. The neutral fraction does not increase monotonically with redshift. It exhibits discontinuities coinciding with the redshifts appended to the node redshift array.

Expected behavior:
Neutral fraction should increase monotonically with redshift.

Details:

• 21cmFAST main branch synced to commit f6bf7ed

sgh25.txt

[steven-murray]

Ah, this is a strange usage. Probably we're just not forcing the node_redshifts to be sorted in InputParameters when we should.

[steven-murray]

I take that back, we are forcing node_redshifts to be decreasing.

[steven-murray]

MWE from @samgagnon:

import numpy as np
import py21cmfast as p21c

import logging, os

#    cache = p21c.OutputCache(args.cache)
cacheconfig = p21c.CacheConfig.off()
cache.list_datasets(
    kind=None,
    inputs=None,
    redshift=None,
    all_seeds=True,
)
p21c.config['HALO_CATALOG_MEM_FACTOR'] = 2.0

min_redshift=5.0
max_redshift=35.0

node_redshifts = p21c.get_logspaced_redshifts(
  min_redshift  = min_redshift,
  max_redshift  = max_redshift,
  z_step_factor = 1.02
)
  
node_redshifts = np.concatenate((node_redshifts, np.array([4.9, 5.7, 6.0, 6.6, 7.0, 7.3])))

inputs = p21c.InputParameters.from_template('./sgh25.toml',\
            random_seed=1,\
            node_redshifts=node_redshifts).evolve_input_structs(N_THREADS=4)

lcn = p21c.RectilinearLightconer.with_equal_cdist_slices(
    min_redshift=inputs.node_redshifts[-1],
    max_redshift=inputs.node_redshifts[0],
    resolution=inputs.simulation_options.cell_size,
    cosmo=inputs.cosmo_params.cosmo,
    quantities=['brightness_temp', 'density', \
                'kinetic_temperature', 'neutral_fraction', 'velocity_z']
)

_, _, _, lc = p21c.run_lightcone(
    lightconer=lcn,
    inputs=inputs,
    write=cacheconfig,
    cache=cache
)