Fix the sorting of gas particles into rays used for launching AGN jets

As we were finalizing some changes to the code before being merging the AGN spin/jet model in (just before Flamingo calibration), @matthieu suggested dealing with the rare and special case where a BH particle runs out of gas particles in one of its hemispheres (relative to the spin axis). The relevant lines here are 339-346 in black_holes_iact.h:

      if (cosine_theta < 0) {
        quantity_to_minimize = -fabsf(cosine_theta) + ray_jet_correction;
        quantity_to_minimize_pos =
            -1e8 * fabsf(cosine_theta) + ray_jet_correction;
      } else {
        quantity_to_minimize = -1e8 * fabsf(cosine_theta) + ray_jet_correction;
        quantity_to_minimize_pos = -fabsf(cosine_theta) + ray_jet_correction;
      }

And lines 371-376 as well:

  /* Loop over rays and do the actual minimization */
  for (int i = 0; i < spinjet_blackhole_number_of_rays; i++) {
    ray_minimise_distance(quantity_to_minimize, bi->rays_jet, arr_size_jet,
                          gas_id, pj->mass);
    ray_minimise_distance(quantity_to_minimize_pos, bi->rays_jet_pos,
                          arr_size_jet, gas_id, pj->mass);

What these lines do is that they minimize (sort) the particles according to how close they are to the spin axis, separately for each hemisphere of the BH. In the ideal case we need only particles with cosine_theta < 0 to populate the ray bi->rays_jet (this is the ray that launches particles in the opposite direction from the BH spin vector), and particles with cosine_theta > 0 to populate bi->rays_jet_pos (equally this launches particles along the spin vector).

The extra lines involving the factor 1e8 I added so that each ray is also populated with particles from not just its own hemisphere, but the entire smoothing kernel (the other hemisphere too). However, the idea is to shove these particles at the very end of the ray and use them only if the main hemisphere is empty. This means multiplying whatever we are sorting (the cosine of the angle in this case) with a very small value, not a very large one. So it should be 1e-8 instead of 1e8, which is the change being implemented in this MR.

With the previous implementation, the particles from the wrong hemisphere were sorted to the very beginning of each of the rays. This means that particles were being launched from the wrong sides of the BH. Each time, two particles were being launched towards one another, and with the BH in between them... I have tested this in an isolated idealized case (with the spin vector pointing along the z-axis), by printing out the jet velocities and heights relative to the BH. Indeed in the current implementation particles launched along the positive z-axis are below the BH, while ones with the velocity pointing in the negative z direction were above the BH. The fix changes this around.

One can also see a pretty clear visual difference:

The top one is now with the fix implemented. The particles make two coherent streams/jets, even though we're talking about only a couple of dozen launching events. The bottom one is the incorrect version. The particles launched into one another get shocked immediately. The BH ends up being surrounded by lots of hot gas, which also suppresses the accretion rate (hence only 2 particle kicks happened here). The gas eventually (probably buoyantly) rises in a random direction.