/******************************************************************************* * This file is part of SWIFT. * Copyright (c) 2020 Josh Borrow (josh@joshborrow.com) * * This program is free software: you can redistribute it and/or modify * it under the terms of the GNU Lesser General Public License as published * by the Free Software Foundation, either version 3 of the License, or * (at your option) any later version. * * This program is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * GNU General Public License for more details. * * You should have received a copy of the GNU Lesser General Public License * along with this program. If not, see . * ******************************************************************************/ #ifndef SWIFT_PARTICLE_SPLITTING_H #define SWIFT_PARTICLE_SPLITTING_H #include "inline.h" #include "io_properties.h" #include "particle_splitting_struct.h" #include /** * @brief Initialise a particle_splitting_data struct * at the start of a run, given an initial * progenitor ID. * * @param splitting_data the uninitialised particle struct. * @param id the ID of the particle (used for progenitor_id) */ __attribute__((always_inline)) INLINE static void particle_splitting_mark_part_as_not_split( struct particle_splitting_data* restrict splitting_data, long long id) { splitting_data->progenitor_id = id; splitting_data->split_tree = 0; splitting_data->split_count = 0; } /** * @brief Updates the binary trees of the particles that * have been split. pi should retain its original ID, * and hence gets the relevant part of the tree set * to zero. * * @param sdi first particle_splitting_data* resulting from * the splitting event. * @param sdj second particle_splitting_data* resulting from * the splitting event. * @param id_i ID of the first particle. * @param id_j ID of the first particle. * @param extra_split_logger File pointer (opened) where to log * the resetting of the split counters. */ __attribute__((always_inline)) INLINE static void particle_splitting_update_binary_tree( struct particle_splitting_data* restrict sdi, struct particle_splitting_data* restrict sdj, const long long id_i, const long long id_j, FILE* extra_split_logger, swift_lock_type* file_lock) { /* Print warnings if we have split these particles more * than the number of times the tree can accommodate. * Warning is only printed once for each particle */ if (sdi->split_count > 0 && sdi->split_count % (8 * sizeof(sdi->split_tree)) == 0) { message( "Warning: Particle (%lld) with progenitor ID %lld with binary tree " "%lld has been split over the maximum %zu times, making its binary " "tree invalid.", id_i, sdi->progenitor_id, sdi->split_tree, sizeof(sdi->split_tree)); /* Are we logging this? */ if (extra_split_logger != NULL) { /* Log the old state before reseting. Use the lock to prevent multiple * threads from writing at the same time. */ lock_lock(file_lock); fprintf(extra_split_logger, " %12d %20lld %20lld %20d %20lld\n", engine_current_step, id_i, sdi->progenitor_id, sdi->split_count, sdi->split_tree); fflush(extra_split_logger); /* Release the lock and continue in parallel */ if (lock_unlock(file_lock) != 0) error("Impossible to unlock particle splitting"); } /* Reset both counters and trees */ sdi->split_tree = 0LL; sdj->split_tree = 0LL; /* Set both particles as having particle i as their progenitor */ sdi->progenitor_id = id_i; sdj->progenitor_id = id_i; } /* Update the binary tree */ sdj->split_tree |= 1LL << sdj->split_count % (8 * sizeof(sdi->split_tree)); /* Increase counters on both; sdi implicitly has a zero * in the relevant spot in its binary tree */ sdj->split_count++; sdi->split_count++; } /** * @brief Specifies the particle-splitting related fields * to write to a dataset. * * @param parts The particle array. * @param xparts The extra particle array. * @param list The list of i/o properties to write. * @param with_cosmology Are we running with cosmology? * * @return Returns the number of fields to write. */ INLINE static int particle_splitting_write_particles(const struct part* parts, const struct xpart* xparts, struct io_props* list, const int with_cosmology) { list[0] = io_make_output_field( "ProgenitorParticleIDs", LONGLONG, 1, UNIT_CONV_NO_UNITS, 0.f, xparts, split_data.progenitor_id, "ID of the progenitor of this particle. If this particle is the result " "of one (or many) splitting events, this ID corresponds to the ID of the " "particle in the initial conditions that its lineage can be traced back " "to. If the particle was never split, this is the same as ParticleIDs."); list[1] = io_make_output_field( "SplitCounts", UINT8, 1, UNIT_CONV_NO_UNITS, 0.f, xparts, split_data.split_count, "Number of times this particle has been split. Note that both particles " "that take part in the splitting have their counter incremented."); list[2] = io_make_output_field( "SplitTrees", LONGLONG, 1, UNIT_CONV_NO_UNITS, 0.f, xparts, split_data.split_tree, "Binary tree describing splitting events. Particles that keep the " "original ID have a value of zero in a splitting event, whereas" "particles given a new ID have a value of one."); return 3; } /** * @brief Specifies which star particle fields to write to a dataset * * @param sparts The star particle array. * @param list The list of i/o properties to write. * * @return Returns the number of fields to write. */ INLINE static int particle_splitting_write_sparticles( const struct spart* sparts, struct io_props* list) { list[0] = io_make_output_field( "ProgenitorParticleIDs", LONGLONG, 1, UNIT_CONV_NO_UNITS, 0.f, sparts, split_data.progenitor_id, "Progenitor ID of the gas particle that became this star. If this " "particle is the result of one (or many) splitting events, this ID " "corresponds to the ID of the particle in the initial conditions that " "its lineage can be traced back to. If the particle was never split, " "this is the same as ParticleIDs."); list[1] = io_make_output_field( "SplitCounts", UINT8, 1, UNIT_CONV_NO_UNITS, 0.f, sparts, split_data.split_count, "Number of times the gas particle that turned into this star particle " "was split. Note that both particles that take part in the splitting " "have their counter incremented."); list[2] = io_make_output_field( "SplitTrees", LONGLONG, 1, UNIT_CONV_NO_UNITS, 0.f, sparts, split_data.split_tree, "Binary tree describing splitting events. Particles that keep the " "original ID have a value of zero in a splitting event, whereas" "particles given a new ID have a value of one."); return 3; } /** * @brief Specifies which black hole particle fields to write to a dataset * * @param bparts The black hole particle array. * @param list The list of i/o properties to write. * * @return Returns the number of fields to write. */ INLINE static int particle_splitting_write_bparticles( const struct bpart* bparts, struct io_props* list) { list[0] = io_make_output_field( "ProgenitorParticleIDs", LONGLONG, 1, UNIT_CONV_NO_UNITS, 0.f, bparts, split_data.progenitor_id, "Progenitor ID of the gas particle that became the seed BH. If this " "particle is the result of one (or many) splitting events, this ID " "corresponds to the ID of the particle in the initial conditions that " "its lineage can be traced back to. If the particle was never split, " "this is the same as ParticleIDs."); list[1] = io_make_output_field( "SplitCounts", UINT8, 1, UNIT_CONV_NO_UNITS, 0.f, bparts, split_data.split_count, "Number of times the gas particle that became this BH seed " "was split. Note that both particles that take part in the splitting " "have their counter incremented."); list[2] = io_make_output_field( "SplitTrees", LONGLONG, 1, UNIT_CONV_NO_UNITS, 0.f, bparts, split_data.split_tree, "Binary tree describing splitting events prior to BH seeding. Particles " "that keep the original ID have a value of zero in a splitting event, " "whereas particles given a new ID have a value of one."); return 3; } /** * @brief Write particle splitting data to the stdout for debugging purposes. * * @param p Particle data. * @param xp Extra particle data. */ __attribute__((always_inline)) INLINE static void particle_splitting_debug_particle(const struct part* p, const struct xpart* xp) { if (xp != NULL) { warning("[PID%lld] particle_splitting_data:", p->id); warning( "[PID%lld] progenitor_id = %lli, split_tree = %lli, " "split_count = %hhu", p->id, xp->split_data.progenitor_id, xp->split_data.split_tree, xp->split_data.split_count); } } #endif /* SWIFT_PARTICLE_SPLITTING_H */