Binding energies and Vmax calculations #132
Merged
Conversation
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
robjmcgibbon
changed the title
Collaborator
Author
|
This is ready for review, In this PR I am reading in the potential energies from HBT, and calculating the kinetic energies in SOAP. Unfortunately this leads to some particles with PE + KE > 0. I think this is due to two differences
If I look at the value of KE_HBT / KE_SOAP for the least bound particle in each subhalo then it is always close to 1 as we would expect. I have tested this on the L0025N0188, and get the following figure. |
Collaborator
Is there a reason why we want to compute kinetic energies in the most bound particle reference frame? Perhaps we can make SOAP and HBT-HERONS consistent with each other. |
robjmcgibbon
added a commit
that referenced
this pull request
Jul 17, 2025
* Record processing time for each halo * Skip density check for satellites * Use written chunks when restarting * Store calculation time for each property * Use lookup table for stellar ages * Add links to footnotes * Remove unneeded sort * Squash sort bug * Non iterative inertia tensors * Catch io errors for chunk files * Move filter to parameter file * Remove VR FOFSubhalo halo type * Remove GroupNr_all * Descriptions for timing information * Sort halos on chunk * Report peak memory usage * Copy values from smaller apertures * Store filter metadata for properties * Fix docs generation * Add script for plotting timings * Added wrong plot_time script * Update plot time script * Plot setup time * Various minor fixes * Wrap COM values * Add Victor symmetric array script * Set min read radius to be optional * Matthieu documentation suggestions * Change assigning of files to ranks * Add MPI compression script * Set compressed properties to be read-only * Add vmax radius to SO properties * Fix velocity scale factors * First attempt at reading in EAGLE * Calculate cell offsets * Working metadata * Working EAGLE conversion script * Add L100 EAGLE * Add license * Add JOSS paper * Fix tests * Create SOAP package * Run formatter * Remove stray file * Update formatter version * Format all files * Minor paper compilation updates * Address pre-review comments (#133) * Use pyproject.toml instead of setup.py * Update packages * List packages * Find packages automatically * Add small DMO example * Test subhalo rank with HBT data * Add scripts for running HBT/VR * Make test helper script MPI-compatible * Generate documentation for small test * Move Rockstar test and read with numpy instead of pandas * Add rockstar test * Run formatter * Move subfind test * Update scripts * Include neutrinos in omega_m * Update dust fraction calculation (#134) * Don't use hydrogen fraction for dust properties * Typo footnote_dust.tex * Update parameter file README (#136) * Update parameter file README * Mention pdf doco * Check each halo type for valid properties (#137) * Have central's HostHaloIndex point to themselves (#138) * Allow for different alias for snipshots (#139) * Adding a number of extra properties (#141) * Add dust properties * Add stellar CoM * Update test * Minor name update * Remove references to centre of potential (#143) * Calculate SOAP index of progenitor/descendant (#142) * Small refactor * Don't update skipped SOAP properties * Format * Calculate and store progenitor/descendant indices * Update handling of named columns * Format * Add footnote * Calculate fof radii (#144) * Script v1 * Create files across ranks * Update docs * Create virtual file * Use relative paths * Format * Add test * Use absolute positions * Proper no wrap * FOF radii * Update SOAP to read new values * Add half-light radii for each GAMA band. (#145) * Defined new half light radii property * Copied get_half_mass_radius to start implementing luminosity weighting. * Create general function to find radius enclosing half of a given quantity * Defined half light radius. * Fix a couple of bugs. * Script to obtain the SOAP property evolution for specified HBT-HERONS subhaloes (#146) * Add script * Fix redshift saving in HDF5 * Update misc README.md * Make the script load the TrackIds and properties from a text file. * Fixing typos * Skip assert check for ranks with no subhalos (#148) * Luminosity-weighted angular momenta and inertia tensors (#135) * Add shorthand to obtain grik luminosities. * Add basic functions to measure disc to total luminosity in different bands. * Make new properties not be computed by default * Add call of unimplemented function that computes light weighted L quantities. * Define function for luminosity weighted L quantities. Not implemented for now. * Add docs on parameters for get_angular_momentum_and_kappa_corot_luminosity_weighted * Make kappa_corot and counterotating mass give a value per GAMMA band * Add counterrotating luminosity array * Handle multiple luminosity bands when computing average angular momentum * Compute projected angular momentum if we do_counterrot_luminosity * We now compute if stars are counterrotating or not for each luminosity band. Expanded description of unyt array shapes to keep track in each step. * fix typo * Basic implementation done. * Add correct return and remove NotImplemented error * Change variables assigned to subhaloes to prevent overwritting mass-weighted calculations * Only return a single property, which includes all GAMMA bands. * Return a mass and light-based D/T. * Add kappa_corot for luminosity bands * change variable name * Add property to save all luminosity-weighted angular momentum * Add property definition in the Aperture properties * Added new properties to the property table * Add properties to the subhalo properties parameter file of COLIBRE * Remove lazy properties that were not needed anymore * Update docstring * Add luminosity-weighted angular momenta to bound subhalo properties * Add properties to the ApertureProperties parameter file section * Change output expected type * Change condition to execute new properties. Mainly because of an annoying inconsistent between how starless subhaloes are handled. Mstar = 0 but Luminosities = None. * Fix some indicies and if condition * Add important clarification * Fix bug in kappa_corot definition. * Formatting alignment * Add luminosity-weighted quantities to SO apertures * Clarification comment * Individual angular momentum for particles is now correct for kinetic energy calculation * Fix GAMA typo * Remove mass term form average angular momentum we keep it in the individual particle angular momentum to obtain the kinetic energy in rotation. * Create copy of inertia tensor to implement luminosity weighting * Update description and inputs of inertia tensor * Mask luminosities if we mask other properties * Add comment... * Create correct shape for eigenval and eigenvec matricies * Shapes for each luminosity band * Add stopping criteria on a per band basis * Define the luminosity-weighted inertia tensor for BoundSubhalos * Compute particle positions for each band. Perhaps unneccessary since they will have the same positions in every band... Check how to improve later on. * Update break condition to be done on a luminosity band basis * Assign luminosity-based weights * is_converged is now a boolean array * Remove comment that was already addressed * Tensor for each band done. Still need to check if it gives me the correct element ordering relative to the default calculation of the inertia tensor. * Re-add removed calculation of inertia tensor eigenvalue/vectors * Create variable to hold the number of provided luminosity bands * Return flattened matricies for each band * Remove bug in stopping criteria for iterative tensor * Added missing np.abs() call It was causing the new function to stop iterating early * Add import of new function * Fix to prevent overflows into negative values for eigenvectors of inertia tensor * Defined luminosity weighted projected inertia tensor * Add call to luminosity-weighted function (undefined yet) * Define luminosity weighted function and pass luminosity arrays. * Correct typo from main branch * Finish implementation of get_projected_inertia_tensor_luminosity_weighted. Directly taken from 3D inertia calculations. Will likely need some restructuring. * Fix incorrect comment statement * Bug: using 3D positions instead of projected... * Change location of axis ratio transpose * Handle q = 0 values when testing convergence. * Handle 0 / 0 division. * Make q be an array from the get go * Only compute tensors for non converged bands * Do not rely on disabling warnings to handle q = 0 cases. * Make q = 0 tensors to have all matrix entries to be equal to 0. * Formatting * Update variable names and parameter files * Make get_inertia_tensor_luminosity_weighted only do computations for non-converged bands More similar to get_projected_inertia_tensor_luminosity_weighted * Fix issue resulting from merging * Make 3D inertia tensors general * Fix typo en error message * Make mass-weighted inertia tensor more self-evident * Make projected inertia tensors general * Two bugs: incorrect function call and units in reduced tensors * Move inertia tensor functions into their own file * Fix unyt array initialisation to specify unit registry * Generalise kappa_corot * Shorten the description of helper functions * Update colibre test * Reformat * Add note on single particle * Remove questions --------- Co-authored-by: robjmcgibbon <robjmcgibbon@gmail.com> * Full SOAP EAGLE catalogues (#150) * Add more EAGLE snapshot properties * Add script to approximate H fractions * Use proper index for SOAP EAGLE * Update scripts * Squash bugs * Define apertures based on BoundSubhalo properties (#147) * Add note on timing flags * Allow for property defined apertures * Add error messages if required radii are disabled * Missing method * Update pdf documentation * Add massive halos for testing * Format * Add contribution guide (#151) * JOSS - Stuart comments 1 (#154) * Merge fof position test * Update docs * Add project scripts & fix tests * Format * Add venv script for strw * Rerun pdf generation * Update README.md * Match subhalos using membership files (#156) * Backbone broken * Fix test in calculate_fof_pos * Run matching * Match to SOAP idx * Check matches are consistent * Output to file * Squash bug when saving * Update docs * Centrals only * Update README * Add __main__ * Run formatter * Require astropy >= 6 * Update README * Generate pdf with all properties * Set lower limit on dust mass fractions (#158) * Set lower limit * Format * Add DM CoM (#159) * Add sbatch script for calculating FoF radii (#162) * Add sbatch script for FOF radii * Skip assert for ranks with no halos * Nicer check * Read new names if present (#163) * Add workflow to check code is formatted (#165) * Add format workflow * Add --check * Run formatter * COLIBRE core excision (#167) * Exclude COLIBRE core * HalfMassRadii H * Binding energies and Vmax calculations (#132) * Read binding energy from HBT * Output TotalBindingEnergy * Vmax for aperture properties * Squash bug and format * Stop ranks hanging when reading binding energy * Read SnapshotIndexOfLastIsolation * Read in potential energy instead of binding * Remove duplicate line * Use vCoM from all particles * Move subhalo_properties check for Ntot * Update COLIBRE test * Remove Ntot warning * Add documentation of SpecificPotentialEnergies * Add parameter file option to read energy * Review changes * Editorial changes (#171) * Editorial changes. * Update paper.bib * A couple of more changes * Update title capitalization * Set read_potential_energy default to false * Add more references * Read sorted HBT (#172) * Add script to recalculate COLIBRE FOF (#166) * Add sbatch script for FOF radii * Skip assert for ranks with no halos * Nicer check * More logging * Add option to copy datasets * Add node requirements --------- Co-authored-by: Victor J. Forouhar Moreno <60343271+VictorForouhar@users.noreply.github.com> Co-authored-by: Ivelina Momcheva <ivelina.momcheva@gmail.com>
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
Potential energy
Depends on SWIFTSIM/HBT-HERONS#76. We want to output the potential energy of each subhalo. HBT now outputs the potential energy of each particle, so I have updated the group membership script to read in those values. In subhalo properties we then sum up those values to get the total binding energy of the subhalo.
TODO
Other
SnapshotIndexOfLastIsolation(lifted directly from the HBT catalogues)subhalo_properties.pyto make sure we have loaded all particles. I've moved this to a slightly different location.