Fix edge cases for saturated PureFluid states / illustrate saturated water properties

[ischoegl]

Changes proposed in this pull request

• fix TQ/PQ assignment at critical point
• add example illustrating saturated water properties / vapor dome
• fix various bugs for edge case calculations (critical point / triple point)
• disable TP setter for saturated mixtures

If applicable, fill in the issue number this pull request is fixing

Fixes #906, fixes #915, fixes #917, fixes #919, fixes #920

Checklist

• There is a clear use-case for this code change
• The commit message has a short title & references relevant issues
• Build passes (scons build & scons test) and unit tests address code coverage
• The pull request is ready for review

Thoughts

I generated a digital version of a saturated steam table for an exam; i am including it here as it could be considered as an example and potentially expanded upon (otherwise I’ll simply drop it).

[codecov]

Codecov Report

Merging #907 (41c0499) into main (d8e62ad) will increase coverage by 0.04%.
The diff coverage is 94.40%.

@@            Coverage Diff             @@
##             main     #907      +/-   ##
==========================================
+ Coverage   71.22%   71.26%   +0.04%     
==========================================
  Files         377      377              
  Lines       46272    46296      +24     
==========================================
+ Hits        32955    32991      +36     
+ Misses      13317    13305      -12     

Impacted Files	Coverage Δ
include/cantera/tpx/Sub.h	89.47% <ø> (ø)
src/tpx/Sub.cpp	89.16% <94.35%> (+2.55%)	⬆️
src/tpx/Water.cpp	97.50% <100.00%> (ø)
src/thermo/PureFluidPhase.cpp	59.72% <0.00%> (+0.92%)	⬆️

---

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update d8e62ad...864caa0. Read the comment docs.

[bryanwweber]

FYI the failures in the multiple-sundials job are fixed by 529b799 in #905

[speth]

Unfortunately, I think fixing this problem is more complicated than just changing the comparisons in the state setters. With this change, the state of the mixture when setting the state to be a saturated state at the critical temperature does not end up at the critical point, and the resulting state depends on the previous state of the mixture. For example:

>>> w = ct.Water()
>>> w.TP = 300, 101325
>>> w.TQ = w.critical_temperature, 0
>>> w.TP
(647.286, 1089571723.7984445) # not even close
>>> w.TQ = w.critical_temperature, 1
>>> w.TP
(647.286, 7623.02783606933) # wrong, and different
>>> w.TP = 645, 101325
>>> w.TQ = w.critical_temperature, 0
(647.286, 22304442.572488125) # closer, but still wrong

[ischoegl]

That’s actually much worse than what I recall for the tests I ran - I saw an insignificant discrepancy depending on whether I set Q equal to 0 or 1, but nothing of that magnitude. Based on what you are showing, this is likely due to having marched along saturated liquid and vapor lines, meaning that I was reasonably close to the correct result, which led me to believe that the simple fix worked as expected. Update: located the bug - see #915 - and added the fix to this PR.

Since this is just a papercut for an edge case, I’d suggest to close both issue and PR - if there’s interest in the example, I could still modify it as it’s not too difficult to circumnavigate. There’s also the question about the proposed change making physical sense, while being mathematically ill defined (there is no unique value for Q).

[ischoegl]

PS: it may still make sense to change the error message, as it clearly doesn’t apply to the critical point itself; I’ll likely push a suggestion before long.

[bryanwweber]

if there’s interest in the example, I could still modify it as it’s not too difficult to circumnavigate.

I think the example is valuable, so I'd be in favor of fixing it up.

There’s also the question about the proposed change making physical sense, while being mathematically ill defined (there is no unique value for Q).

This was also my concern about the proposed change that I never got around to writing up.

[ischoegl]

@speth / @bryanwweber ... I replaced the 'fix' by just changing the error message, and updated the example.

On a separate note, I think that this is the last pending PR from my side - unless there's interest in resurrecting #888 (a simplification of .gitignore that has no other technical contribution). #904 is something that I likely cannot push further, and the two 'stubs' depend on other work (or clarifications).

PS: I haven't forgotten about Cantera/cantera-website#104 and Cantera/cantera-jupyter#27, which I'll get to one of these days.

Update: Adding a unit test surfaced a previously unknown issue for the PQ setter (the exception was actually not reachable), which is now fixed. Edit: another issue surfaced with #915.

[ischoegl]

Added fixes for #915.

This should open the option to allow for TQ/PQ convenience setters at the critical point after all (with the mathematical caveat that the value of Q would not matter). At this point, I am leaning against it.

[ischoegl]

I added a couple of tests to the suite. There is one last edge case left for discussion: at the critical point, Q can be read, although it cannot be specified in setters and a unique value cannot be defined. Nevermind: this is a larger issue, to be discussed in #916.

Work on this PR is done, except for review.

[speth]

There’s also the question about the proposed change making physical sense, while being mathematically ill defined (there is no unique value for Q).

I'm not convinced that it's mathematically ill defined just because any value of Q corresponds to the same state. I think it's just that there's no longer a natural inverse (that is, a way to determine Q from the intrinsic state) based on the physical definition of Q. That doesn't preclude defining Q in a consistent way at the critical point or even outside the vapor dome.

Of course, we wouldn't want to add this complication unless it were useful, but I think the vapordome.py example makes a case for being able to set the state to the Tcrit and any value of Q. Otherwise, you're forced to jump through hoops in order to traverse the whole saturation curve.

[ischoegl]

@speth ... the convenience argument is what made me file the issue and this PR in the first place, so it’s not difficult to persuade me here 😉 ... I’ll push an update shortly.

Regarding definitions, it’s no longer possible to differentiate liquid from vapor at the critical point, so it actually becomes impossible to define Q. It’s a different discussion (#916) though and won’t affect the implementation here.

[ischoegl]

@speth ... TQ/PQ setters are now enabled for the critical point.

I found quite a few irregularities of various solvers for edge cases while working on this PR, and fixed whatever popped up (it didn’t quite extend to #605 though). I also added unit tests to ensure proper behavior.

[ischoegl]

@speth ... thanks for the review, I believe I addressed everything.

Since triple point properties are currently not broken out (Python's current min_temp is actually T_triple, and there isn't an equivalent min_pressure), I'd suggest leaving a pre-calculation for another time. If you agree, I will open a feature request on Cantera/enhancements ...

PS: I can rebase once #946 is merged, so CI will give expected results.

[ischoegl]

close/reopen to trigger checks

[speth]

I think this looks fine, and your suggestions for a few fine points that can resolved later are reasonable. My only suggestion is to flag a couple of these things in the code with @todo so they're easier to find later, although keeping track of them in the related Github issue is probably also a good idea.

[ischoegl]

I think this looks fine, and your suggestions for a few fine points that can resolved later are reasonable. My only suggestion is to flag a couple of these things in the code with @todo so they're easier to find later, although keeping track of them in the related Github issue is probably also a good idea.

@speth ... thanks! Your suggestions are straightforward (as long as there isn't a preference for upper/lowercase - I used @TODO).