Add non-relativistic radiation transport solvers (explicit and implicit) and cosmic ray transport #492

felker · 2023-04-13T19:14:23Z

Just opening this draft PR to start parsing the massive amount of changes and create an area for discussion.

Would definitely need some Wiki documentation before merging. At a minimum:
- New User Guide page "Newtonian radiation"
- https://github.com/PrincetonUniversity/athena/wiki/The-Input-File
- https://github.com/PrincetonUniversity/athena/wiki/Configuring
Consider excising separating particles, BlockFFT, thermal conduction from this PR? https://github.com/PrincetonUniversity/athena/tree/rad_newtonian/src/thermal_conduction
tst/multi_group/ needs to be reorganized and cleaned up.
- This would be the first time we are considering adding .ipynb notebooks to version control; do we want to do that?
- Remove the many duplicates of athena_read.py

CHANGELOG.md

CONTRIBUTING.md

Makefile.in

README.md

configure.py

yanfeij · 2023-04-14T04:03:52Z

@felker How to add equations in the documentation easily.

munan · 2023-04-14T12:01:28Z

I am working on merging the chemistry_scalar branch, which should be done in the next ~2 weeks before the Athena++ conference. There is also a radiation class there, although it's mostly just a container. See src/radiation on the branch.

c-white · 2023-04-14T12:52:15Z

@felker How to add equations in the documentation easily.

@yanfeij I just learned that as of very recently, github brought back mathjax support! So you can enclose basic latex commands in $...$ or $$...$$ and it should work.

yanfeij · 2023-04-14T14:04:18Z

I am working on merging the chemistry_scalar branch, which should be done in the next ~2 weeks before the Athena++ conference. There is also a radiation class there, although it's mostly just a container. See src/radiation on the branch.

I am afraid it will have exactly the same name as what I am using. Looks like you are just following six fixed rays, which can be a very special case of general radiation transport.

yanfeij · 2023-04-14T14:04:41Z

@felker How to add equations in the documentation easily.

@yanfeij I just learned that as of very recently, github brought back mathjax support! So you can enclose basic latex commands in $...$ or $$...$$ and it should work.

It works! Thank you.

munan · 2023-04-14T14:53:18Z

@yanfeij Yes, the six-ray integrator is a very special case. I am hoping the container Radiation class is general enough to be shared among different radiation modules. Perhaps we should discuss this. See src/radiation/radiation.hpp (or.cpp) in chemistry_scalar branch.

yanfeij · 2023-04-14T15:02:51Z

@yanfeij Yes, the six-ray integrator is a very special case. I am hoping the container Radiation class is general enough to be shared among different radiation modules. Perhaps we should discuss this. See src/radiation/radiation.hpp (or.cpp) in chemistry_scalar branch.

There are many more steps to do for general radiation transport. The implicit radiation transport is also using a completely different integrator. The data array structure and boundary values I need are also different. You can also take a look at the rad_newtonian branch to get an idea.

c-white · 2023-04-14T15:19:20Z

Just as a heads up, after this is merged I will work on merging the gr_rad branch as well. This will probably be more difficult, as there is a lot of overlap in files and functions but not details of what the code is doing.

Maybe we should decide among the following cases soon:

All radiation -- whether written by @munan or @yanfeij or me or anyone else) will eventually live in the same radiation module (i.e., the same class Radiation and the same directory src/radiation/). This might make for a confusingly large module.
Each radiation implementation will be separate, unless they already coexist (maybe some of @yanfeij 's do). There might be some duplicated code, but at least those interested in one use case would not have to sift through code pertaining to vastly different cases.
Something in between.

yanfeij · 2023-04-14T15:38:50Z

Just as a heads up, after this is merged I will work on merging the gr_rad branch as well. This will probably be more difficult, as there is a lot of overlap in files and functions but not details of what the code is doing.

Maybe we should decide among the following cases soon:

All radiation -- whether written by @munan or @yanfeij or me or anyone else) will eventually live in the same radiation module (i.e., the same class Radiation and the same directory src/radiation/). This might make for a confusingly large module.

Each radiation implementation will be separate, unless they already coexist (maybe some of @yanfeij 's do). There might be some duplicated code, but at least those interested in one use case would not have to sift through code pertaining to vastly different cases.

Something in between.

The two modules I have share some common data structures. But I think it will be a lot of work to make all these "radiation" to use the same integrator or data arrays. Also we do not expect to turn on any two of these all at once. I think it is easier to keep each one separate. But do give them some specific names, not just radiation. For example, six_ray_tracing, gr_radiation or something like that.

… into rad_newtonian merge remote changes

yanfeij · 2023-07-11T04:04:51Z

OK. My branch has a long history.

…

On Tue, Jul 11, 2023 at 11:55 AM Kyle Gerard Felker < ***@***.***> wrote: ***@***.**** commented on this pull request. ------------------------------ On tst/regression/scripts/tests/shearingbox/mri3d.py <#492 (comment)> : looks like a dirty merge from master. @tomo-ono <https://github.com/tomo-ono> removed the test in Dec 2020 in 8271040 <8271040> (not sure why though?) We should also remove https://github.com/PrincetonUniversity/athena/blob/rad_newtonian/tst/regression/data/mhd_mri_3d.hst too — Reply to this email directly, view it on GitHub <#492 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ACCWEJYIZV2DFJMDBJONPVTXPTFCBANCNFSM6AAAAAAW5P23YY> . You are receiving this because you were mentioned.Message ID: ***@***.***>

-- Best, Yan-Fei Jiang Associate Research Scientist, Flatiron Institute, 160 Fifth Avenue New York, NY 10010

felker · 2023-07-11T04:07:52Z

@tomidakn might get a chance to take a look at this PR on Wednesday--- if not, we merge immediately on Thursday

felker · 2023-07-12T04:47:51Z

retest this

tomidakn · 2023-07-12T16:55:56Z

I've reviewed about half. Will do the rest tomorrow.

tomidakn

Thank you very much for the huge effort. I quickly reviewed the code, and I just have some comments, mostly cosmetic. I think some files (matrix utils) can be merged, and also some function names can be more self-explanatory.

src/coordinates/cylindrical.cpp

src/coordinates/spherical_polar.cpp

tomidakn · 2023-07-12T16:33:56Z

src/nr_radiation/frequencygrid.cpp

+// frequency is scaled with kT_0/h
+// using fitting formula to return \int_0^nu_min and \int_0^nu_max
+
+Real NRRadiation::FitBlackBody(Real nu_t) {


The function names "FitBlackBody" and "BlackBodySpec" sound confusing to me. I prefer to explicitly indicate that they are integration of the Planck function.

src/nr_radiation/implicit/rad_iteration.cpp

src/nr_radiation/integrators/rad_source.cpp

tomidakn · 2023-07-13T08:13:59Z

src/outputs/outputs.cpp

+  }// End (RADIATION_ENABLED)
+
+
+  if (CR_ENABLED) {


How about "adding cr" data set?

can you clarify? like another user option for the output block like outputN/variable=cr?

Yes. It may convenient to group CR data sets (as well as radiation data sets).

tomidakn · 2023-07-13T09:08:41Z

src/task_list/im_rad_task_list.hpp

+struct IMRadTask {
+  TaskID task_id;      //!> encodes task using bit positions in IMRadTaskNames
+  TaskID dependency;   //!> encodes dependencies to other tasks using IMRadTaskNames
+  TaskStatus (IMRadTaskList::*TaskFunc)(MeshBlock *); //!> ptr to a task
+};
+
+//----------------------------------------------------------------------------------------
+//! \class MultigridTaskList
+//! \brief data and function definitions for IMRadTaskList class
+
+class IMRadTaskList {


Can we derive the TaskList class instead of redefining these codes? Also, the comment must be updated.

I just noticed that MultigridTaskList is not derived either?

True, but it is intentional for a few reasons. MG works on Multigrid objects, not MeshBlocks. This is for future extension to apply different parallelizations for MG and the main part. Maybe I can still derive it from the TaskList, though.

For this radiation, I think it can be done similarly to STS Task List.

I think the IMRadTaskList has a similar argument as the MultigridTaskList then, since the function signature is also different. It does not have a notion of a time integrator stage:

TaskStatus (IMRadTaskList::*TaskFunc)(MeshBlock *); //!> ptr to a task

The SuperTimeStepTaskList does use int stage in a similar fashion to the main task list, so it isnt quite the same as the radiation task list.

Anyway, let's continue this discussion on #537

tomidakn · 2023-07-13T10:54:42Z

src/utils/gauleg.cpp

+
+
+//--------------------------------------------------------------------------------------
+//! \fn  void Gauleg(int n, Real x1, Real x2,  AthenaArray<Real> &x,


We do not severely limit the length of function names. I think it is better to use names that are not shortened in general. For example, GaussLegendre().

src/utils/inverse_matrix.cpp

… into rad_newtonian merge remote changes

tomidakn · 2023-07-18T03:55:35Z

@yanfeij, when you have time, can you add some practical usage of the radiation transport and CR transport modules in the Wiki, including outputs and analyses? For example, Brief description of the test codes is helpful for users.

yanfeij · 2023-07-18T09:12:59Z

Will do after my travel this week.

…

On Tue, Jul 18, 2023 at 11:55 AM Kengo TOMIDA ***@***.***> wrote: @yanfeij <https://github.com/yanfeij>, when you have time, can you add some practical usage of the radiation transport and CR transport modules in the Wiki, including outputs and analyses? For example, Brief description of the test codes is helpful for users. — Reply to this email directly, view it on GitHub <#492 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ACCWEJ3XPZCO5SDAOZAUMB3XQYCMFANCNFSM6AAAAAAW5P23YY> . You are receiving this because you were mentioned.Message ID: ***@***.***>

-- Best, Yan-Fei Jiang Associate Research Scientist, Flatiron Institute, 160 Fifth Avenue New York, NY 10010

tomidakn · 2023-07-19T04:05:04Z

Will do after my travel this week.

Thanks!!

Merge radiation with master

f1ed2b0