Rebased splice potential

[AnatoleStorck]

Here is my code to spliced the potential, rebased onto the current master branch.

[apontzen]

Thanks, this looks great. I've added some very cosmetic comments, but if you could run over these that'd be much appreciated.

Also, could you add a test potential splicing scenario, so that we catch any later regression? e.g. you could make a copy of and adapt this test: https://github.com/pynbody/genetIC/tree/master/genetIC/tests/test_24_splicing

[apontzen]

Thanks, changes look good - are you OK to add a test as described above?

[AnatoleStorck]

Absolutely! I will work on this sometime in the next few days.

…

On Mon, Sep 30, 2024 at 20:14 Andrew Pontzen ***@***.***> wrote: Thanks, changes look good - are you OK to add a test as described above? — Reply to this email directly, view it on GitHub <#127 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AUELSFNWVJCJZKQUSWNBEEDZZGIHHAVCNFSM6AAAAABPDH634OVHI2DSMVQWIX3LMV43OSLTON2WKQ3PNVWWK3TUHMZDGOBTHA3DKMRQHE> . You are receiving this because you authored the thread.Message ID: ***@***.***>

[apontzen]

Shouldn't the result be deterministic to within a certain numerical tolerance? If so, I think we should adjust the tolerance rather than force it to run on a single thread (which is perhaps not a realistic test).

[AnatoleStorck]

I tried this up to a relative tolerance of 1e-7, with 93% elements matching at that point. The problem is it needs to do 1e4 iterations so when the tester outputs the log, its all there is in the console. I will test and set it to 1e-8 which should match all elements, but outputting the log will have to be done differently.

[apontzen]

I'm not sure I entirely follow - for this test, can't we just ignore what's in the log and examine the final result?

[AnatoleStorck]

So instead of trying to set arbitrarily lower tolerances (which were never reached as the minimizer reached 10*res**3 iterations anyways), I've reduced the splicing region to 5 Mpc/h on a 50 Mpc/h box. This makes the number of iterations needed to reach a standard tolerance of 1e-5 around 15, which is not enough for nondeterministic deviations to accumulate. Not only that, but a smaller splicing regions is more in line with the current science case of potential splicing, dark matter halo Lagrangian patches.

[apontzen]

Perfect, thanks.

Sorry for asking one more thing but I just realised that all the existing method names in the ICGenerator class (with one mistaken exception) are camelCased whereas you have used some_underscores -- would you mind renaming your methods to be consistently camel-cased?

I promise this is the last thing then we can merge.

[apontzen]

This seems ready to go. To understand for the documentation -- why do you recommend ZELDOVICH_GRADIENT_FOURIER_SPACE is disabled when splicing potential?

[AnatoleStorck]

The decision to disable ZELDOVICH_GRADIENT_FOURIER_SPACE was made pretty early on in the development of potential splicing by @cphyc . I think it came down to it being more advantageous to calculate the displacement field in real space using finite differences than in Fourier space. Corentin can probably motivate the decision better than I can.

[AnatoleStorck]

Perhaps adding a warning in the splicePotential function like in line 325 of linearmodification.hpp for compiling with ZELDOVICH_GRADIENT_FOURIER_SPACE could be helpful.

[cphyc]

The motivation to disable ZELDOVICH_GRADIENT_FOURIER_SPACE is to force the velocities (and hence the density) to be obtained using a local finite-difference approach. That way, once the potential is set up to machine precision in the spliced region, so are the velocities and densities (minus a few pixels on the edge).
If gradients are computed in Fourier space, the differencing isn't local anymore and some information from the outside can leak into the spliced region when going from potential to velocities.

[apontzen]

That makes sense, but I don't think it will be a huge issue for most potential users (unless you have tests that show otherwise)?

I'd suggest we update the warning message to be something more like

You are splicing potential with ZELDOVICH_GRADIENT_FOURIER_SPACE. Due to non-locality, the density field may not exactly match within the spliced region. Consider recompiling with ZELDOVICH_GRADIENT_FOURIER_SPACE disabled if precise density field matching is important to you.

appropriately split across lines...

Does that sounds reasonable?

[AnatoleStorck]

Well if you are interested in keeping all local linear information about the evolution of your patch the same, then it would be important. Otherwise you wouldn't know if the effect on your evolved patch is from the environment or on the changes made in the initial velocity and density fields while taking the gradient in fourier space. Its not something Corentin and I have quantified but maybe doing a quick check in my free time might be possible. In any case, I'm pretty happy with the changes and suggestions made, thank you Andrew!

[cphyc]

This looks good to me! Thanks @AnatoleStorck for the work! Once https://arxiv.org/abs/2409.13010 is accepted and this PR merged, it'd be great to add this to the manual. @apontzen what is the procedure for this?

[apontzen]

@cphyc Agree we need to document it.

At the moment the manual is very 'manual'; it's just a latex document.

Ideally we should turn it into markdown or something like that and host it in a repository that uses Jekyll or something, hosted in this repository. Do you feel like doing that? Or would you prefer just to edit a latex file and send back to me?

[apontzen]

@cphyc @AnatoleStorck quick ping re the documentation discussion above - we can add something very basic if you like, but we should definitely describe the basics at a minimum. Let me know! It will be good to get this merged and released shortly.