photoelectric heating - Githubissues

chongchonghe commented 1 month ago

Description

In this PR, I implement photoelectric heating within the current gas-radiation-dust interaction framework. PE can be enabled by setting enable_photoelectric_heating = true in ISM_Traits. The prerequisite is nGroups > 1 and enable_dust_gas_thermal_coupling_model = true.

Added two test problems, RadiationMarshakDustPE-coupled and RadiationMarshakDustPE-decoupled, to test photoelectric heating in cases where gas and dust are well coupled and decoupled. These two tests share the same problem code test_radiation_marshak_dust_and_PE.cpp but use different runtime parameters.

Here is the code structure:

- AddSourceTermsMultiGroup
  - outer loop:
    - compute `v_times_F_term` 
    - solve gas-dust-radiation energy exchange:
      - if constexpr (!enable_dust_gas_thermal_coupling_model_)
        - `SolveGasRadiationEnergyExchange`:
          - The matter-radiation coupling scheme from the multigroup paper
      - else if constexpr (!enable_photoelectric_heating_)
        - `SolveGasDustRadiationEnergyExchange()`:
          - The gas-dust-radiation coupling scheme proposed in #733, which has two cases: the well-coupled case and decoupled case
      - else
        - `SolveGasDustRadiationEnergyExchangeWithPE()`:
          - Based on the gas-dust-radiation coupling scheme, but with modified Jacobian in the well-coupled case to include PE heating. In the decoupled case, PE heating is added via a forward-Euler update using the updated radiation energy density after the gas-dust-radiation step.
    - compute updated `v_times_F_term` 
    - if `v_times_F_term` has converged, break
  - `UpdateFlux()`

Related issues

A modification to #733 to include PE heating.

Checklist

Before this pull request can be reviewed, all of these tasks should be completed. Denote completed tasks with an x inside the square brackets [ ] in the Markdown source below:

[x] I have added a description (see above).
[x] I have added a link to any related issues see (see above).
[x] I have read the Contributing Guide.
[x] I have added tests for any new physics that this PR adds to the code.
[x] I have tested this PR on my local computer and all tests pass.
[x] I have manually triggered the GPU tests with the magic comment /azp run.
[x] I have requested a reviewer for this PR.

psharda commented 1 month ago

PE heating will be included in chemistry. Do you still want to have a separate implementation here in radiation?

markkrumholz commented 1 month ago

PE heating will be included in chemistry. Do you still want to have a separate implementation here in radiation?

Yes -- @chongchonghe and I discussed this, and while we will include this with chemistry, there are also uses cases for a standalone implementation that does not require chemistry. For example, the method implemented in Bate & Keto (2015) does photoelectric heating without chemistry, and ours is loosely based on that. The main use case here is that if we want to simulate the atomic ISM without resolving the molecular phase, we can do the partition of the gas into hot and cold atomic phases with this implementation.

chongchonghe commented 1 month ago

PE heating will be included in chemistry. Do you still want to have a separate implementation here in radiation?

and speed is another big factor. I bet this implementation of multigroup RHD + dust + PE is at least 3 times faster than using the VODE solver in microphysics.

chongchonghe commented 1 month ago

/azp run

azure-pipelines[bot] commented 1 month ago

Azure Pipelines successfully started running 2 pipeline(s).

BenWibking commented 1 month ago

This has the same code maintainability issues that I thought we were trying to avoid by using Microphysics for the radiation update. I would prefer to discuss this on the telecon tomorrow before abandoning our previous plan to use Microphysics.

However, regardless of whether we use Microphysics or not, in principle, I think RadiationSystem should be agnostic to the type of non-LTE effects one wants to model. That is, the code should be structured such that implementation details of those effects should be handled in a separate class.

chongchonghe commented 1 month ago

This has the same code maintainability issues that I thought we were trying to avoid by using Microphysics for the radiation update. I would prefer to discuss this on the telecon tomorrow before abandoning our previous plan to use Microphysics.

However, regardless of whether we use Microphysics or not, in principle, I think RadiationSystem should be agnostic to the type of non-LTE effects one wants to model. That is, the code should be structured such that implementation details of those effects should be handled in a separate class.

This is simple to do. Looking at the code structure I wrote above, it should be easy to put SolveGasDustRadiationEnergyExchange and SolveGasDustRadiationEnergyExchangeWithPE into a separate file. I'll create a issue to remind me to do it later. For this PR, I'll keep it here because this way it's easier to review.

BenWibking commented 1 month ago

This has the same code maintainability issues that I thought we were trying to avoid by using Microphysics for the radiation update. I would prefer to discuss this on the telecon tomorrow before abandoning our previous plan to use Microphysics. However, regardless of whether we use Microphysics or not, in principle, I think RadiationSystem should be agnostic to the type of non-LTE effects one wants to model. That is, the code should be structured such that implementation details of those effects should be handled in a separate class.

This is simple to do. Looking at the code structure I wrote above, it should be easy to put SolveGasDustRadiationEnergyExchange and SolveGasDustRadiationEnergyExchangeWithPE into a separate file. I'll create a issue to remind me to do it later. For this PR, I'll keep it here because this way it's easier to review.

I think it would actually be easier to review if those were in a separate file.

chongchonghe commented 1 month ago

@BenWibking ROCm failed due to 'Disk quota exceeded'. Can you fix that?

markkrumholz commented 1 month ago

Will wait to review this until the opacity simplification PR I just merged is synced with this, then come back and look at it.

chongchonghe commented 1 month ago

OK. I have merged in the dev branch and also moved everything about dust and PE heating into a separate file source_terms_imperfectly_coupled_dust.hpp. Now, source_terms_multi_group.hpp is a lot simpler.

Ready for review and merge.

azure-pipelines[bot] commented 1 month ago

Azure Pipelines successfully started running 2 pipeline(s).

markkrumholz commented 1 month ago

I'm looking through this, and it looks like the Jacobian solver with decoupled dust + radiation and photoelectric heating is still sitting in radiation_system.hpp. I thought the plan was to pull this out into a separate file, so not only would the custom source terms sit in different files, the custom Jacabians would too, and we wouldn't wind up with a bunch of bespoke, model-specific Jacobian solvers sitting in that file. Am I misremembering?

chongchonghe commented 1 month ago

I'm looking through this, and it looks like the Jacobian solver with decoupled dust + radiation and photoelectric heating is still sitting in radiation_system.hpp. I thought the plan was to pull this out into a separate file, so not only would the custom source terms sit in different files, the custom Jacabians would too, and we wouldn't wind up with a bunch of bespoke, model-specific Jacobian solvers sitting in that file. Am I misremembering?

Yes, you are right. I'll make that change.

markkrumholz commented 1 month ago

I'm looking through this, and it looks like the Jacobian solver with decoupled dust + radiation and photoelectric heating is still sitting in radiation_system.hpp. I thought the plan was to pull this out into a separate file, so not only would the custom source terms sit in different files, the custom Jacabians would too, and we wouldn't wind up with a bunch of bespoke, model-specific Jacobian solvers sitting in that file. Am I misremembering?

Yes, you are right. I'll make that change.

Great. Just ping me when I should look again.

chongchonghe commented 1 month ago

@markkrumholz OK, all comments addressed.

chongchonghe commented 1 month ago

/azp run

azure-pipelines[bot] commented 1 month ago

Azure Pipelines successfully started running 2 pipeline(s).

markkrumholz commented 1 month ago

OK, another stylistic question. The design here seems a bit different than our normal paradigm, in that you have the various source and exchange models protected by ifdef brackets, and you forward-declare the various exchange and solver functions is radiation_system.hpp so that even if you turn off a particular exchange model by not #define'ing the appropriate flag to turn it on, you don't get compilation errors thanks to the forward declaration. This is contrary to our style thus far of mostly using if constexpr's and templates to manage conditional compilation. The paradigm that I would think would work would be remove the forward declarations from radiation_system.hpp and the protection ifdef's from the various model files, and just #include the various model files in radiation_system.hpp, but then select which model to use at compile time with a trait in the problem file. Is there a reason not to use this paradigm here?

chongchonghe commented 1 month ago

OK, another stylistic question. The design here seems a bit different than our normal paradigm, in that you have the various source and exchange models protected by ifdef brackets, and you forward-declare the various exchange and solver functions is radiation_system.hpp so that even if you turn off a particular exchange model by not #define'ing the appropriate flag to turn it on, you don't get compilation errors thanks to the forward declaration. This is contrary to our style thus far of mostly using if constexpr's and templates to manage conditional compilation. The paradigm that I would think would work would be remove the forward declarations from radiation_system.hpp and the protection ifdef's from the various model files, and just #include the various model files in radiation_system.hpp, but then select which model to use at compile time with a trait in the problem file. Is there a reason not to use this paradigm here?

I also tried another paradigm and I just pushed it. In this paradigm, there are two main files in src/radiation: radiation_system.hpp and radiation_dust_system.hpp. The latter file contains all the definition of the functions used solely in dust or dust+PE models, and it includes radiation_system.hpp. Now, in most tests where there is no dust model, we include radiation_system.hpp and those functions for dust is never read by the compiler. In tests with dust, we include radiation_dust_system.hpp in the header. Does this resolve your concern?

Both of the two hpp files have to be protected by ifdef in this new paradigm.

chongchonghe commented 1 month ago

/azp run

azure-pipelines[bot] commented 1 month ago

Azure Pipelines successfully started running 2 pipeline(s).

markkrumholz commented 1 month ago

Thanks, that structure is much cleaner. I left a few minor stylistic comments (mostly about your comments!), but once those are resolved I'm happy to approve.

chongchonghe commented 1 month ago

OK. All done.