Skip to content

Commit 0ef4292

Browse files
simulkadesimulkade
committed
refactoring the code to object oriented matlab and octave
1 parent 76b0d07 commit 0ef4292

29 files changed

Lines changed: 2148 additions & 705 deletions

‎Classes/@BIP/BIP.m‎

Lines changed: 81 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,81 @@
1+
classdef BIP
2+
% BIP Binary Interaction Parameters container.
3+
%
4+
% bip = BIP(n) creates an n-component BIP object with all parameter
5+
% matrices initialised to zero. Set individual matrices after
6+
% construction to provide non-zero interaction parameters.
7+
%
8+
% Properties (all n-by-n matrices unless noted):
9+
% EOScons, EOStdep - EOS kij parameters
10+
% NRTLcons, NRTLtdep, NRTLtdep2 - NRTL: A + B*T + C*T^2
11+
% NRTLtdepm1, NRTLtdeplog - NRTL: D/T + E*ln(T) terms
12+
% NRTLalfa - NRTL non-randomness parameter
13+
% Wilsoncons, Wilsontdep - Wilson parameters
14+
% UNIQUACcons, UNIQUACtdep - UNIQUAC interaction parameters
15+
% UNIQUACR (1-by-n) - UNIQUAC volume parameters
16+
% UNIQUACQ (1-by-n) - UNIQUAC surface parameters
17+
18+
%{
19+
Copyright (c) 2012, 2013, Ali Akbar Eftekhari
20+
All rights reserved.
21+
22+
Redistribution and use in source and binary forms, with or
23+
without modification, are permitted provided that the following
24+
conditions are met:
25+
26+
* Redistributions of source code must retain the above copyright notice,
27+
this list of conditions and the following disclaimer.
28+
* Redistributions in binary form must reproduce the above
29+
copyright notice, this list of conditions and the following
30+
disclaimer in the documentation and/or other materials provided
31+
with the distribution.
32+
33+
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
34+
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,
35+
THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
36+
PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR
37+
CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
38+
EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
39+
PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,
40+
OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
41+
LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
42+
NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
43+
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
44+
%}
45+
46+
properties
47+
EOScons
48+
EOStdep
49+
NRTLcons
50+
NRTLtdep
51+
NRTLtdep2
52+
NRTLtdepm1
53+
NRTLtdeplog
54+
NRTLalfa
55+
Wilsoncons
56+
Wilsontdep
57+
UNIQUACcons
58+
UNIQUACtdep
59+
UNIQUACR
60+
UNIQUACQ
61+
end
62+
63+
methods
64+
function obj = BIP(n)
65+
obj.EOScons = zeros(n);
66+
obj.EOStdep = zeros(n);
67+
obj.NRTLcons = zeros(n);
68+
obj.NRTLtdep = zeros(n);
69+
obj.NRTLtdep2 = zeros(n);
70+
obj.NRTLtdepm1 = zeros(n);
71+
obj.NRTLtdeplog = zeros(n);
72+
obj.NRTLalfa = zeros(n);
73+
obj.Wilsoncons = zeros(n);
74+
obj.Wilsontdep = zeros(n);
75+
obj.UNIQUACcons = zeros(n);
76+
obj.UNIQUACtdep = zeros(n);
77+
obj.UNIQUACR = zeros(1, n);
78+
obj.UNIQUACQ = zeros(1, n);
79+
end
80+
end
81+
end

‎Classes/@Component/Component.m‎

Lines changed: 198 additions & 86 deletions
Original file line numberDiff line numberDiff line change
@@ -1,7 +1,47 @@
11
classdef Component
2-
%UNTITLED Summary of this class goes here
3-
% Detailed explanation goes here
4-
2+
% Component Pure-component thermodynamic properties.
3+
%
4+
% Typical usage: load from database
5+
% comp = Component.fromDatabase('CH4');
6+
% comps = Component.fromDatabaseArray({'CH4', 'C2H6'});
7+
%
8+
% Or construct directly:
9+
% comp = Component(name, formula, MW, Tc, Pc, Vc, Zc,
10+
% acentric_factor, Psat_eq, Psat_coefs,
11+
% PsatTrange, Psatrange,
12+
% dh_vap_eq, dh_vap_coefs, dh_vap_Trange, dh_vap_range,
13+
% cp_liq_eq, cp_liq_coefs, cp_liq_Trange, cp_liq_range,
14+
% cp_ig_eq, cp_ig_coefs, cp_ig_Trange, cp_ig_range,
15+
% dhf_ig, dgf_ig, ds_ig, dh_comb)
16+
17+
%{
18+
Copyright (c) 2012, 2013, Ali Akbar Eftekhari
19+
All rights reserved.
20+
21+
Redistribution and use in source and binary forms, with or
22+
without modification, are permitted provided that the following
23+
conditions are met:
24+
25+
* Redistributions of source code must retain the above copyright notice,
26+
this list of conditions and the following disclaimer.
27+
* Redistributions in binary form must reproduce the above
28+
copyright notice, this list of conditions and the following
29+
disclaimer in the documentation and/or other materials provided
30+
with the distribution.
31+
32+
THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
33+
AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO,
34+
THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
35+
PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR
36+
CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
37+
EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
38+
PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA,
39+
OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
40+
LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
41+
NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
42+
SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
43+
%}
44+
545
properties
646
name
747
formula
@@ -32,98 +72,170 @@
3272
ds_ig
3373
dh_comb
3474
end
35-
75+
3676
methods
37-
function obj = Component(formula, ...
38-
MW, ...
39-
Tc, ...
40-
Pc, ...
41-
Vc, ...
42-
Zc, ...
43-
acentric_factor, ...
44-
Psat_eq, ...
45-
Psat_coefs, ...
46-
PsatTrange, ...
47-
Psatrange, ...
48-
dh_vap_eq, ...
49-
dh_vap_coefs, ...
50-
dh_vap_Trange, ...
51-
dh_vap_range, ...
52-
cp_liq_eq, ...
53-
cp_liq_coefs, ...
54-
cp_liq_Trange, ...
55-
cp_liq_range, ...
56-
cp_ig_eq, ...
57-
cp_ig_coefs, ...
58-
cp_ig_Trange, ...
59-
cp_ig_range, ...
60-
dhf_ig, ...
61-
dgf_ig, ...
62-
ds_ig, ...
63-
dh_comb)
64-
%Component Construct an instance of this class
65-
% Creates a Component object; currently it can be created from
66-
% a database from the Perry's handbook
67-
obj.name = name;
68-
obj.formula = formula;
69-
obj.MW = MW;
70-
obj.Tc = Tc;
71-
obj.Pc = Pc;
72-
obj.Vc = Vc;
73-
obj.Zc = Zc;
77+
function obj = Component(name, formula, ...
78+
MW, Tc, Pc, Vc, Zc, acentric_factor, ...
79+
Psat_eq, Psat_coefs, PsatTrange, Psatrange, ...
80+
dh_vap_eq, dh_vap_coefs, dh_vap_Trange, dh_vap_range, ...
81+
cp_liq_eq, cp_liq_coefs, cp_liq_Trange, cp_liq_range, ...
82+
cp_ig_eq, cp_ig_coefs, cp_ig_Trange, cp_ig_range, ...
83+
dhf_ig, dgf_ig, ds_ig, dh_comb)
84+
% Zero-argument call is used by Octave/MATLAB to pre-allocate
85+
% object arrays (e.g. obj(n) = Component(...)).
86+
if nargin == 0; return; end
87+
obj.name = name;
88+
obj.formula = formula;
89+
obj.MW = MW;
90+
obj.Tc = Tc;
91+
obj.Pc = Pc;
92+
obj.Vc = Vc;
93+
obj.Zc = Zc;
7494
obj.acentric_factor = acentric_factor;
75-
obj.Psat_eq = Psat_eq;
76-
obj.Psat_coefs = Psat_coefs;
77-
obj.PsatTrange = PsatTrange;
78-
obj.Psatrange = Psatrange;
79-
obj.dh_vap_eq = dh_vap_eq;
80-
obj.dh_vap_coefs = dh_vap_coefs;
81-
obj.dh_vap_Trange = dh_vap_Trange;
82-
obj.dh_vap_range = dh_vap_range;
83-
obj.cp_liq_eq = cp_liq_eq;
84-
obj.cp_liq_coefs = cp_liq_coefs;
85-
obj.cp_liq_Trange = cp_liq_Trange;
86-
obj.cp_liq_range = cp_liq_range;
87-
obj.cp_ig_eq = cp_ig_eq;
88-
obj.cp_ig_coefs = cp_ig_coefs;
89-
obj.cp_ig_Trange = cp_ig_Trange;
90-
obj.cp_ig_range = cp_ig_range;
91-
obj.dhf_ig = dhf_ig;
92-
obj.dgf_ig = dgf_ig;
93-
obj.ds_ig = ds_ig;
94-
obj.dh_comb = dh_comb;
95+
obj.Psat_eq = Psat_eq;
96+
obj.Psat_coefs = Psat_coefs;
97+
obj.PsatTrange = PsatTrange;
98+
obj.Psatrange = Psatrange;
99+
obj.dh_vap_eq = dh_vap_eq;
100+
obj.dh_vap_coefs = dh_vap_coefs;
101+
obj.dh_vap_Trange = dh_vap_Trange;
102+
obj.dh_vap_range = dh_vap_range;
103+
obj.cp_liq_eq = cp_liq_eq;
104+
obj.cp_liq_coefs = cp_liq_coefs;
105+
obj.cp_liq_Trange = cp_liq_Trange;
106+
obj.cp_liq_range = cp_liq_range;
107+
obj.cp_ig_eq = cp_ig_eq;
108+
obj.cp_ig_coefs = cp_ig_coefs;
109+
obj.cp_ig_Trange = cp_ig_Trange;
110+
obj.cp_ig_range = cp_ig_range;
111+
obj.dhf_ig = dhf_ig;
112+
obj.dgf_ig = dgf_ig;
113+
obj.ds_ig = ds_ig;
114+
obj.dh_comb = dh_comb;
95115
end
96-
97-
function p_sat = vapor_pressure(obj,T)
98-
%vapor_pressure(T)
99-
% p_sat = component.vapor_pressure(T)
100-
% T is temperature in K
101-
% p_sat is pressure in Pa
116+
117+
function p_sat = vapor_pressure(obj, T)
118+
% vapor_pressure(T) returns p_sat [Pa] at temperature T [K]
102119
p_sat = obj.Psat_eq(T, obj.Psat_coefs);
103120
end
104-
105-
function dh_v = dh_vap(obj,T)
106-
%dh_vap(T)
107-
% dh_v [J/mol] = component.dh_vap(T)
108-
% T is temperature in K
121+
122+
function dh_v = dh_vap(obj, T)
123+
% dh_vap(T) returns enthalpy of vaporisation [J/mol] at T [K]
109124
dh_v = obj.dh_vap_eq(T/obj.Tc, obj.dh_vap_coefs);
110125
end
111-
112-
function cp_l = cp_liq(obj,T)
113-
%cp_liq(T)
114-
% cp_l [J/(mol.K)] = component.cp_liq(T)
115-
% liquid phase specific heat capacity
116-
% T is temperature in K
126+
127+
function cp_l = cp_liq(obj, T)
128+
% cp_liq(T) returns liquid heat capacity [J/(mol·K)] at T [K]
117129
cp_l = obj.cp_liq_eq(T/obj.Tc, obj.cp_liq_coefs);
118130
end
119-
120-
function cp = cp_ig(obj,T)
121-
%cp_ig(T)
122-
% Ideal gas heat capacity
123-
% cp [J/(mol.K)] = component.cp_ig(T)
124-
% T is temperature in K
125-
cp = obj.cp_ig(T, obj.cp_ig_coefs);
131+
132+
function cp = cp_ig(obj, T)
133+
% cp_ig(T) returns ideal-gas heat capacity [J/(mol·K)] at T [K]
134+
cp = obj.cp_ig_eq(T, obj.cp_ig_coefs);
126135
end
127136
end
128-
end
129137

138+
methods (Static)
139+
function comp = fromDatabase(name_or_formula)
140+
% fromDatabase Load a single Component from puredata.mat.
141+
% comp = Component.fromDatabase('CH4')
142+
% comp = Component.fromDatabase('methane')
143+
comps = Component.fromDatabaseArray({name_or_formula});
144+
if isempty(comps)
145+
error('Component:notFound', ...
146+
'Component "%s" not found in database.', name_or_formula);
147+
end
148+
comp = comps(1);
149+
end
150+
151+
function [components, flag] = fromDatabaseArray(name_formula_cell)
152+
% fromDatabaseArray Load multiple Components from puredata.mat.
153+
% [comps, flag] = Component.fromDatabaseArray({'CH4','C2H6'})
154+
% flag contains indices of names that were not found.
155+
load puredata.mat
156+
N = length(name_formula_cell);
157+
data_index = zeros(1, N);
158+
for n = 1:N
159+
idx = find(strcmpi(component_formula, name_formula_cell(n)), 1);
160+
if ~isempty(idx)
161+
data_index(n) = idx;
162+
else
163+
idx = find(strcmpi(component_name, name_formula_cell(n)), 1);
164+
if ~isempty(idx)
165+
data_index(n) = idx;
166+
end
167+
end
168+
end
169+
flag = find(data_index == 0);
170+
data_index = data_index(data_index > 0);
171+
Nfound = length(data_index);
172+
if Nfound == 0
173+
components = [];
174+
return
175+
end
176+
% Pre-allocate by constructing the last element first so MATLAB
177+
% knows the array type, then fill in order.
178+
i = Nfound;
179+
ci = data_index(i);
180+
if liquid_heat_capacity_equation_num(ci) == 1
181+
cp_liq_fn = liquid_heat_capacity_equation_1;
182+
else
183+
cp_liq_fn = liquid_heat_capacity_equation_2;
184+
end
185+
if ideal_gas_heat_capacity_equation_num(ci) == 1
186+
cp_ig_fn = ideal_gas_heat_capacity_equation_1;
187+
elseif ideal_gas_heat_capacity_equation_num(ci) == 2
188+
cp_ig_fn = ideal_gas_heat_capacity_equation_2;
189+
else
190+
cp_ig_fn = ideal_gas_heat_capacity_equation_3;
191+
end
192+
components(Nfound) = Component( ...
193+
component_name{ci}, component_formula{ci}, ...
194+
Molecular_weight_data(ci), ...
195+
critical_temperature_data(ci), critical_pressure_data(ci), ...
196+
critical_volume_data(ci), critical_compressibility_factor(ci), ...
197+
acentric_factor(ci), ...
198+
vapor_pressure_equation, vapor_pressure_coefs(ci,:), ...
199+
vapor_pressure_T_range(ci,:), vapor_pressure_p_range(ci,:), ...
200+
dh_vaporization_equation, dh_vaporization_coefs(ci,:), ...
201+
dh_vaporization_T_range(ci,:), dh_vaporization_dhv_range(ci,:), ...
202+
cp_liq_fn, liquid_heat_capacity_coef(ci,:), ...
203+
liquid_heat_capacity_T_range(ci,:), liquid_heat_capacity_cp_range(ci,:), ...
204+
cp_ig_fn, ideal_gas_heat_capacity_coefs(ci,:), ...
205+
ideal_gas_heat_capacity_T_range(ci,:), ideal_gas_heat_capacity_cp_range(ci,:), ...
206+
dh_formation_ig(ci), dg_formation_ig(ci), ...
207+
ds_ideal_gas(ci), dh_combustion(ci));
208+
for n = 1:Nfound-1
209+
ci = data_index(n);
210+
if liquid_heat_capacity_equation_num(ci) == 1
211+
cp_liq_fn = liquid_heat_capacity_equation_1;
212+
else
213+
cp_liq_fn = liquid_heat_capacity_equation_2;
214+
end
215+
if ideal_gas_heat_capacity_equation_num(ci) == 1
216+
cp_ig_fn = ideal_gas_heat_capacity_equation_1;
217+
elseif ideal_gas_heat_capacity_equation_num(ci) == 2
218+
cp_ig_fn = ideal_gas_heat_capacity_equation_2;
219+
else
220+
cp_ig_fn = ideal_gas_heat_capacity_equation_3;
221+
end
222+
components(n) = Component( ...
223+
component_name{ci}, component_formula{ci}, ...
224+
Molecular_weight_data(ci), ...
225+
critical_temperature_data(ci), critical_pressure_data(ci), ...
226+
critical_volume_data(ci), critical_compressibility_factor(ci), ...
227+
acentric_factor(ci), ...
228+
vapor_pressure_equation, vapor_pressure_coefs(ci,:), ...
229+
vapor_pressure_T_range(ci,:), vapor_pressure_p_range(ci,:), ...
230+
dh_vaporization_equation, dh_vaporization_coefs(ci,:), ...
231+
dh_vaporization_T_range(ci,:), dh_vaporization_dhv_range(ci,:), ...
232+
cp_liq_fn, liquid_heat_capacity_coef(ci,:), ...
233+
liquid_heat_capacity_T_range(ci,:), liquid_heat_capacity_cp_range(ci,:), ...
234+
cp_ig_fn, ideal_gas_heat_capacity_coefs(ci,:), ...
235+
ideal_gas_heat_capacity_T_range(ci,:), ideal_gas_heat_capacity_cp_range(ci,:), ...
236+
dh_formation_ig(ci), dg_formation_ig(ci), ...
237+
ds_ideal_gas(ci), dh_combustion(ci));
238+
end
239+
end
240+
end
241+
end

0 commit comments

Comments
 (0)