[newchem-cpp] transcribe: solve_rate_cool_g to C++

[mabruzzo]

This PR was originally proposed as brittonsmith#24

---

This needs to be reviewed after PRs #275, #276, #277, and #315

---

High-level Summary

This PR transcribes solve_rate_cool_g from Fortran to C++ (the helper functions within solve_rate_cool_g.F were not transcribed). I think the following table is somewhat instructive about the significance of this PR

	# of Lines	# of args	# of local variables
original subroutine:	~1700	469	265
transcribed function:	~400	7	38

Overview

At its core, this PR is a straight-forward one-to-one from transcription from Fortran to C++. This is plainly visible in the first couple of commits. All of the subsequent commits correspond to efforts to essentially "polish" the code

All of this polishing was primarily motivated by a desire to minimize unnecessary future work and this routine's status as the lowest common ancestor of the call-graphs of every non-trivial Fortran subroutine. In more detail:

• the arguments of every subroutines further down the call stack are composed of a subset of solve_rate_cool_g's arguments and hundreds of local variables
• a key feature of the transcription tools is that during transcription from Fortran to C++, it can be directed to replace a set of arguments that are read from a single struct with a single argument referring to the struct and modify all references to the original argument appropriately (this can be done for multiple argument-sets)
• Consequently, I was strongly incentivized to try to aggregate the local variables into structs based on how they are used in other subroutines.¹

Polishing the code also made it easier to understand the general control flow (and how the routines are all inter-related), but that's honestly a secondary motivation

Important

To be clear:

• I took great effort to perfectly preserve all of the logic. I not only confirmed that each commit was atomic and passed every test, but also made sure not to refactor or reorder any of the logic. All mathematical operations, array modifications, etc. should be perfectly preserved. There are a small handful of locations where control-flow constructs (namely if/else and goto) have superficial changes (the control-flow is preserved) that are described down below.²
• When looking at the original routine, you should easily be able to locate the corresponding logic in the new function
• All comments have been preserved (or improved!)

More About Polishing:

• I aggregated local variables into a series of structs. I spent a lot of time figuring out how to do this in a way that will make transcription of future subroutines a lot easier. This took a significant amount of time
• I put all calls to untranslated fortran subroutines behind wrapper function.
  • This helped make the purpose of certain local variables a lot more clear and helped me make decisions on how to best aggregate local variables
  • the header for this file also explains how it will also help with the transcription of fortran routines that are called in multiple locations
• I relocated chunks of logic from within the subcycle loop into helper functions.
  • This made the subcycle loop a lot easier to read (and helped me understand how to better aggregate the local variables).
  • they're named: enforce_max_heatcool_subcycle_dt_, select_chem_scheme_update_masks_, calc_Heq_div_dHeqdt_, & set_subcycle_dt_from_chemistry_scheme_
  • these functions all have descriptive docstrings and there are explanatory comments just before each function-call (the logic has not changed -- this was just a relocation)

Important new types:

• IndexRange:
  • this is the subject of PR [newchem-cpp] tweak index_helper.h for use with transcribed routines (and general C++ compatibility) #315
  • I think this will play a central role in running Grackle on GPUs (happy to expand on that)
• InternalGrUnits:
  • this aggregates all of the unit variables commonly passed around the Fortran routines.
  • Perhaps more importantly, I have copied all of the duplicated logic related to constructing these values directly into a function that will construct this type directly from an instance of code_units. I was EXTREMELY careful to make sure I preserved the original logic and make sure we used the right version of all of the constants (The value of pi differs between the C and Fortran cases and the default precision various constants differ)
  • we can worry about better reconciling this with code_units in the future. For no, this simplifies a lot of transcription
• SpeciesRateSolverScratchBuf:
  • the main thing to note here is that this struct aggregates local-variable arrays as well as structs of arrays. This is intended for organizational purposes inside of solve_rate_cool_g (i.e. it provides a clear distinction between local variables required for cooling and variables needed for species chemistry)

Places where logic changed/relocated:

• convert the goto statement into a break statement
• the handling of ierr has been tweaked slightly. (the actual logic is all the same, but now we directly return the value)
• some of the unit-calculation logic has been moved to the functions associated with the new InternalGrUnitstype. But the logic has been perfectly preserved.
• the outer for-loop (over j,k index-pairs) was converted to use the index_helper machinery used in our other existing C routines (the logic is exactly equivalent and this was necessary to use InternalGrUnits)
• the logic where we use the radiation field to modify the iteration mask has been relocated to the coupled_rt_modify_itmask_ local function. There was a minor tweak: I moved the statement checked by both inner if-statements to the outer if-statement

Assorted Closing Thoughts

The scope of this PR definitely ballooned, but I think it's worth it. All of the commits were logically ordered and have a somewhat descriptive commit message (the ordering got a little more chaotic near the end as I kept finding things that I forgot to address).

If it would make things easier, I could divide this PR into smaller parts. But you would need to be okay with merging very messy looking code into the branch.

Note

PRs for other transcribed routines should all be a lot simpler to review. As noted earlier, there were unique and strong incentives to do a lot of "polishing" in this PR.

---

Brief Thoughts for the future:

• If we are willing to embrace a little more C++ I think we could further improve organization
• Now that I finally have a complete understanding of how all of the logic works (I've never understood been able to keep track of all of it at once before), I have a good idea of how we should pursue adding GPU support. (I think there's a "right way" to do it in a performance portable manner. At some point soon, I plan to make an Issue explaining how to do it)

Footnotes

1. The groupings of local variables don't have to be perfect, we can always adjust them later. But, the process of aggregating them will never be easier than it is right now ↩

2. Additionally, there were a couple spots that had large "if-conditions" where I pre-computed parts of the conditions and temporarily stored them within local variables to simply improve legibility. The logic did not change ↩

[brittonsmith]

It can wait, but I think the answer is yes. We have discussed this before and I think we both agree that carrying a metal-related parameter and also checking for the existence of the metal density field is not something we should be doing. If the metal field is needed and not provided, we should error out.

[brittonsmith]

I believe the reasoning for this predates me. I genuinely wonder why we would want to do this. If I understand correctly, we could temporarily find ourselves in a state with short subcycle timesteps, get through it, and then be forced to very slowly increase the timestep again. This seems like it just makes a bad situation worse.

I'll leave this comment only to keep a record since this isn't really the place to interrogate whether this should be here.

[brittonsmith]

Just thinking out loud now, but it may be useful at a later date to understand the frequency with which various timestep limiters are coming into play and evaluating how useful they truly are.

[brittonsmith]

Let's not change this now since it may require a gold standard update, but I am now fairly certain it should be like this metallicity is already in units of Zsun.

Suggested change

	const double min_metallicity = 1.e-9 / my_chemistry->SolarMetalFractionByMass;
	const double min_metallicity = 1.e-9;

[brittonsmith]

In thinking about it, to me having multiple itmasks that need to coordinate is making our lives more difficult. What about a single itmask that can have more than two values. For example:
0: do nothing
1: do gauss-seidel
2: do newton-raphson v1
3: do newton-raphson v2

[brittonsmith]

Maybe this is brittle or dangerous, but I wonder if all this should be inside an itmask[i] check. If we're here, there should always be at least one true itmask value still and we can save the time of looping over the whole array. This is for the future, but putting this here to make a record of it.

[brittonsmith]

This was a joy to read. I've left some notes for possible future changes, but I would rather pull this in first and do those later. I've written down all these notes locally, too. This looks great, let's take the leap.