WIP: Merge black hole (and potentially some other) updates from Yannick's development branch into master

This merge request integrates the changes made in my (primarily) black hole development branch. There are some more experimental features in there that we likely do not want to include, so feel free to pick (and adapt) what seems appropriate (a copy of this development branch is kept at yb-bhdev).

The main elements are:

• Options to fine-adjust many elements of the black hole model via YAML parameters
• Option to use black hole subgrid masses from (idealised) initial conditions
• Option to use "multi-phase" black hole accretion approach
• Option to prescribe a repositioning velocity for black holes
• Slight re-design of the black hole // gas interactions, so that gas neighbour velocities are always recorded in the black hole frame
• Some extra black hole outputs (which, as discussed previously, we may not want to take over to master)
• Crude option to use different levels of output (e.g. snap- and snipshots), and modify the output list upon restarting from restart files
• Option to easily block entire particle types from outputs
• Some (hopefully) clarifications and improvements to the comments

Also (added late):

• Changed the black hole output 'CumulativeNumberSeeds' to 'CumulativeNumberOfSeeds'. I kept getting this one wrong, since it is not quite consistent with the general convention of using unabbreviated, descriptive names. Feel free to ignore if there was a good reason for the original name.

So that branch does not compile on cosma or my laptop.

I think we want to keep the multi-level output out of it for now. The full mechanism is processing rapidly.

Can you tell me what error messages you get on compilation? It seems to run fine for me on cosma or my laptop.

I get many hundred lines of errors. It seems OK when I compile with the full EAGLE model, but if I switch that off (for instance for DMO runs) then it all breaks. For instance ./configure breaks.

I can deal with it and extract the changes that make sense in master if that's helpful.

added 1 commit

• fc320191 - Fixed incompatible function signatures in non-EAGLE black hole model

Compare with previous version

Ah, I see what the problem was -- I had changed a few of the black hole related function signatures, but hadn't mirrored these changes in the 'Default' black hole model (I had never used this branch for runs without EAGLE, which is why I hadn't noticed it before).

If it's straightforward, could you please check whether it compiles correctly now?

As for selecting what (if anything) should go into master, feel free to either pick the relevant bits out yourself or let me know what you want and I can trim the merge branch to that. Whatever is easiest.

Ah good. All works now. I'll have a more detailed look and probably fish out directly the changes that do not conflict with the broader snipshot update.

added 1 commit

• 50e94d1b - Applied code formatting tool

Compare with previous version

changed the description

I have done a scan of the changes. Here are a few questions for you:

• In the multi-phase Bondi, why use mj * wi and not pj->rho as the density of the particle entering the Bondi formula?

• In the multi-phase Bondi, is there ever a scenario where the denominator can be 0? Or is the check there only for debugging purposes?

• In black_holes_prepare_feedback(), shouldn't we check that bh_props->multi_phase_bondi is set and not just that the accretion rate is > 0?

• When recording the number of reposition attempts, shouldn't that be done only for BHs with a mass below the max mass for reposition?

• Is it more readable if the repositioning based on upsilon is activated via an on/off switch and then upsilon is used? I think I'd make the parameter files more human readable.

• The number of time-steps should not be increased in black_holes_init_bpart(). This function can be called multiple times per time-step if a BH did not find the appropriate number of neighbours at the first attempt.

• Could you update the parameter file of the EAGLE_ICs/EAGLE_25 case with the new fields set such that we reproduce the "old" behaviour in this case? I can then cross-check that as well

added 2 commits

• cf75f878 - Corrected counting of BH time steps, improved multi-phase branch
• f535055d - Merge branch 'yb-bhdev' into yb-bhdev-merge

Compare with previous version

• In the multi-phase Bondi, why use mj * wi and not pj->rho as the density: --> My reasoning here is that we want the (weighted) contribution of the particle to the gas density at the black hole (mj * wi), rather than the (already smoothed) density at the location of the gas particle, which would then (presumably) have to be weighted by another factor of wi to account for the distance from the black hole. So I think this is closer to the intention behind the Bondi formula. Let me know if you think I'm wrong.

• Zero denominator in multi-phase Bondi: --> no, I don't think this should ever happen, as in the 'standard' case (even if, by chance, the velocity difference is zero, the sound speed should always be strictly positive). It's a sanity check.

• Check for bh_props->multi_phase_bondi in black_holes_prepare_feedback: --> definitely, changed already. For some reason, when I wrote this I seem to have been under the impression that the props were not available here, hence the workaround. Thanks for spotting this.

• Number of reposition attempts: --> I had done it deliberately so that it would count all possible cases, including BHs above the mass limit (which is easy to spot in the analysis anyway). The idea was to get a feel for how commonly there were any eligible gas particles within the BH kernel.

• Upsilon-based repositioning switch: --> you're probably right, it was just the lazy option. I'll change that, also in the analogous case of setting the velocity-to-sound-speed threshold ratio.

• Good catch about the black_holes_init_bpart() -- I've moved the counter to black_holes_prepare_feedback() now. This may well have been the reason why I had so far seen that black holes even attempted repositioning in rarely more than ~30% of "time steps"...

• I'll update the EAGLE_25 parameter file together with the two switches discussed above, good point.

Great, thanks. Sorry for the level of scrutiny. I am trying to minimize the number of times we re-run the first wave of calibration...

Let me think about the rho in the formula. At first sight, I am not convinced the current version is the ideal choice but I may be wrong. In any case, the difference should be small. I'll upload some notes when I am done.

added 2 commits

• 3d16084d - Small updates to the black hole routines (incl. minimum velocity threshold…
• 4fa75ff7 - Merge branch 'yb-bhdev-merge' of https://gitlab.cosma.dur.ac.uk/swift/swiftsim into yb-bhdev-merge

Compare with previous version

No problem, thanks for checking it carefully.

I've updated the black hole parameters now, including two explicit switches for enabling a velocity cut in repositioning (with_reposition_velocity_threshold), and for enabling the prescribed-speed repositioning scheme (set_reposition_speed). I've also fixed a minor bug in that the option to set a fixed minimum velocity threshold appears to not have actually been used so far in this branch.

Both the 'main' parameter example and EAGLE_25 parameter files are updated. I have slightly re-arranged the ordering of the parameters to correspond roughly to the order/grouping in which they are used. In the EAGLE_25 YAML file, I omitted the three non-applicable parameters (use_subbgrid_mass_from_ics, reposition_coefficient_upsilon, and reposition_exponent_xi).

Note that, for idealised simulations with BHs in the ICs, I now by default read BH subgrid masses, rather than initialising them to the particle masses. I think this is the better approach in terms of what people would expect, but if you think overwriting them by default is better, let me know.

So, I now think I agree with your expression:

I agree, reading the subgird masses is better if they are present in the ICs.

I'll have a second look at the updates later tonight.

Good, glad you think my approach to the density is correct. Thanks for checking it!

About the subgrid masses: I guess the one downside of this is that there may not be a SubgridMasses data set in the ICs, in which case the subgrid masses are probably initialised to garbage. My take on this would be that anyone using black hole ICs that do not contain subgrid masses and not instructing the code to initialise them to the particle masses deserves what they get... Maybe a warning in this case would be best, provided there is a way to catch the situation that an optional input field was not found?

No hurry on the updates, whenever you have time is good.

Now that we know these changes are not problematic, I'll have another look.

Was there something else you wanted to add or change?

The only thing I could think of was a test for whether the Subgrid masses are actually present in the ICs when the option to initialise them to particle masses is not set in the YAML file.

One crude check would be to test whether the subgrid masses are all > 0, which could be done in black_holes_first_init_bpart(). I cannot immediately think of any situation when one would want to have a zero-subgrid-mass black hole, but maybe you can?