Ionbox refactor #427

daviesje · 2024-09-16T11:40:38Z

This is a coarse modularisation of the ionised box calculation which separates the large function into a number of smaller steps. This should make future development easier. I have tested that we get the same results under a few variations of parameters, and it has not slowed down.

Key differences:

constants are calculated at the beginning and held in a struct which is passed into each function.
Radii are explicitly calculated at the beginning of the run, instead of on-the-fly in a while loop
the RNG for partial ionisations has been explicitly removed (it only affected the !USE_MASS_DEPENDENT_ZETA case)

Further modularisation is possible, but best left to another PR

…stant mfp

…into ionbox-refactor

codecov · 2024-11-11T13:19:06Z

Codecov Report

Attention: Patch coverage is 62.06897% with 11 lines in your changes missing coverage. Please review.

Project coverage is 79.51%. Comparing base (4071f4e) to head (246d89c).
Report is 38 commits behind head on v4-prep.

Files with missing lines	Patch %	Lines
src/py21cmfast/wrapper/structs.py	20.00%	8 Missing ⚠️
src/py21cmfast/wrapper/globals.py	33.33%	1 Missing and 1 partial ⚠️
src/py21cmfast/drivers/single_field.py	80.00%	1 Missing ⚠️

Additional details and impacted files

@@             Coverage Diff             @@
##           v4-prep     #427      +/-   ##
===========================================
+ Coverage    74.11%   79.51%   +5.40%     
===========================================
  Files           23       23              
  Lines         3747     3764      +17     
  Branches       641      644       +3     
===========================================
+ Hits          2777     2993     +216     
+ Misses         742      549     -193     
+ Partials       228      222       -6

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

steven-murray

Nice work @daviesje, this is going to make things much easier!

steven-murray · 2024-11-12T11:24:16Z

src/py21cmfast/src/InitialConditions.c

+//  This should be very rare, and even when it occurs, it will usually use a different initialiser,
+//  as well as be applied to different cells/halos, so I can't forsee any issues.
+//Just in case I'm not using it for the initial conditions, which is okay since the slower version is only used once
+void seed_rng_threads_fast(gsl_rng * rng_arr[], unsigned long long int seed){


Let's add a line that just checks if any of the seeds are the same and re-samples if so.

steven-murray · 2024-11-12T11:26:31Z

src/py21cmfast/src/InitialConditions.h

+//TODO: these seeding functions should probably be somewhere else
+//  Possibly make an rng.c/h


I agree with this

steven-murray · 2024-11-12T11:30:06Z

src/py21cmfast/src/Stochasticity.c

@@ -299,6 +299,46 @@ int stoc_halo_sample(struct HaloSamplingConstants *hs_constants, gsl_rng * rng,
    return 0;
 }

+//same as stoc_halo_sample but with a mass tolerance for comparison with Ivan
+int stoc_halo_sample_tol(struct HaloSamplingConstants *hs_constants, gsl_rng * rng, int *n_halo_out, float *M_out){


I think this function needs to be documented a bit to say what it does and why it's here.

steven-murray · 2024-11-12T11:36:40Z

src/py21cmfast/src/IonisationBox.c

+//Making a dummy previous box which has the required fields for the first snapshot:
+void setup_first_z_prevbox(IonizedBox *previous_ionize_box, PerturbedField *previous_perturb, int n_radii){
+    LOG_DEBUG("first redshift, do some initialization");
+    int i,j,k;
+    unsigned long long int ct;

-#pragma omp parallel shared(box) private(ct) num_threads(user_params->N_THREADS)
+    //z_re_box is used in all cases
+    previous_ionize_box->z_re_box    = (float *) calloc(HII_TOT_NUM_PIXELS, sizeof(float));
+    #pragma omp parallel shared(previous_ionize_box) private(i,j,k) num_threads(user_params_global->N_THREADS)
    {
-#pragma omp for
-        for (ct=0; ct<HII_TOT_NUM_PIXELS; ct++) {
-            box->z_re_box[ct] = -1.0;
+        #pragma omp for
+        for (i=0; i<user_params_global->HII_DIM; i++){
+            for (j=0; j<user_params_global->HII_DIM; j++){
+                for (k=0; k<HII_D_PARA; k++){
+                    previous_ionize_box->z_re_box[HII_R_INDEX(i, j, k)] = -1.0;
+                }
+            }
        }


is this not something that's better done from the python wrapper?

steven-murray · 2024-11-12T11:46:37Z

src/py21cmfast/src/IonisationBox.c

+    //TODO: moved from inside the R loop, either place verification in python wrapper or write a function to verify backend parameters
+    if (global_params.FIND_BUBBLE_ALGORITHM != 2 && global_params.FIND_BUBBLE_ALGORITHM != 1){ // center method
+        LOG_ERROR("Incorrect choice of find bubble algorithm: %i",
+                    global_params.FIND_BUBBLE_ALGORITHM);


This should be done in the GlobalOptions validator. I think we should be writing the C functions essentially without validation, assuming all has been validated externally.

steven-murray · 2024-11-12T12:30:07Z

src/py21cmfast/src/IonisationBox.c

+    //boxes which aren't guaranteed to have every element assigned to need to be initialised
+    //TODO: they should be zero'd in the wrapper, right?
+    if(flag_options->INHOMO_RECO) {
+        if(INIT_RECOMBINATIONS) {
+            init_MHR();
+            INIT_RECOMBINATIONS=0;
        }

-        // remember to add the factor of VOLUME/TOT_NUM_PIXELS when converting from
-        //  real space to k-space
-        // Note: we will leave off factor of VOLUME, in anticipation of the inverse FFT below
-#pragma omp parallel shared(deltax_unfiltered,xe_unfiltered,N_rec_unfiltered,prev_deltax_unfiltered,\
-                            log10_Mturnover_unfiltered,log10_Mturnover_MINI_unfiltered,stars_unfiltered) \
-                    private(ct) num_threads(user_params->N_THREADS)
+        #pragma omp parallel shared(box) private(ct) num_threads(user_params->N_THREADS)
        {
-#pragma omp for
-            for (ct=0; ct<HII_KSPACE_NUM_PIXELS; ct++){
-                deltax_unfiltered[ct] /= (HII_TOT_NUM_PIXELS+0.0);
-                if(flag_options->USE_TS_FLUCT) { xe_unfiltered[ct] /= (double)HII_TOT_NUM_PIXELS; }
-                if (recomb_filter_flag){ N_rec_unfiltered[ct] /= (double)HII_TOT_NUM_PIXELS; }
-                if(flag_options->USE_HALO_FIELD) {
-                    stars_unfiltered[ct] /= (double)HII_TOT_NUM_PIXELS;
-                    sfr_unfiltered[ct] /= (double)HII_TOT_NUM_PIXELS;
-                }
-                if(flag_options->USE_MINI_HALOS){
-                    prev_deltax_unfiltered[ct]          /= (HII_TOT_NUM_PIXELS+0.0);
-                    log10_Mturnover_unfiltered[ct]      /= (HII_TOT_NUM_PIXELS+0.0);
-                    log10_Mturnover_MINI_unfiltered[ct] /= (HII_TOT_NUM_PIXELS+0.0);
-                }
+            #pragma omp for
+            for (ct=0; ct<HII_TOT_NUM_PIXELS; ct++) {
+                box->Gamma12_box[ct] = 0.0;
+                box->MFP_box[ct] = 0.0;


Yes, this is all wrapper stuff.

steven-murray · 2024-11-12T12:31:51Z

src/py21cmfast/src/IonisationBox.c

+    // perform a very rudimentary check to see if we are underresolved and not using the linear approx
+    if ((user_params->BOX_LEN > user_params->DIM) && !(global_params.EVOLVE_DENSITY_LINEARLY)){
+        LOG_WARNING("Resolution is likely too low for accurate evolved density fields\n It Is recommended \
+                    that you either increase the resolution (DIM/BOX_LEN) or set the EVOLVE_DENSITY_LINEARLY flag to 1\n");
+    }


I feel like this also should be in the wrapper.

steven-murray

Looks great, thanks @daviesje !

steven-murray · 2024-11-13T15:10:46Z

@daviesje most recent changes look good -- I didn't understand the switch from true to false but I assume it was for a good reason. The increase in the halo buffer seems perhaps a little arbitrary, though not really a problem.

The two failing tests with nomdz -- these segfault here but pass locally? Can you run with valgrind or something?

daviesje · 2024-11-13T15:19:14Z

@daviesje most recent changes look good -- I didn't understand the switch from true to false but I assume it was for a good reason.

This was a bug introduced in the refactor where the wrong interpolation tables were being initialised (the previous redshift instead of the current).

The increase in the halo buffer seems perhaps a little arbitrary, though not really a problem.

This is entirely arbitrary, but doesn't effect the test results. I just need to make sure that it can hold all the halos that the tests produce, I thought about 500MB should be enough.

The two failing tests with nomdz -- these segfault here but pass locally? Can you run with valgrind or something?

This has just been fixed, it was an error on my part incorrectly setting interpolation table limits when USE_MASS_DEPENDENT_ZETA=False. For some reason running tests with -v made it not segfault, although the results would surely have been garbage

daviesje added 11 commits September 11, 2024 15:16

WIP filter tests

5abd202

WIP filters 0,3,4 test as expected

b7cd6a1

fix test comparison cases

23eb13c

fix bug in mfp filter, set expected mean

9efea6d

WIP refactor of ionbox

5dc7339

WIP installs

9e9c171

bugfixes

206ec93

make IC rng consistent with master, use faster version per snapshot

a7b797e

don't calculate log10turn in every cell unless we need to

2b3c4cf

WIP

08bcae9

remove partial ion rng model

67ecfd0

daviesje linked an issue Sep 16, 2024 that may be closed by this pull request

Refactor of ComputeIonizedBox #360

Open

daviesje added 9 commits September 26, 2024 10:25

Merge branch 'v4-prep' into ionbox-refactor

05eecd3

fix mcrit setting on minihalos

e7968de

fix minihalos?

b665ff1

add median scaling flag, set mfp filter to use filtered density & con…

a1c9df5

…stant mfp

reduce debug logs, fix recombination dz bug

d26be59

add tolerance sampler for Nikolic Comparison

b35cadb

add constant sfr for compatibility

3c56a38

comment cleanup

7d7ba1b

merge v4-prep

56b425b

daviesje mentioned this pull request Nov 4, 2024

Ternary Option for USE_INTERPOLATION_TABLES #437

Open

Merge branch 'v4-prep' into ionbox-refactor

6f0239e

daviesje requested a review from qyx268 November 7, 2024 17:49

daviesje added 5 commits November 8, 2024 16:35

fixing tolerance sampling

d148995

fix merge conflict typo

53ee836

Merge branch 'v4-prep' into ionbox-refactor

b9a4b90

fix merge typo

fdce42b

Merge branch 'ionbox-refactor' of https://github.com/21cmfast/21cmFAST …

083dee0

…into ionbox-refactor

steven-murray reviewed Nov 12, 2024

View reviewed changes

daviesje added 4 commits November 12, 2024 14:53

move rng to separate file, no duplicate seeds

3af6101

add new files

dd5d89f

move validations to wrapper

fb48534

no more re-allocating first ionbox, fix some bugs

6326f68

daviesje marked this pull request as ready for review November 12, 2024 18:02

steven-murray approved these changes Nov 12, 2024

View reviewed changes

daviesje added 4 commits November 12, 2024 21:20

fix binning and increase sampling on halo tests

5ba0de5

fix previous nion bug, make tests for segfaults

d584482

fix sampler tests

4e01ba6

rebalance halo sampler test buffer size

76d594b

fix ionbox fgtrm segfault

862cf3f

fix xray_source turnover bug, improve box summary printing

246d89c

steven-murray approved these changes Nov 13, 2024

View reviewed changes

daviesje merged commit d974de9 into v4-prep Nov 13, 2024
6 of 12 checks passed

daviesje deleted the ionbox-refactor branch November 13, 2024 17:27

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Ionbox refactor #427

Ionbox refactor #427

daviesje commented Sep 16, 2024

codecov bot commented Nov 11, 2024 •

edited

Loading

steven-murray left a comment

steven-murray Nov 12, 2024

steven-murray Nov 12, 2024

steven-murray Nov 12, 2024

steven-murray Nov 12, 2024

steven-murray Nov 12, 2024

steven-murray Nov 12, 2024

steven-murray Nov 12, 2024

steven-murray Nov 12, 2024

steven-murray left a comment

steven-murray commented Nov 13, 2024

daviesje commented Nov 13, 2024 •

edited

Loading

		//TODO: these seeding functions should probably be somewhere else
		// Possibly make an rng.c/h

Ionbox refactor #427

Ionbox refactor #427

Conversation

daviesje commented Sep 16, 2024

codecov bot commented Nov 11, 2024 • edited Loading

Codecov Report

steven-murray left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

steven-murray left a comment

Choose a reason for hiding this comment

steven-murray commented Nov 13, 2024

daviesje commented Nov 13, 2024 • edited Loading

codecov bot commented Nov 11, 2024 •

edited

Loading

daviesje commented Nov 13, 2024 •

edited

Loading