Matlab Refprop: Incorrect solutions for mixture working fluid & Unit system issue

[Yili-Zhang]

Hi all,

The Matlab Refprop installed in my windows computer works well for pure working fluids. However, it never returns correct solutions when the working fluids is mixture. In addition, it returns different solutions if I run exact same code twice when the working fluid is mixture. For example, When the working fluid is: hFld='Methanol;0.5;Ethanol;0.5'; The first law efficiency it returns is same as when the working fluid is pure: hFld='Methanol'; And if I run the code with mixture fluids twice, the first law efficiency it returns are totally different. Just for the readers' reference, here is the complete code:

clear
clc

%% Set Up - Call RefProp in MatLAB
% Refprop MatLAB setup instruction web: https://github.com/usnistgov/REFPROP-wrappers/tree/master/wrappers/MATLAB
[v,e] = pyversion; system([e,' -m pip install --user -U ctREFPROP']);
RP = py.ctREFPROP.ctREFPROP.REFPROPFunctionLibrary('C:\Program Files (x86)\REFPROP');

% Example 1 (https://nbviewer.jupyter.org/github/usnistgov/REFPROP-wrappers/blob/master/wrappers/python/notebooks/Tutorial.ipynb)
% hFld='water';
% r = RP.REFPROPdll(hFld,'PQ','PHASE',int8(0),int8(0),int8(0),101.325,0.0,{1.0});
% o = double(r.Output);
% o(1)

% Example 2
% SI_WITH_C = RP.GETENUMdll(int8(3), "SI WITH C").iEnum;
% p_Pa = .101325;
% Q = 0.0;
% r = RP.REFPROPdll('Water','PQ','T',SI_WITH_C,int8(0),int8(0),p_Pa,Q,{1.0});
% o=double(r.Output);
% o(1)

%% Basic Cycle - Pure Fluids
MASS_SI=RP.GETENUMdll(int8(0), "MASS SI").iEnum; % Unit System: Mass SI

% Given Condition
hFld='Methanol;0.5;Ethanol;0.5'; % Working fluid
Q_dot_in=1.6e5; % [kW] Total heat into the secondary cycle
eat_t=0.85; % Turbine efficiency
eat_p=0.75; % Pump efficiency

T_max=574.15; % [K] Maximum temperature allowed in the secondary cycle
T_min=303.15; % [K] Minimum temperature allowed in the secondary cycle
P_max=12.5; % [Mpa] Maximum pressure allowed in the secondary cycle

% state 1 properties
T_1=T_max; % Temperature at state 1
P_1=P_max; % Pressure at state 1
H=RP.REFPROPdll(hFld,'TP','H',MASS_SI,int8(0),int8(0),T_1,P_1,{1.0});
H_num=double(H.Output);
H_1=H_num(1); % Enthalpy at state 1
S=RP.REFPROPdll(hFld,'TP','S',MASS_SI,int8(0),int8(0),T_1,P_1,{1.0});
S_num=double(S.Output);
S_1=S_num(1); % Entropy at state 1
Q=RP.REFPROPdll(hFld,'TP','QMASS',MASS_SI,int8(0),int8(0),T_1,P_1,{1.0});
Q_num=double(Q.Output);
Q_1=Q_num(1); % Quality at state 1
PHASE=RP.REFPROPdll(hFld,'TP','PHASE',MASS_SI,int8(0),int8(0),T_1,P_1,{1.0});
PHASE_1=convertCharsToStrings(char(PHASE.hUnits));
V=RP.REFPROPdll(hFld,'TP','V',MASS_SI,int8(0),int8(0),T_1,P_1,{1.0});
V_num=double(H.Output);
V_1=H_num(1);

% State 3 properties
T_3=T_min; % Temperature at state 3
Q_3=0; % Quality at state 3(Quality after condensing)
P=RP.REFPROPdll(hFld,'TQ','P',MASS_SI,int8(0),int8(0),T_3,Q_3,{1.0});
P_num=double(P.Output);
P_3=P_num(1); % Pressure at state 3
H=RP.REFPROPdll(hFld,'TQ','H',MASS_SI,int8(0),int8(0),T_3,Q_3,{1.0});
H_num=double(H.Output);
H_3=H_num(1); % Enthalpy at state 3
S=RP.REFPROPdll(hFld,'TQ','S',MASS_SI,int8(0),int8(0),T_3,Q_3,{1.0});
S_num=double(S.Output);
S_3=S_num(1); % Entropy at state 3
PHASE=RP.REFPROPdll(hFld,'TQ','PHASE',MASS_SI,int8(0),int8(0),T_3,Q_3,{1.0});
PHASE_3=convertCharsToStrings(char(PHASE.hUnits));
V=RP.REFPROPdll(hFld,'TQ','V',MASS_SI,int8(0),int8(0),T_3,Q_3,{1.0});
V_num=double(H.Output);
V_3=H_num(1);

% State 2 properties
P_2=P_3; % Pressure at state 2
S_2s=S_1;
H=RP.REFPROPdll(hFld,'PS','H',MASS_SI,int8(0),int8(0),P_2,S_2s,{1.0});
H_num=double(H.Output);
H_2s=H_num(1); % Isentropic enthalpy at state 2
H_2=eat_t*(H_2s-H_1)+H_1; % Enthalpy at state 2
T=RP.REFPROPdll(hFld,'PH','T',MASS_SI,int8(0),int8(0),P_2,H_2,{1.0});
T_num=double(T.Output);
T_2=T_num(1); % Temperature at state 2
S=RP.REFPROPdll(hFld,'PH','S',MASS_SI,int8(0),int8(0),P_2,H_2,{1.0});
S_num=double(S.Output);
S_2=S_num(1); % Entropy at state 2
Q=RP.REFPROPdll(hFld,'PH','QMASS',MASS_SI,int8(0),int8(0),P_2,H_2,{1.0});
Q_num=double(Q.Output);
Q_2=Q_num(1); % Quality at state 2
PHASE=RP.REFPROPdll(hFld,'PH','PHASE',MASS_SI,int8(0),int8(0),P_2,H_2,{1.0});
PHASE_2=convertCharsToStrings(char(PHASE.hUnits));
V=RP.REFPROPdll(hFld,'PH','V',MASS_SI,int8(0),int8(0),P_2,H_2,{1.0});
V_num=double(H.Output);
V_2=H_num(1);

% State 4 Properties
P_4=P_1; % Pressure at state 4
S_4s=S_3;
H=RP.REFPROPdll(hFld,'PS','H',MASS_SI,int8(0),int8(0),P_4,S_4s,{1.0});
H_num=double(H.Output);
H_4s=H_num(1); % Isentropic enthalpy at state 4
H_4=(H_4s-H_3)/eat_p+H_3; % Enthalpy at state 4
T=RP.REFPROPdll(hFld,'PH','T',MASS_SI,int8(0),int8(0),P_4,H_4,{1.0});
T_num=double(T.Output);
T_4=T_num(1); % Temperature at state 4
S=RP.REFPROPdll(hFld,'PH','S',MASS_SI,int8(0),int8(0),P_4,H_4,{1.0});
S_num=double(S.Output);
S_4=S_num(1); % Entropy at state 4
Q=RP.REFPROPdll(hFld,'PH','QMASS',MASS_SI,int8(0),int8(0),P_4,H_4,{1.0});
Q_num=double(Q.Output);
Q_4=Q_num(1); % Quality at state 4
PHASE=RP.REFPROPdll(hFld,'PH','PHASE',MASS_SI,int8(0),int8(0),P_4,H_4,{1.0});
PHASE_4=convertCharsToStrings(char(PHASE.hUnits));
V=RP.REFPROPdll(hFld,'PH','V',MASS_SI,int8(0),int8(0),P_4,H_4,{1.0});
V_num=double(H.Output);
V_4=H_num(1);

% Cycle Performance
m_dot=Q_dot_in/(H_1-H_4); % [kg/s] Secondary cycle mass flow rate
w_dot_out=m_dot*(H_1-H_2); % [kW] Turbine power out
w_dot_in=m_dot*(H_4-H_3); % [kW] Pump power in
w_dot_net=w_dot_out-w_dot_in; % [kW] Net power out of the cycle
Q_dot_out=m_dot*(H_2-H_3);  % [KW] Heat out of the cycle
eta_I=w_dot_net/Q_dot_in; % 1st law efficiency

% Export Data into Excel
Fluid=[hFld;"N/A";"N/A";"N/A"];
State=["state 1"; "state 2"; "state 3"; "state 4"];
T=[T_1; T_2; T_3; T_4];
P=[P_1; P_2; P_3; P_4];
H=[H_1; H_2; H_3; H_4];
S=[S_1; S_2; S_3; S_4];
Q=[Q_1; Q_2; Q_3; Q_4];
PHASE=[PHASE_1; PHASE_2; PHASE_3; PHASE_4];
v=[V_1;V_2;V_3;V_4];
M_dot=[m_dot;"N/A";"N/A";"N/A"];
Eta_I=[eta_I;"N/A";"N/A";"N/A"];
Table=table(Fluid,State,T,P,H,S,Q,PHASE,M_dot,Eta_I);
Table(1:4,:)

% filename='PureFluidProp.xlsx';
% writetable(Table,filename,'sheet',1,'Range','D1')

%% Comments
% 1. There is a bug on the Refprop High-level API website
% (https://refprop-docs.readthedocs.io/en/latest/DLL/high_level.html#f/_/REFPROPdll):
% When iUnits=0, it actually matches with the MASS SI unit system.
% The unit systems actually matches with the Excel dataset.'

Moreover, by following this API website: https://refprop-docs.readthedocs.io/en/latest/DLL/highlevel.html#f//REFPROPdll, the MASS SI unit system's iUnits number is 0 instead of 2 shown in the website. Actually, after a few test, it turned out that only when the iUnits number is "0", the unit system is valid.

Looking forward to hearing suggestions!

Thanks a lot,

-Yili

[ianhbell]

Nice post; you're really getting the hang of it! Can you please fill out the entire template in the future? I made your example formatted nicer, see my edits to see the changes.

Some notes:

1) We have had some issues with caching with mixture compositions, and you may have run into another bug with that. Usually the problem arises when the composition is on a mass-fraction, but in your case, you are using a mole fraction, so I don't think it is the same bug, but it could be another related bug. REFPROP is maybe getting confused because you passed a composition in the string and in the slot for z.

2) I note that you trust that the ierr codes are always zero, but you should check after every single call to REFPROPdll to ensure that ierr is either zero or less than 100. Maybe you have an error in one of the calls? Particularly the first call?

3) Do you believe the mixture model for methanol + ethanol? I'm not totally sure that I do. This mixture should be quite ideal, but the model for this system is probably too simple to yield highly accurate predictions of mixture behavior. So watch out, and be cautious with any conclusions you make based on this model.

?Methanol/Ethanol                                          [MEO/ETO]
?T.M. Blackham and E.W. Lemmon, NIST (2010)
  c41aa690/3d5f67e0
    XR0      1.             0.999          1.             1.012          0.             0.             0. 0. 0. 0. 0. 0.
    TC1     -0.013455912    0.55791296    -0.92403449     2.6625274     -2.372638       0.897343       0. 0. 0. 0. 0. 0.
    VC1      0.0015400039   0.00099950601  0.016948491   -0.054484426    0.077497449   -0.033556216    0. 0. 0. 0. 0. 0.

[Yili-Zhang]

Hi ,

Thank you very much for your quick reply!

By " fill out the entire template", do you mean fill out:

REFPROPdll(hFld, hIn, hOut, iUnits, iMass, iFlag, a, b, z, Output, hUnits, iUCode, x, y, x3, q, ierr, herr, hFld_length, hIn_length, hOut_length, hUnits_length, herr_length) ?

If that were the case, I would love to do that. However, there is not any example about it (filling out all the parameter), and it would be very challenging to me for I don't even know what data type to be used...Therefore, if you could provide me a couple examples. it would be very helpful!

Furthermore, in term of your notes, I am mostly confused:

1. I don't quite understand what you mean by "REFPROP is maybe getting confused because you passed a composition in the string and in the slot for z". If my current pass for z were incorrect. what could be the potentially correct ones? May you give me some hints or examples?

2. For checking with the ierr, again, I need a complete example for filling out the entire subroutine in order to be able to edit on that.

3. The reason why our team chose NIST is based on lots of paper reading in the thermodynamic fields, and almost all the mixture fluids are tested with NIST RefProp. So, we assumes this software is the authority in this area. However, combined with your statements, I am confused about the level of reliability of this software?

Thanks,

-Yili

[ianhbell]

1. Nope, by template, I mean that when you start a new issue on github, there are a set of common pieces of information we need, like operating system, version of REFPROP, how you are calling it, what errors you get, etc. You just deleted the entire template, rather than filling in the relevant information.

2. When you call the REFPROPdll function, a structure (class) is returned, and you can access the value of ierr, see the example here: https://nbviewer.jupyter.org/github/usnistgov/REFPROP-wrappers/blob/master/wrappers/python/notebooks/Tutorial.ipynb#Use

3. A clarification: NIST is the National Institute of Standards and Technology, REFPROP is the software tool that you are using. Mixture models are empirical, and fit to experimental data. For some mixtures, there is very limited experimental data available, and therefore, we can only build as reliable a model as the data will allow. I don't know what data were used to build the mixture model for methanol + ethanol, but you are considering a mixture model of two problematic equations of state, which is never a good starting point.

[Yili-Zhang]

Hi Lan,

Thanks for your reply. I have been busy and sorry for the late reply.

I recently encountered another problem with RefProp in MatLAB, and the debug leads to the source of issue as the reference point of enthalpies in the model are not consistent between different fluids.

From the https://refprop-docs.readthedocs.io/en/latest/DLL/highlevel.html#f//REFPROPdll website, I found the "SETREF" command under the "subroutine RefPROPdll". However, I have hard time on programming it specifically in MatLAB. Will you give me a explicit example, please?

Thanks, -Yili

[ianhbell]

Good question! I added a section to the Python tutorial: https://github.com/usnistgov/REFPROP-wrappers/blob/master/wrappers/python/notebooks/Tutorial.ipynb (see the section at the end on reference states)

[EricLemmon]

Has this issue been resolved? If so, could you close this thread?

[Yili-Zhang]

Hi Eric,

Unfortunately, no. Sorry about it. We might need more time.

-Yili

On Mon, May 13, 2019 at 7:20 AM Eric Lemmon notifications@github.com wrote:

Has this issue been resolved? If so, could you close this thread?

— You are receiving this because you authored the thread. Reply to this email directly, view it on GitHub https://github.com/usnistgov/REFPROP-wrappers/issues/109?email_source=notifications&email_token=AHEBAGLWF6V4OREM3PDI67DPVFTJTA5CNFSM4HAR2MRKYY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODVIJBMQ#issuecomment-491819186, or mute the thread https://github.com/notifications/unsubscribe-auth/AHEBAGJVE6RV6H4E64WGU2LPVFTJTANCNFSM4HAR2MRA .