Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[bugfix] H nuclei are being double-counted in the Arepo frontend #4211

Merged
merged 3 commits into from
Nov 17, 2022

Conversation

jzuhone
Copy link
Contributor

@jzuhone jzuhone commented Nov 16, 2022

PR Summary

yt computes "nuclei density" fields by summing up ion number density fields for the same nuclei with different ionization states. The Arepo frontend currently has "H", "H_p0", and "H_p1" all added to the species_names attribute in ArepoFieldInfo, which is then used to compute the set up the species fields and compute the fields for specific abundances.

The "H_p0" and "H_p1" do not belong in this list for two reasons:

  1. Their respective associated fields are handled later on in the method ArepoFieldInfo.setup_gas_particle_fields, and do not need to be handled by setup_species_fields.
  2. Putting them in species_names results in the number of hydrogen nuclei being double-counted

This PR removes those ionization species from the species_names attribute and also removes an alias between H_nuclei_density and H_number_density which was supposed to work (and maybe did at one point) but is no longer working. It is also unnecessary with the above fix.

PR Checklist

  • New features are documented, with docstrings and narrative docs
  • Adds a test for any bugs fixed. Adds tests for new features.

@jzuhone
Copy link
Contributor Author

jzuhone commented Nov 16, 2022

@matthewturk or @neutrinoceros this is a fix which should be merged pretty quickly if possible

@neutrinoceros neutrinoceros added this to the 4.1.2 milestone Nov 16, 2022
@neutrinoceros
Copy link
Member

This seems sound. I'd like to get a deeper look at this tomorrow to make sure I understand what's going on. In particular I want to check that removing the alias does indeed not count as a breaking change, as you suggested.

@neutrinoceros
Copy link
Member

So, I find that there's only one sample dataset that this can be tested with, namely "TNGHalo"
Running this script

import yt

ds = yt.load_sample("TNGHalo")

ad = ds.all_data()

print(ad["gas", "H_number_density"])
print(ad["gas", "H_nuclei_density"])

I get

...
[0.06318717 0.09067829 0.1182974  ... 0.02359772 0.00351616 0.06363773] cm**(-3)
[0.12637434 0.18135657 0.2365948  ... 0.04719544 0.00703232 0.12727546] cm**(-3)

I see that indeed, the "alias" field ("gas", "H_nuclei_density") is incorrectly aliased.
I can also confirm, with this minimal reproducer, that your patch fixes the bug. To my surprise, it also does break my script (I was expecting ("gas", "H_nuclei_density") to simply become unavailable), so this patch is indeed eligible for back porting to yt 4.1.x (yay!).
Could you please add a small test case ? it doesn't have to be much more complicated that my script, changing the print for some assert statements.

@neutrinoceros neutrinoceros added the code frontends Things related to specific frontends label Nov 17, 2022
@jzuhone
Copy link
Contributor Author

jzuhone commented Nov 17, 2022

Grrr stupid python 3.9 test

@neutrinoceros
Copy link
Member

(Soon we'll get rid of these with #4122)

@neutrinoceros neutrinoceros merged commit 6906db5 into yt-project:main Nov 17, 2022
meeseeksmachine pushed a commit to meeseeksmachine/yt that referenced this pull request Nov 17, 2022
neutrinoceros added a commit that referenced this pull request Nov 17, 2022
…1-on-yt-4.1.x

Backport PR #4211 on branch yt-4.1.x ([bugfix] H nuclei are being double-counted in the Arepo frontend)
@jzuhone jzuhone deleted the arepo_H_fix branch April 10, 2023 13:19
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
bug code frontends Things related to specific frontends
Projects
None yet
Development

Successfully merging this pull request may close these issues.

2 participants