-
Notifications
You must be signed in to change notification settings - Fork 32
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
fixed dual-energy formalism synchronization. #356
Conversation
I expect @evaneschneider is going to want my opinion about this... here it is. As mentioned, the dual energy condition is important in preventing numerical errors when the kinetic energy is much greater than the thermal energy, but it can also suppress shock heating if it is not calibrated correctly. Currently, Cholla applies the eta_1 condition when computing the pressure from the conserved variables and the eta_2 condition when synchronizing the advected thermal energy with Note that the current values of eta_1 and eta_2 (defined here) are:
The value of eta_2 has been calibrated to produce gas temperatures in cosmological simulations that match the expectation from combining the IGM temperature and the temperature of the virialized gas that collapses into dark matter halos. If you are interested, a discussion is presented in sections 2.3 and 3.4 of Villasenor et al. 2021. Since the If you decide to proceed with this change, I recommend leaving the current condition (only the eta_2 condition) when running cosmological simulations. @mabruzzo, is there a particular reason for you to suggest this change? Have you observed unexpected behavior somewhere? |
@bvillasen Thanks for your comments! A little more background: I asked @mabruzzo to make this change because we've been running some galaxy sims with supernova lately and seeing incorrect results inside the supernova bubbles -- in particular, with dual energy turned on and the current version of Select_Internal_Energy, the gas inside the bubbles ends up being cold, rather than hot. I found that when I turned dual energy off, the problem was fixed, and realized that the issue is the way that the total energy synchronization after the selection is currently applied. I believe that in the original formulation of DE in Cholla, both eta criteria had to be met in order for the value of the total energy to be updated with internal energy value at the end of the function. As currently written, ONLY the eta_2 check is applied, which means that in any circumstance other than when the eta_2 check is met, the value of the total energy is updated to contain the advected internal energy instead of the conservatively calculated value. This I think is definitely wrong, since it applies all over the place if there are no nearby shocks, and has nothing to do with whether you actually need the additional precision (as Matthew noted above). It also gave rise to the cold gas in bubbles we saw above. To summarize, I actually think the "or" statement in this corrected version should be an "and" -- you only update the total energy with the advected internal energy if both criteria are met, and you use the conservatively calculated value of the internal energy in all other cases. That is the version that I have in my branch (where things seem to be behaving correctly now), but I haven't done extensive tests, or looked at the cosmological case, which I imagine could be different. That said, as Matthew points out, the discussion of dual energy in Enzo doesn't actually talk about synchronization in the context of eta_2. So it's possible that a better solution is to not synchronize the total energy with the internal energy at the end of the step -- basically let the advected internal energy do its thing at all times unless the eta_2 criterion is met, and just use eta_1 to figure out when to use that quantity elsewhere in the code, as @bvillasen currently has it written. My understanding from @mabruzzo is that is more similar to what Enzo actually does, and it would avoid the problem we are currently having with the total energy being updated to contain the wrong value of the internal energy. Happy to hear more thoughts on this if anyone has them! |
Hi Bruno, thanks for reaching out! We definitely value your opinions on this! I've actually got some experience with the dual-energy formalism (Several years back, I needed to implement it in the VL+CT integrator I had implemented in Enzo-E1 and I've spoken at length about the topic with Greg Bryan). When I came across your paper a little while back I found it very insightful! (Calibrating the value of Unexpected behaviorBefore saying anything else, let me address:
@EvanSchneider and I have been trying to debug issues with supernovae bubbles in galaxy-scale simulations. When running these simulations we were finding that the cavities formed by supernovae were filled with lots of cool gas. A little while back Evan realized that the dual-energy formalism was making the problem worse. More recently Evan realized that the application of Highlighting subtle differences
So your recollection is exactly right! As we've both established, the description in the Enzo-paper talks about tracking the the separately-advected internal energy, However, there is a subtle, important difference between the description in that paper and what Cholla actually does (and for that matter what Enzo does). The difference comes down to "synchronization". In the Enzo-paper, the only discussion of "synchronization" is relates to updating the value of
Now, in practice, there is a strong temptation to somehow update
Discussion of your concerns
This is a very good point! You make a compelling case that this may probably change the temperature structure of the cosmological simulation. I certainly don't want to do that. With that said, this highlights a potential issue with how well the current choice of
Thus, if the proposed change makes a big difference in the temperature of halos, this suggests that our calibrated
Where to go from hereSo I think we have 3 options:
B. We can go with the approach historically used in Enzo's ppm solver (variation number 3 from above). This doesn't let C. We could go with the approach described in the Enzo-paper (i.e. we never overwrite the value of I'm a little agnostic here. On the one hand, I like choices 2 or 3 a lot in principle (and they are somewhat battle-tested). On the other hand, the fact that Enzo-Classic doesn't really use choices 2 or 3 and hasn't for the last 12 years mutes my feelings on this EDIT: Changed enumeration of options from numbers to letters. Footnotes
|
I'll update my comment based on @mabruzzo's more extensive and better articulated explanation: |
IMO, we should do whatever is the most physically motivated. I agree that “and” seems preferable to “or”, but we need to check the behavior of the cosmological simulation test before we know whether this should be used generally or if the behavior should be optional.On Dec 11, 2023, at 6:47 AM, Evan Schneider ***@***.***> wrote:
I'll update my comment based on @mabruzzo more extensive and better articulated explanation:
The reason I prefer "and" to "or" in this check is that my thinking on dual energy in general is that it should only be necessary in a scenario where we lack precision in the internal energy calculation due to truncation error. So, if the eta_1 criterion isn't met, we should not need it (and definitely should not change the value of the total energy based on the non-conservative update). I guess that is most similar in philosophy to @mabruzzo's option 3, since it essentially imposes quite a strict set of requirements before you are allowed to overwrite the value in E. That said, if @bvillasen thinks that abandoning the synchronization at the end of the function (and therefore, only applying the eta_2 check to adjust the value of e) would make this fix more compatible with the cosmological sims, I'm happy for us to go that route.
—Reply to this email directly, view it on GitHub, or unsubscribe.You are receiving this because you are subscribed to this thread.Message ID: ***@***.***>
|
About the current implementation: @evaneschneider I think you and I are conceptually on the same page, but I still think "or" is correct (it's possible I am misunderstanding the code). In the currently proposed change: if ((U_total / E > eta_1) or (U_total / Emax > eta_2)) {
U = U_total;
} else {
U = U_advected;
} We are saying here, that we overwrite
If we made that an "and" statement, then its possible that we won't overwrite the |
Just a note - I realized that the enumerated options don't match up with the enumerated variations of the dual-energy formalism. Therefore I just changed the enumerated options to use letters instead of numbers (to try to avoid confusion) |
Ah yes, I see. I agree! |
Could we have a meeting to discuss this issue?On Dec 11, 2023, at 7:33 AM, Evan Schneider ***@***.***> wrote:
About the current implementation: @evaneschneider I think you and I are conceptually on the same page, but I still think "or" is correct (it's possible I am misunderstanding the code).
In the currently proposed change:
if ((U_total / E > eta_1) or (U_total / Emax > eta_2)) {
U = U_total;
} else {
U = U_advected;
}
We are saying here, that we overwrite internal_energy with the value computed from total_energy if EITHER:
E−v2//E>η1 -- the total energy has enough precision to make the shock-heating unnecessary
ρ(E−v2/2)/max(ρE)>η2 -- there's possibly nearby shock-heating that must be included
If we made that an "and" statement, then its possible that we won't overwrite the internal_energy in cases where there is strong shock heating (in that scenario the precision of the internal-energy computed from total_energy would be low)
Ah yes, I see. I agree!
—Reply to this email directly, view it on GitHub, or unsubscribe.You are receiving this because you commented.Message ID: ***@***.***>
|
Sorry for bombarding you with messages -- Just to clarify:
Again, I'm fine with anything. Maybe scheduling a meeting would be good. |
HI @mabruzzo, can we get this one finalized today so I can merge it in? |
Background
A brief summary of the problem (for future generations), the dual-energy formalism (see the end of section 4.1.1 of the Enzo-paper for a summary) is parameterized by quantities that$\eta_1$ and $\eta_2$
about$\eta_1$ : The whole point of tracking the internal-energy field is to have extra precision for the thermal energy when the kinetic energy is a very large fraction of the total energy). Essentially, $\eta_1$ governs when this this extra precision becomes necessary. We don't don't bother using the internal-energy field unless $((Etot-KE) / Etot) < \eta_1$ , where $\eta_1$ is usually 0.001.
about$\eta_2$ : you can think of this as a synchronization parameter that tells us when we should updating the internal energy-field field to use the value computed from the total energy field. This is most relevant for including the effects of shock-heating.
It's worth mentioning that in most descriptions of the dual-energy formalism, only$\eta_2$ is discussed in the context of synchronization. However, $\eta_1$ becomes relevant in Cholla during synchronization of the internal energy because immediately after we synchronize the internal energy field, we use recompute the total energy from the kinetic energy and the total energy (If we didn't recompute the total energy - we would not need to consider $\eta_1$ at this point).
What I changed
The logic I needed to change was choosing which energy is assigned to the Internal Energy field by the
Select_Internal_Energy_3D
kernel (& in the 2D and 1D versions) since the subsequent kernel that gets called synchronizes total energy with gas energy.if-statement
:if (U_total / Emax > eta_2) {
if ((U_total / E > eta_1) or (U_total / Emax > eta_2)) {