Adiabatic Flash Calculations using the Peng Robinson Equation of State
Authors: Housam Binous and Ahmed Bellagi
Adiabatic Flash Calculations using the Peng Robinson Equation of State
Authors: Housam Binous and Ahmed Bellagi
Authors: Housam Binous and Ahmed Bellagi
H . Binous, Applications of the Peng - Robinson Equation of State using Mathematica®, Chemical Engineering Education, 42 (1), pp . 47 - 51, 2008.
Acknowledgement :
This problem was presented by Henley and Seader in Equilibrium-Stage Separation Operations in Chemical Engineering, Wiley, 1981.
This problem was presented by Henley and Seader in Equilibrium-Stage Separation Operations in Chemical Engineering, Wiley, 1981.
In[]:=
Off[General::"spell1"]
Reduced Pressure and Temperature
In[]:=
Tre[i_]:=T/Tc[i]
In[]:=
Pre[i_]:=P/Pc[i]
Equilibrium constants using the ϕ-ϕ method and the Peng Robinson Equation of State. Binary interaction parameters for the Peng Robinson EOS were taken from HYSYS 3.2.
In[]:=
ap[i_]:=0.45724a[i]R^2Tc[i]^2/Pc[i]
In[]:=
Ap[i_]:=ap[i]P/(RT)^2;
In[]:=
bp[i_]:=0.07780RTc[i]/Pc[i]
In[]:=
Bp[i_]:=Pbp[i]/(RT)
In[]:=
ϕ
v
Z
v
B
v
Z
v
B
v
A
v
2
B
v
4
∑
j=1
A
v
B
v
Z
v
2
B
v
Z
v
2
B
v
In[]:=
ϕ
l
Z
l
B
l
Z
l
B
l
A
l
2
B
l
4
∑
j=1
A
l
B
l
Z
l
2
B
l
Z
l
2
B
l
In[]:=
a[i_]:=(1+m[i](1-Tre[i]^0.5))^2
In[]:=
m[i_]:=0.37464+1.54226ω[i]-0.26992ω[i]^2
In[]:=
A
v
4
∑
i=1
4
∑
j=1
In[]:=
Ab[i_,j_]=(Ap[i]Ap[j])^0.5(1-k[i,j]);
In[]:=
B
v
4
∑
i=1
In[]:=
A
l
4
∑
i=1
4
∑
j=1
In[]:=
B
l
4
∑
i=1
In[]:=
K[i_]:=[i]/[i];
ϕ
l
ϕ
v
In[]:=
k[1,1]=0;k[2,2]=0;k[3,3]=0;k[4,4]=0;k[2,1]=k[1,2];k[3,1]=k[1,3];k[4,1]=k[1,4];k[2,3]=k[3,2];k[3,4]=k[4,3];k[2,4]=k[4,2];k[1,2]=0.202000007033348;k[1,3]=0.285100013017654;k[1,4]=0.285100013017654;k[3,2]=3.9999000728130310^-2;k[4,2]=6.4900003373622910^-2;k[4,3]=9.5191010041162410^-4;
In[]:=
R=1.987;
Critical Pressure and Temperature and acentric factor for Hydrogen, Methane, Benzene and Toluene
In[]:=
Pc[1]=190.8;Tc[1]=59.7;ω[1]=0.0;
In[]:=
Pc[2]=673.1;Tc[2]=343.9;ω[2]=0.0;
In[]:=
Pc[3]=714.2;Tc[3]=1012.7;ω[3]=0.2116;
In[]:=
Pc[4]=587.8;Tc[4]=1069.1;ω[4]=0.2415;
Feed Pressure and temperature are 485 psia and 100°F. Feed is composed of 0.41% Hydrogen, 5.71% Methane, 70.97% Benzene and 22.91% Toluene. Feed pressure is decreased adiabatically from 485 psia to 165 psia.
In[]:=
P=165;
In[]:=
z[1]:=0.0041;z[2]=0.0571;z[3]=0.7097;z[4]=0.2291;
In[]:=
z[1]+z[2]+z[3]+z[4]
Out[]=
1.
Adiabatic Flash Calculations using Enthalpy balance with Q=0.
In[]:=
Eq[i_]:=y[i]K[i]x[i]
In[]:=
V=F;L=(1-)F;
ϕ
vap
ϕ
vap
In[]:=
EQ[i_]:=x[i]z[i]/(1+(K[i]-1));
ϕ
vap
In[]:=
AC[1,1]=29.088;AC[2,1]=-0.19210^-2;AC[3,1]=0.410^-5;AC[4,1]=-0.8710^-9;AC[1,2]=19.875;AC[2,2]=5.02110^-2;AC[3,2]=1.26810^-5;AC[4,2]=-11.00410^-9;AC[1,3]=-36.193;AC[2,3]=48.44410^-2;AC[3,3]=-31.54810^-5;AC[4,3]=77.57310^-9;AC[1,4]=-34.364;AC[2,4]=55.88710^-2;AC[3,4]=-34.43510^-5;AC[4,4]=80.33510^-9;
In[]:=
HL1=AC[1,i]((T-459.67-32)/1.8+273.15)x[i]+AC[2,i]((T-459.67-32)/1.8+273.15)^2/2x[i]+AC[3,i]((T-459.67-32)/1.8+273.15)^3/3x[i]+AC[4,i]((T-459.67-32)/1.8+273.15)^4/4x[i]-AC[1,i](298.15)x[i]+AC[2,i](298.15)^22x[i]+AC[3,i](298.15)^33x[i]+AC[4,i](298.15)^44x[i];
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
In[]:=
HL2=1.987T(-1)+121.987TPLog+1++1-(TD[(1.987T)^2/P,{T,1}]-(1.987T)^2/P);
Z
l
2
B
l
Z
l
2
B
l
Z
l
2
B
l
A
l
A
l
In[]:=
HL2=%/.P165;
In[]:=
HL2=%/9.48610^4/453.593;
In[]:=
HL=HL1+HL2;
In[]:=
HV1=AC[1,i]((T-459.67-32)/1.8+273.15)y[i]+AC[2,i]((T-459.67-32)/1.8+273.15)^2/2y[i]+AC[3,i]((T-459.67-32)/1.8+273.15)^3/3y[i]+AC[4,i]((T-459.67-32)/1.8+273.15)^4/4y[i]-AC[1,i](298.15)y[i]+AC[2,i](298.15)^22y[i]+AC[3,i](298.15)^33y[i]+AC[4,i](298.15)^44y[i];
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
In[]:=
HV2=1.987T(-1)+121.987TPLog+1++1-(TD[(1.987T)^2/P,{T,1}]-(1.987T)^2/P);
Z
v
2
B
v
Z
v
2
B
v
Z
v
2
B
v
A
v
A
v
In[]:=
HV2=%/.P165;
In[]:=
HV2=%/9.48610^4/453.593;
In[]:=
HV=HV1+HV2;
In[]:=
(100-32)/1.8+273.15
Out[]=
310.928
In[]:=
apf[i_]:=0.45724af[i]R^2Tc[i]^2/Pc[i]
In[]:=
Apf[i_]:=apf[i]Pf/(RT)^2;
In[]:=
bpf[i_]:=0.07780RTc[i]/Pc[i]
In[]:=
Bpf[i_]:=Pfbpf[i]/(RT)
In[]:=
af[i_]:=(1+m[i](1-Tre[i]^0.5))^2
In[]:=
m[i_]:=0.37464+1.54226ω[i]-0.26992ω[i]^2
In[]:=
AF=z[i]z[j]AFb[i,j];
4
∑
i=1
4
∑
j=1
In[]:=
AFb[i_,j_]=(Apf[i]Apf[j])^0.5(1-k[i,j]);
In[]:=
BF=z[i]Bpf[i];
4
∑
i=1
In[]:=
solZF=Solve[^3+(-1+BF)^2+(AF-3BF^2-2BF)-AFBF+BF^2+BF^30/.T100+459.67/.Pf485,]//Chop
Z
F
Z
F
Z
F
Z
F
Out[]=
{{0.116934},{0.392537-0.962899},{0.392537+0.962899}}
Z
F
Z
F
Z
F
In[]:=
HF1=AC[1,i](310.928)z[i]+AC[2,i](310.928)^2/2z[i]+AC[3,i](310.928)^3/3z[i]+AC[4,i](310.928)^4/4z[i]-AC[1,i](298.15)z[i]+AC[2,i](298.15)^22z[i]+AC[3,i](298.15)^33z[i]+AC[4,i](298.15)^44z[i]
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
4
∑
i=1
Out[]=
1102.76
In[]:=
HF2=1.987T(-1)+12+1++1-
Z
F
2
BF1.987TPfLogZ
F
2
BFZ
F
2
BF(TD[AF(1.987T)^2/Pf,{T,1}]-AF(1.987T)^2/Pf);In[]:=
HF2=%/.Pf485/.T->(100+459.67)/.solZF[[1]];
In[]:=
HF2=%/9.48610^4/453.593
Out[]=
-31016.
In[]:=
HF=HF1+HF2
Out[]=
-29913.3
In[]:=
sol=FindRootEq[1],Eq[2],Eq[3],y[i]1,EQ[1],EQ[2],EQ[3],x[i]1,0,^3+(-1+)^2+(-3^2-2)-+^2+^30,^3+(-1+)^2+(-3^2-2)-+^2+^30,HFHV+(1-)HL,{x[1],0.0004},{x[2],0.03},{x[3],0.7},{x[4],0.25},{,0.01},{y[1],0.1},{y[2],0.85},{y[3],0.02},{y[4],0.001},{,0.05},{,0.95},{T,540},MaxIterations1000//Chop
4
∑
i=1
4
∑
i=1
4
∑
i=1
z[i](K[i]-1)
1+(K[i]-1)
ϕ
vap
Z
l
B
l
Z
l
Z
l
A
l
B
l
B
l
A
l
B
l
B
l
B
l
Z
v
B
v
Z
v
Z
v
A
v
B
v
B
v
A
v
B
v
B
v
B
v
ϕ
vap
ϕ
vap
ϕ
vap
Z
l
Z
v
Results of the iterative technique of Henley and Seader:
Equilibrium liquid mole fractions: { 0.000413642, 0.0280863, 0.734359, 0.237141 } ,
Equilibrium vapor mole fractions: { 0.108406, 0.869565, 0.0197101, 0.00231884 } ,
Only 3.44% of the feed is vaporized.
All these results compare very well with the approach using FindRoot presented here.
Equilibrium liquid mole fractions: { 0.000413642, 0.0280863, 0.734359, 0.237141 } ,
Equilibrium vapor mole fractions: { 0.108406, 0.869565, 0.0197101, 0.00231884 } ,
Only 3.44% of the feed is vaporized.
All these results compare very well with the approach using FindRoot presented here.