chmarti1 / PYroMat

PYroMat thermodynamic properties in Python
http://pyromat.org
Other
71 stars 13 forks source link

Numerics iteration boundary issue in ._d() #44

Closed jranalli closed 1 year ago

jranalli commented 2 years ago

Similar to issue #43, initial guesses for certain cases end up being unacceptable for generating the initial iterations. A few test cases for _d() are:

sub = pm.get('mp.CH4')
print(sub.d(p=0.11707696, T=624.375)) # T supercritical, p low
print(sub.d(p=9999.0, T=624.375)) # T&p supercritical

What I found here is that at high temperatures and various pressures, you get insane pressures at the high density boundary and it screws up the initiation of the iteration. You need better guesses of the possible limiting densities.

Temperature above the dome, and pressure below can use the critical density as a lower limit. The lower limit for both supercritical is harder to come by and I don't have a perfect suggestion there.

chmarti1 commented 2 years ago

It's starting to look like a specialized bracketing function could be a good idea instead of requiring each code segment that calls _hybrid1() to define its own bracket. It's clear that there are enough nuances to the bracketing problem that it is not efficient to solve them one-by-one. This is now my first priority for v2.2.2.

chmarti1 commented 2 years ago

Addendum... This is actually another instance of Issue #46. The default upper density for methane is 1600 kg/m3. However, densities that high are so far out of methane's practical range that pressure cannot be reasonably extrapolated that far. p() passes an inflection and becomes negative - hence the bracket problem.

This problem vanishes if

>>> ch4 = pm.get('mp.CH4')
>>> ch4.data['dlim'] = [0, 1000]
>>> ch4.d(p=.11707696, T=624.375)
array([0.03617947])
>>> ch4.d(p=9999, T=624.375)
array([502.41168232])

The bracketing function idea is still not a bad one, but establishing better density boundaries is a much easier fix.

chmarti1 commented 1 year ago

This is now resolved with new density limits in version 2.2.2.