Skip to content

[newchem-cpp] Factor out gas property calculation#529

Merged
brittonsmith merged 34 commits into
grackle-project:newchem-cppfrom
mabruzzo:ncc/gas-props2
May 12, 2026
Merged

[newchem-cpp] Factor out gas property calculation#529
brittonsmith merged 34 commits into
grackle-project:newchem-cppfrom
mabruzzo:ncc/gas-props2

Conversation

@mabruzzo
Copy link
Copy Markdown
Collaborator

@mabruzzo mabruzzo commented Apr 16, 2026

To be reviewed after #528 has been merged


This PR factors out a bunch of calculations pertaining to basic gas properties. I had started working on this 3 months ago before getting distracted. Most of these calculations are now performed by the basic_gas_props function

mabruzzo added 26 commits April 19, 2026 09:40
This switches to using GRIMPL_NAMESPACE_DECL and puts all the helper
function into the namespace
The result is definitely not bitwise identical, but that's not going to
make tests fail (since the logic is not tested)
@mabruzzo
Copy link
Copy Markdown
Collaborator Author

The force-push was the result of a rebase atop the most recent version of PR #528; I choose to rebase (and not merge) since nobody else has pulled in this branch yet

@mabruzzo mabruzzo marked this pull request as ready for review April 30, 2026 13:28
Comment thread src/clib/gas_props.hpp
Comment on lines +297 to +306
#else
// in this branch, we roughly approximate the adjustment
// -> it's equivalent to a single iteration of fixed point iteration
// -> the resulting tgas is technically not self-consistent (i.e.
// applying the adjustment again would yield a different answer)
double adjusted_gamma = chemistry_T_detail::variable_gamma(
tgas[i], nH2, nother, my_chemistry->Gamma);
tgas[i] =
tgas[i] * (adjusted_gamma - 1.) / (my_chemistry->Gamma - 1.);
#endif
Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If this is basically just one iteration of self_consistent_Tgas, can we refactor this to just use that with n_iter=1? Or, I wonder if it's worth just making a runtime parameter out of n_iter with default of 1 and also making another runtime parameter of the tolerance, which is currently 1e-3 (which I suspect is completely arbitrary).

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I mean there is a difference here: using self_consistent_Tgas we apply a temperature floor based on the minimum temperature value in the shared 1d temperature grid. This leads to a few minor considerations:

  1. do we actually want to use the minimum temperature of the shared temperature grid as a floor? If so then we should absolutely document that (if somebody adjust the minimum temperature it could have significant repercussions)
  2. we're already applying that floor later when we compute ln(T). We should stop doing this twice
  3. should we be applying a similar floor for primordial_chemistry == 0?

These questions started to pop up in the subsequent ln(T) PR.


More broadly, I'm ok with us computing a more self-consistent Temperature (in fact it's probably a good thing), but I think fixed point iteration is probably a mistake (I don't think there's a guarantee that it converges at all let alone how fast it converges). If we are going to do self-consistent temperature calculations, I think we should be doing something more robust. I've thought about this in the past, and I actually think Newton's method with analytic derivatives is the way to go -- I describe this in #549 (it's surprisingly tractable and I think it would be significantly faster than fixed point iteration).

Copy link
Copy Markdown
Collaborator Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Let me know how you want to proceed

Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If you have a more robust solution planned I think we should leave this as is and move on. I don't like having the ifdef is the main thing. It's not at all clear to me how much of a difference this makes. I'm more interested in moving forward with everything else.

@brittonsmith brittonsmith merged commit 5e0e568 into grackle-project:newchem-cpp May 12, 2026
5 checks passed
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

2 participants