-
Notifications
You must be signed in to change notification settings - Fork 8
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
Fix Bug in PTE Solver and Several Code Quality Improvements #383
Conversation
…ll sie and density
…roduces read herrings
I'm still getting negative volume fractions despite the fact that I believe the code should be catching this. I'll add a unit test here shortly that should demonstrate this and then work on whatever fix is needed to prevent this. |
… and temporarily comment out vfrac checks
I think I got all the bugs worked out with running the test on device. The test should be a good platform to add tests that used to break the PTE solver (or still do while we're figuring out a solution). And it can be a bit of a template for using the new |
…CMakeLists and other '*_USE_*' or '*_BUILD_*' variables are in the main project
great work @jhp-lanl! This not only fixes a long-standing bug in the solver, but it uses the latest best practices for portability and adds improvements to the configure/build system. In addition, this provides a tool to explore problematic PTE inputs separate from the context of a host code! Assuming this passes CI and pending review by those that understand the physics better than I, this is ready! 👍 |
I think I'm still having some issues with the CI on gitlab unfortunately. I thought I had fixed these issues on my build when executing on device, but I'll keep working to reproduce |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Asum *= uscale / Tnorm; | ||
rhoBsum /= Tnorm; | ||
PORTABLE_REQUIRE(rhoBsum > 0., "Non-positive energy density"); | ||
PORTABLE_REQUIRE(rhoBsum > 0., "Non-positive Ideal PTE energy density guess"); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Is this only for ideal gas? Because non-ideal gas EOS's can have non-positive energies and it's totally expected/fine.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The next line is
Tideal = uscale / rhoBsum / Tnorm;
This would imply a negative temperature guess, so no, this code cannot admit a negative energy density. I think this is probably part of why it's not useful outside of gas mixtures.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ah ok. I just want to make sure this is only used for ideal gas, as tabulated EOSs definitely support negative energy.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
yep! I think this is a good reason not to use this guess for tabular EOS. You could also imagine a reacting system where the reaction products is represented by a shifted ideal gas with a negative energy. I think this guess isn't really used, but I figured it might be good to put guardrails on it.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This function is called in the Init
for PTESolverRhoU
which I believe is supposed to handle non-ideal EOS.
Notice that TryIdealPTE
is commented out from PTESolverRhoT
with the following comment:
// Leave this in for now, but comment out because I'm not sure it's a good idea
// TryIdealPTE(this);
// Set the current guess for the equilibrium temperature. Note that this is already
// scaled.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I unresolved this conversation because I was concerned that PTESolverRhoU
now no longer works for non-ideal EOS wherein u can be less than 0. Any thoughts here?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I suggest that we floor the energy/temperature inside the solver, as it's only used for guesses
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I unresolved this conversation because I was concerned that
PTESolverRhoU
now no longer works for non-ideal EOS wherein u can be less than 0. Any thoughts here?
If PTESolverRhoU
relies on a negative temperature guess, I'd argue that it never worked for non-ideal EOS that can have negative energies and that the guess infrastructure should be re-worked for that solver.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ah I misread the above. I forgot that TryIdealPTE()
is separate from GetIdealPTE()
. So these error checks being located in GetIdealPTE()
is a mistake since TryIdealPTE()
will check the results of GetIdealPTE()
and see if the results are sane.
Even though these changes are highlighted just because I changed the message and not because they were added here, I'll move them to the TryIdealPTE()
function.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Actually, I ended up getting rid of all the GetIdealPTE
asserts and instead I added a check in TryIdealPTE()
to make sure the temperature guess is actually positive.
@Yurlungur: Jeff (@jhp-lanl) and I discussed the timeline for this PR. Given the timeline for a release into EAP, we are planning on merging this no later than EOD Monday and beginning the release process immediately following. This gives some, but not much time for review but we really can't wait any longer than that. If fixes are needed beyond what is caught on Monday, they will need to go in after the release. |
Yes understood. I believe @pdmullen at least plans to review Monday. We shouldn't wait for @jdolence . |
@Yurlungur I'd really like @jdolence's eyes on it as some point for a sanity check if at all possible, even if it's after the merge. Sometimes the original author has knowledge that is hard to obtain by reading alone. |
That said, I still need to fix the CI. I'll work on that this afternoon. |
I agree. My concern is I don't know when I will be able to get him to look. We were visiting him in Ann Arbor this week and I asked him to look, though. So there's hope. |
I think his input will be valuable whenever we receive it. |
CI Fixed! Array was out of bounds, and I just bit the bullet and enabled proxy variables for all my terminal sessions rather than limiting them to interactive ones. |
Asum *= uscale / Tnorm; | ||
rhoBsum /= Tnorm; | ||
PORTABLE_REQUIRE(rhoBsum > 0., "Non-positive energy density"); | ||
PORTABLE_REQUIRE(rhoBsum > 0., "Non-positive Ideal PTE energy density guess"); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This function is called in the Init
for PTESolverRhoU
which I believe is supposed to handle non-ideal EOS.
Notice that TryIdealPTE
is commented out from PTESolverRhoT
with the following comment:
// Leave this in for now, but comment out because I'm not sure it's a good idea
// TryIdealPTE(this);
// Set the current guess for the equilibrium temperature. Note that this is already
// scaled.
…e check to make sure the guess is valid
I think this is ready to go now |
👍 Merge when ready |
Merging... @dholladay00 I think we're ready for a release |
PR Summary
This PR is dependent on #382
This PR makes a few main changes:
The previousScaleDx
algorithm was written assuming that only one guard (small volume fraction, small temperature, minimum density, etc.) would activate at once. I believe I'm running into a test case where multiple guards are activating so they end up overwriting each other.Actually, since the newly-shrunk
scale
is used in each successiveif
block, they have the effect of only decreasing scale. I'll add a comment to this effect.else
block.ScaleDx
function. First, there were missing parentheses. Second, the logic was incorrect when the volume fraction of a material has already exceeded the imposed maximum. This maximum occurs at the minimum pressure on the EOS, effectively guarding against the EOS entering unphysical (thermodynamically unstable) states. However, it's possible for this maximum to change as the temperature guess changes such that an old volume fraction is above this maximum. In that case, the logic of the function was increasing the step to get back into the physical regime, but possibly at the cost of violating other checks.TestUpdate()
function has been changed so that it's a bit more clear when it's also updating the state versus when it's supposed to cache the previous state.PR Checklist
make format
command after configuring withcmake
.If preparing for a new release, in addition please check the following: