CalebBell / chemicals

chemicals: Chemical database of Chemical Engineering Design Library (ChEDL)
MIT License
186 stars 36 forks source link

Incorrect calculation of density when using saturation temperature as input #30

Closed thoerbr closed 3 years ago

thoerbr commented 3 years ago

The IAPWS functions return an erroneous density value when calculatiog with the saturation temperature.

import numpy as np
from chemicals import iapws95_Psat, iapws95_Tsat, iapws95_rho

pressure = 558497.3519367648
tsat = iapws95_Tsat(pressure)
iapws95_rho(tsat, pressure) # density wrong
iapws95_rho(np.floor(tsat), np.floor(pressure)) # rounding down leads to proper density

I guess this error is invoked due to some boundary issues with IAPWS. If it is too much of an effort to make it work at the saturation temperature exactly, an error message with a hint would be welcome.

Thanks for the great package anyway!

CalebBell commented 3 years ago

Hi thoerbr, There will always be a transition between liquid and vapor, because that is what every chemical does. Whether the calculation returns a liquid or vapor density at the exact saturation temperature is left up the the solver - it will converge to one or another. I confirmed this is exactly what is happening with the original two points you posted. This is expected behavior; and the iapws95_Tsat solver will never be as tight as iapws95_rho either.

My solvers are quite precise. When you know for sure which density you are looking for, it is up to you to make sure you are providing conditions which have the correct solution. In this case you can increase the pressure a trivial amount.

iapws95_rho(tsat, pressure*(1+1e-14))

Sincerely, Caleb

thoerbr commented 3 years ago

Makes sense. I didn't think about the fact, that above saturation temperature the library would return gas phase properties. Thanks for the quick answer.