Photoheating from radiative transfer

[brittonsmith]

Originally reported by Yves Revaz (Bitbucket: revaz, GitHub: revaz)

---

Hi all,

Currently, when using the heating rates from radiative transfer solutions, only the heating due to HI is included through the variable photogamma:

in cool1d_multi_g.F
! Photoheating from radiative transfer
edot(i) = edot(i) + real(ipiht, DKIND) * photogamma(i,j,k) / coolunit * HI(i,j,k) /dom

However, when computing the radiative heating due to the UB background through the rates provided by the cloudy table, grackle also include the He heating:

in cool1d_multi_g.F
! Photoionization heating
edot(i) = edot(i) + real(ipiht, DKIND)*(
& piHI *HI (i,j,k) ! pi of HI
& + piHeI *HeI (i,j,k)0.25_DKIND ! pi of HeI
& + piHeIIHeII(i,j,k)*0.25_DKIND ! pi of HeII
& )/dom

For the sake of consistency, in the radiative transfert mode, I suggest to let the possibility to provide piHI, piHeI and piHeII, replacing photogamma.

Doing so will still be compatible with the current version, if setting piHeI=0 and piHeII=0.

I think there is very few lines to change and if this is fine, I can certainly to the modifications.

yves

---

• Bitbucket: https://bitbucket.org/grackle/grackle/issue/15

[brittonsmith]

Original comment by Britton Smith (Bitbucket: brittonsmith, GitHub: brittonsmith)

---

Hi Yves,

Like many things in grackle, this has been carried over from how things are done in Enzo. The modification you propose would be very welcome! I definitely encourage you to go ahead with it. Please, feel free to contact the users list or post to this issue if you have any questions. Thanks for taking this on!

Britton

[brittonsmith]

Original comment by Yves Revaz (Bitbucket: revaz, GitHub: revaz)

---

Hi Britton,

Thanks for the very motivating answer.
I went deeper into the code and have other suggestions still related to the rates linked to the radiative transfer.
My motivation is that the ionizing/heating rates provided by the cloudy tables should be replaced by ionizing/heating rates
provided by the radiative transver variables (at the execption that the rates provided by the cloudy tables cans be shielded by different
methods). Implementing this "symetry" will be very very usefull to my opinion.

Thus, I suggest to rename the RT variables to be more consistent with the names addopted for the cloudy tables. This will
strongly improve the code clarity:

variables in the C wrapper:

RT_HI_ionization_rate, (no modification)
RT_HeI_ionization_rate, (no modification)
RT_HeII_ionization_rate, (no modification)
RT_H2_dissociation_rate, (no modification)
RT_heating_rate, ->RT_HI_heating_rate
->RT_HeI_heating_rate
->RT_HeII_heating_rate

variables in the FORTRAN routines:

kphHI, ->RT_k24
kphHeI, ->RT_k26
kphHeII, ->RT_k25
kdissH2I, ->RT_k31
photogamma, ->RT_piHI
->RT_piHeI
->RT_piHeII

By the way, another improvement would be to sum RT_k24, RT_k26, RT_k25 directly
with the shielded values : k24shield, k26shield, k25shield. This will remove a lots of

if (iradtrans .eq. 1)

slightly speeding the code, as it is done for

k31shield(i) = k31 + kdissH2I(i,j,k)

and also improve the code clarity.

I can certainly implemente those suggestions but as it implies quite a lots of changes,
I would like again you opinions before investing time on that.

Cheers,

yves

[brittonsmith]

Original comment by Britton Smith (Bitbucket: brittonsmith, GitHub: brittonsmith)

---

Hi Yves,

When you say that you would like to replace the ionizing/heating rates provided by the Cloudy tables with the radiative transfer variables, do you mean that you would like to completely remove the use of the Cloudy rates? I believe running Grackle with the Cloudy tables is the most popular use right now, so I don't think that should happen. However, if you just mean that you would like to add the option of providing H/He/He+ heating rates from the radiative transfer, then I think this would be a very nice addition.

As for renaming the variables in the C wrapper, I agree with this in principle, but such a change would break backward compatibility, which we would like very much to preserve. I think it would be fine to add RT_H*_heating_rate variables and allow for both the current use case as well as providing them individually.

For renaming the internal FORTRAN variables, this I think would be ok. We usually try not to do too much variable renaming as it can cause issues with merging in pull requests, but I think this would be fine as it would serve to clarify things for developers.

With regard to summing the RT and shielded rates, this should probably stay as is since running with radiative transfer is optional and I think not that widely used to begin with. When radiative transfer is not used, the kdissH2I fields will be NULL pointers. If you have something different in mind, please correct me.

In conclusion, I definitely encourage you to issue a pull request for the first two of your proposed changes as long as we can preserve backward compatibility. You are also welcome to issue multiple PRs if you think it would be simpler. Finally, if you'd like to get more opinions on any proposed changes, I would suggest you write to the mailing list as more people will see it.

Thanks!
Britton

[brittonsmith]

Original comment by Yves Revaz (Bitbucket: revaz, GitHub: revaz)

---

Hi Britton,
Thanks for your mail which precisely answers my questions.

• Of course my aim was not to remove the use of the Cloudy tables, but just provide H/He/He+ heating rates, the same way they are provided by the Cloudy tables.
• Concerning the C wrapper, I agree with the backward compatibility preservation. I will take care of that.
• Ok for the FORTRAN renaming variables.
• Summing the RT and shielded rates is definitively not a necessity as it has no impact on the physics. Following your advise, I won't change that. It will save me a lot of times as this was the most tricky part to do.

Cheers,
yves

[brittonsmith]

Original comment by Britton Smith (Bitbucket: brittonsmith, GitHub: brittonsmith)

---

Fantastic! Looking forward to your pull request(s)!