“HumidAirProps” violates the first law of thermodynamic for isentropic process

[thnttt]

For an isentropic compression process, the relation dh=vdp can be found in any text book. When I calculate the compression work for high pressure humid air, the enthalpy obtained with function “HumidAirProps” seems to produce incorrect results. Below are the codes I used in Mathematica 9.0.1, Win7 x64.

{P1, T1, [Phi]1} = {5000, 303.15, 1};
P2 = 10000;
{W, V1, H1, S1} =
{HumidAirProps["W", "P", P1, "T", T1, "R", [Phi]1],
HumidAirProps["Vha", "P", P1, "T", T1, "R", [Phi]1],
HumidAirProps["Hha", "P", P1, "T", T1, "R", [Phi]1],
HumidAirProps["S", "P", P1, "T", T1, "R", [Phi]1]};
{T2, V2, H2, S2} =
{HumidAirProps["T", "P", P2, "S", S1, "W", W],
HumidAirProps["Vha", "P", P2, "S", S1, "W", W],
HumidAirProps["Hha", "P", P2, "S", S1, "W", W],
S1};

Ws1 = NIntegrate[
HumidAirProps["Vha", "P", p, "S", S1, "W", W], {p, P1, P2}];
Ws2 = H2 - H1;

In[415]:= {Ws1,Ws2}
Out[415]= {67.7215,70.9963}

"Ws1" is calculated by numerical integration along the isentropic path, while "Ws2" is the enthalpy difference between the starting and ending point of the compression process.
"Ws1" should equal "Ws2". However, this isn't true for any value of relative humidity [Phi]1. Besides, "Ws2" does not vary continuously at relative humidity [Phi]1=0.
The large difference between these methods confuses me a lot, and it will lead to large error in the prediction of compressor efficiency.

I tested the function "Props" for dry air, and got the right result "Ws1"="Ws2".
So I wonder If there is a bug in "HumidAirProps" or the method described in ASHRAE RP-1485.

Thanks!

[ibell]

Do you end up with the same entropies at the ends of the compression
processes when you integrate? If the process is truly isentropic, I
suppose you should. Or your numerical integration isn't very accurate. I
can't test Mathematica, so if you redid this simple example in python or
MATLAB I could help you.

On Thu, Apr 24, 2014 at 4:09 PM, thnttt [email protected] wrote:

For an isentropic compression process, the relation dh=vdp can be found in
any text book. When I calculate the compression work for high pressure
humid air, the enthalpy obtained with function “HumidAirProps” seems to
produce incorrect results. Below are the codes I used in Mathematica 9.0.1,
Win7 x64.

{P1, T1, [Phi]1} = {5000, 303.15, 1};
P2 = 10000;
{W, V1, H1, S1} =
{HumidAirProps["W", "P", P1, "T", T1, "R", [Phi]1],
HumidAirProps["Vha", "P", P1, "T", T1, "R", [Phi]1],
HumidAirProps["Hha", "P", P1, "T", T1, "R", [Phi]1],
HumidAirProps["S", "P", P1, "T", T1, "R", [Phi]1]};
{T2, V2, H2, S2} =
{HumidAirProps["T", "P", P2, "S", S1, "W", W],
HumidAirProps["Vha", "P", P2, "S", S1, "W", W],
HumidAirProps["Hha", "P", P2, "S", S1, "W", W],
S1};

Ws1 = NIntegrate[
HumidAirProps["Vha", "P", p, "S", S1, "W", W], {p, P1, P2}];
Ws2 = H2 - H1;

In[415]:= {Ws1,Ws2}
Out[415]= {67.7215,70.9963}

"Ws1" is calculated by numerical integration along the isentropic path,
while "Ws2" is the enthalpy difference between the starting and ending
point of the compression process.
"Ws1" should equal "Ws2". However, this isn't true for any value of
relative humidity [Phi]1. Besides, "Ws2" does not vary continuously at
relative humidity [Phi]1=0.
The large difference between these methods confuses me a lot, and it will
lead to large error in the prediction of compressor efficiency.

I tested the function "Props" for dry air, and got the right result
"Ws1"="Ws2".
So I wonder If there is a bug in "HumidAirProps" or the method described
in ASHRAE RP-1485.

Thanks!

—
Reply to this email directly or view it on GitHubhttps://github.com//issues/211
.

[ibell]

Another thought, the enthalpies and entropies are per kg of dry air,
perhaps that is throwing things off too

On Thu, Apr 24, 2014 at 4:14 PM, Ian Bell [email protected] wrote:

Do you end up with the same entropies at the ends of the compression
processes when you integrate? If the process is truly isentropic, I
suppose you should. Or your numerical integration isn't very accurate. I
can't test Mathematica, so if you redid this simple example in python or
MATLAB I could help you.

On Thu, Apr 24, 2014 at 4:09 PM, thnttt [email protected] wrote:

For an isentropic compression process, the relation dh=vdp can be found
in any text book. When I calculate the compression work for high pressure
humid air, the enthalpy obtained with function “HumidAirProps” seems to
produce incorrect results. Below are the codes I used in Mathematica 9.0.1,
Win7 x64.

{P1, T1, [Phi]1} = {5000, 303.15, 1};
P2 = 10000;
{W, V1, H1, S1} =
{HumidAirProps["W", "P", P1, "T", T1, "R", [Phi]1],
HumidAirProps["Vha", "P", P1, "T", T1, "R", [Phi]1],
HumidAirProps["Hha", "P", P1, "T", T1, "R", [Phi]1],
HumidAirProps["S", "P", P1, "T", T1, "R", [Phi]1]};
{T2, V2, H2, S2} =
{HumidAirProps["T", "P", P2, "S", S1, "W", W],
HumidAirProps["Vha", "P", P2, "S", S1, "W", W],
HumidAirProps["Hha", "P", P2, "S", S1, "W", W],
S1};

Ws1 = NIntegrate[
HumidAirProps["Vha", "P", p, "S", S1, "W", W], {p, P1, P2}];
Ws2 = H2 - H1;

In[415]:= {Ws1,Ws2}
Out[415]= {67.7215,70.9963}

"Ws1" is calculated by numerical integration along the isentropic path,
while "Ws2" is the enthalpy difference between the starting and ending
point of the compression process.
"Ws1" should equal "Ws2". However, this isn't true for any value of
relative humidity [Phi]1. Besides, "Ws2" does not vary continuously at
relative humidity [Phi]1=0.
The large difference between these methods confuses me a lot, and it will
lead to large error in the prediction of compressor efficiency.

I tested the function "Props" for dry air, and got the right result
"Ws1"="Ws2".
So I wonder If there is a bug in "HumidAirProps" or the method described
in ASHRAE RP-1485.

Thanks!

—
Reply to this email directly or view it on GitHubhttps://github.com//issues/211
.

[thnttt]

Thank you for your fast reply.

In my codes above, the compression process has the same entropy, and the parameters “Hha” and “Vha” are used to get the enthalpy and volume at per kg of humid air.

It seems “Sha” for entropy is invalid, so I have to use “S” per kg of dry air. I don’t think this is the cause of my problem.

Sorry I don’t have Python or MATLAB installed. Instead, I can explain my calculation procedure:
At the compressor inlet surface, pressure “P1” and temperature “T1” of humid air with relative humidity “[Phi]1” are known; so humidity ratio “W”, entropy “S1” and enthalpy “Hha1” can be obtained by “HumidAirProps”; the compressor outlet conditions are calculated with a given pressure “P2”, and the same humidity ratio “W” and entropy “S1” as those of inlet.

Let's discuss the accuracy of the numerical integration for compressor work later. Here is another interesting phenomenon I found, and maybe it should be solved first.

Below are the compressor parameters at two different relative humidity “[Phi]1”(0 and 0.0001).
Notice the severe enthalpy jump at compressor outlet (Hha, from 85.4 to 90.97). Why?

[ibell]

Lets start off with an even smaller W (say 1e-10) - though I am not sure
that is the cause. You are adding such a small amount of water, it
shouldn't be a problem... It looks like maybe it is calculating the wrong
drybulb temperature, I'll have to see why. If you then use the drybulb
temperature you obtain + p + W, do you get the right S again?

Also, Ill add Sha as an input for consistency's sake.

On Fri, Apr 25, 2014 at 4:15 AM, thnttt [email protected] wrote:

Thank you for your fast reply.

In my codes above, the compression process has the same entropy, and the
parameters “Hha” and “Vha” are used to get the enthalpy and volume at per
kg of humid air.

It seems “Sha” for entropy is invalid, so I have to use “S” per kg of dry
air. I don’t think this is the cause of my problem.

Sorry I don’t have Python or MATLAB installed. Instead, I can explain my
calculation procedure:
At the compressor inlet surface, pressure “P1” and temperature “T1” of
humid air with relative humidity “[Phi]1” are known; so humidity ratio “W”,
entropy “S1” and enthalpy “Hha1” can be obtained by “HumidAirProps”; the
compressor outlet conditions are calculated with a given pressure “P2”, and
the same humidity ratio “W” and entropy “S1” as those of inlet.

Let's discuss the accuracy of the numerical integration for compressor
work later. Here is another interesting phenomenon I found, and maybe it
should be solved first.

Below are the compressor parameters at two different relative humidity
“[Phi]1”(0 and 0.0001).
Notice the severe enthalpy jump at compressor outlet (Hha, from 85.4 to
90.97). Why?

[image: 11]https://cloud.githubusercontent.com/assets/7394920/2796962/cfed3f02-cc1d-11e3-9cbc-db64d215e3a1.jpg
[image: 22]https://cloud.githubusercontent.com/assets/7394920/2796963/cfeeb9d6-cc1d-11e3-9477-c51b5482cbfd.jpg

—
Reply to this email directly or view it on GitHubhttps://github.com//issues/211#issuecomment-41353132
.

[thnttt]

I think the problem I encountered in the isentropic compression process may be caused by the strange behavior of the Entropy subroutine in "HumidAirProps".

For the dry air at any certain condition (W=0), even adding a very small amount of water vapor (W=1e-10), "HumidAirProps" will show a severe drop of entropy (especially at high pressure and low temperature conditions).

Another interesting thing happens when W is a negative value (W = -1e-10):

In[139]:= HumidAirProps["S","P",100,"T",380,"W",-10.^-10]
Out[139]= Indeterminate

In[140]:= HumidAirProps["H","P",100,"T",380,"W",-10.^-10]
Out[140]= 107.706

The function of enthalpy and volume are continuous at W=0 and have real values when W<0.
I know W<0 is meaningless. I just wonder why the Entropy function is so special.

[ibell]

Interesting. I'll have to take a look and see what is going on. Thanks
for the very clear description of the problem.

On Fri, Apr 25, 2014 at 1:27 PM, thnttt [email protected] wrote:

I think the problem I encountered in the isentropic compression process
may be caused by the strange behavior of the Entropy subroutine in
"HumidAirProps".

For the dry air at any certain condition (W=0), even adding a very small
amount of water vapor (W=1e-10), "HumidAirProps" will show a severe drop of
entropy (especially at high pressure and low temperature conditions).

[image: 33]https://cloud.githubusercontent.com/assets/7394920/2799983/80af98d4-cc6a-11e3-9ac6-264487759db5.jpg
[image: 44]https://cloud.githubusercontent.com/assets/7394920/2799984/80b06bf6-cc6a-11e3-902f-a12505ec6b59.jpg

Another interesting thing happens when W is a negative value (W = -1e-10):

In[139]:= HumidAirProps["S","P",100,"T",380,"W",-10.^-10]
Out[139]= Indeterminate

In[140]:= HumidAirProps["H","P",100,"T",380,"W",-10.^-10]
Out[140]= 107.706

The function of enthalpy and volume are continuous at W=0 and have real
values when W<0.
I know W<0 is meaningless. I just wonder why the Entropy function is so
special.

—
Reply to this email directly or view it on GitHubhttps://github.com//issues/211#issuecomment-41382728
.

[thnttt]

Is there any good news about this problem?

I still need to calculate the isentropic compression work, but the strange results yielded by "HumidAirProps" worries me a lot. REFPROP doesn't support humid air, and I don't think the ideal gas assumption is valid for humid air above 5MPa.

So please help me! Thank you very much!

[ibell]

Sorry at the moment, no, I am working hard on a new version of CoolProp and
haven't had the chance to take a look yet. Can you do a plot for fixed T,P
with varied W? I'm curious to see whether it is truly discontinuous or just
a very sharp change.

On Wed, May 21, 2014 at 10:43 AM, thnttt [email protected] wrote:

Is there any good news about this problem?

I still need to calculate the isentropic compression work, but the strange
results yielded by "HumidAirProps" worries me a lot. REFPROP doesn't
support humid air, and I don't think the ideal gas assumption is valid for
humid air above 5MPa.

So please help me! Thank you very much!

—
Reply to this email directly or view it on GitHubhttps://github.com//issues/211#issuecomment-43728981
.

[thnttt]

Thank you for your attention! Here are the entropy and enthalpy varying with W.

When using "Props" for dry air, the two ways to calculate the isentropic compression work (described above) produce consistent results. However, "HumidAirProps" gives incorrect result even when W=0.

[ibell]

Interesting. The HAProps definition of the enthalpy is not the same as for
the pure fluid - you can see that by reading
http://rp.ashrae.biz/page/ASHRAE-D-RP-1485-20091216.pdf . It is not
entirely surprising behavior really I guess. Perhaps there is still an
error here that needs to be remedied though

On Wed, May 21, 2014 at 1:27 PM, thnttt [email protected] wrote:

Thank you for your attention! Here are the entropy and enthalpy varying
with W.

[image: 1]https://cloud.githubusercontent.com/assets/7394920/3039133/4b0f13c6-e0d8-11e3-8204-bc89a5536fdf.jpg
[image: 2]https://cloud.githubusercontent.com/assets/7394920/3039134/4b348ef8-e0d8-11e3-95d7-7db9ce7347d7.jpg

When using "Props" for dry air, the two ways to calculate the isentropic
compression work (described above) produce consistent results. However,
"HumidAirProps" gives incorrect result even when W=0.

—
Reply to this email directly or view it on GitHubhttps://github.com//issues/211#issuecomment-43742782
.

[ibell]

The basic problem is that if you read that document, they have terms that
are log(psi_w), which is -infinity at psi_w = 0, psi_w being the mole
fraction of water vapor. We have a conditional in our code (
https://github.com/CoolProp/CoolProp/blob/master/src/HumidAirProp.cpp#L774-L778)
to avoid this problem. I have emailed the authors before but we didn't
come to a nice solution.

On Wed, May 21, 2014 at 1:39 PM, Ian Bell [email protected] wrote:

Interesting. The HAProps definition of the enthalpy is not the same as
for the pure fluid - you can see that by reading
http://rp.ashrae.biz/page/ASHRAE-D-RP-1485-20091216.pdf . It is not
entirely surprising behavior really I guess. Perhaps there is still an
error here that needs to be remedied though

On Wed, May 21, 2014 at 1:27 PM, thnttt [email protected] wrote:

Thank you for your attention! Here are the entropy and enthalpy varying
with W.

[image: 1]https://cloud.githubusercontent.com/assets/7394920/3039133/4b0f13c6-e0d8-11e3-8204-bc89a5536fdf.jpg
[image: 2]https://cloud.githubusercontent.com/assets/7394920/3039134/4b348ef8-e0d8-11e3-95d7-7db9ce7347d7.jpg

When using "Props" for dry air, the two ways to calculate the isentropic
compression work (described above) produce consistent results. However,
"HumidAirProps" gives incorrect result even when W=0.

—
Reply to this email directly or view it on GitHubhttps://github.com//issues/211#issuecomment-43742782
.

[thnttt]

Thanks. Now I know the discontinuity at W=0 is an inherent character of the Entropy model.

Here is the performance data of a compressor stage calculated with HAProps :

Notice the abnormal isentropic efficiency, which should be less than the polytropic efficiency.
Would you give me some advice on how to obtain a reasonable isentropic efficiency for high pressure humid air. Thank you!