We reproduce the results obtained for by shindo et al. 1985 for the cross-flow configuration and the ternary system (NH3, N2 and H2). The solution strategy is based on finite differences technique.
Table 1 of Ref. 1 lists the values of the permeabilities of H2, N2 and NH2. Feed is composed of 45 mole percent NH3, 25 mole % N2 and 30 mole % H2. The pressure ratio, g, is equal to 0.13.
In the first part, we find the stage cut (given in Table 2 of Ref. 1) when the total dimensionless area is equal to unity. While in the second part, we find the dimensionless total area (given in Table 3 of Ref. 1) when the stage cut is equal to 0.5.
Table 1 of Ref. 1 lists the values of the permeabilities of H2, N2 and NH2. Feed is composed of 45 mole percent NH3, 25 mole % N2 and 30 mole % H2. The pressure ratio, g, is equal to 0.13.
In the first part, we find the stage cut (given in Table 2 of Ref. 1) when the total dimensionless area is equal to unity. While in the second part, we find the dimensionless total area (given in Table 3 of Ref. 1) when the stage cut is equal to 0.5.
Feed composition
x
fA
x
fB
x
fC
Permeability ratios
qA=1;qB=0.317;qC=0.065;
Permeate flow rate
q
P
q
f
Feed flow rate
q
f
Pressure ratio
p
l
p
h
Total dimensionless area set equal to one
st=1;
In[]:=
eq1i=/()qAst(-/)/(qBst(-/));eq2i=/()qAst(-/)/(qCst(-/));eq3i=++1;soli=FindRoot[{eq1i,eq2i,eq3i},{,0.2},{,0.4},{,0.1}](*molefractionofthepermeateattheentranceofthemodule*)
q
P
y
PAi
q
P
y
PBi
x
fA
p
l
p
h
y
PAi
x
fB
p
l
p
h
y
PBi
q
P
y
PAi
q
P
y
PCi
x
fA
p
l
p
h
y
PAi
x
fC
p
l
p
h
y
PCi
y
PAi
y
PBi
y
PCi
y
PAi
y
PBi
y
PCi
Out[]=
{0.791261,0.165131,0.0436083}
y
PAi
y
PBi
y
PCi
Assume area or st and fetch for theta and permeate and retentate composition at outlet
In[]:=
st=.
Finite differences allow us to discretize the governing equations given in paper by Shindo et al 1985
In[]:=
Nx=20;hs=1/20.;xgrid=hsRange[0,Nx];valuesxA=Table[xA[i],{i,1,Nx+1}];valuesyA=Table[yA[i],{i,1,Nx+1}];valuesxB=Table[xB[i],{i,1,Nx+1}];valuesyB=Table[yB[i],{i,1,Nx+1}];valuesq=Table[q[i],{i,1,Nx+1}];xAp=1/stNDSolve`FiniteDifferenceDerivative[Derivative[1],xgrid,valuesxA,"DifferenceOrder"4];yAp=1/stNDSolve`FiniteDifferenceDerivative[Derivative[1],xgrid,valuesyA,"DifferenceOrder"4];xBp=1/stNDSolve`FiniteDifferenceDerivative[Derivative[1],xgrid,valuesxB,"DifferenceOrder"4];yBp=1/stNDSolve`FiniteDifferenceDerivative[Derivative[1],xgrid,valuesyB,"DifferenceOrder"4];qp=1/stNDSolve`FiniteDifferenceDerivative[Derivative[1],xgrid,valuesq,"DifferenceOrder"4];
The nodes are equally spaced in finite difference approach
In[]:=
xgrid
Out[]=
{0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.}
The governing equations y stands for permeate side mol fraction and x for retentate side mol fraction
In[]:=
eq1[i_]:=-yA[i]qp[[i]]qA(xA[i]-/yA[i])
p
l
p
h
In[]:=
eq2[i_]:=-yB[i]qp[[i]]qB(xB[i]-/yB[i])
p
l
p
h
In[]:=
eq3[i_]:=-yC[i]qp[[i]]qC(xC[i]-/yC[i])
p
l
p
h
In[]:=
eq4[i_]:=yA[i]qp[[i]]q[i]xAp[[i]]+qp[[i]]xA[i]
In[]:=
eq5[i_]:=yB[i]qp[[i]]q[i]xBp[[i]]+qp[[i]]xB[i]
Summation rules
yC[i_]:=1-yA[i]-yB[i]xC[i_]:=1-xA[i]-xB[i]
In[]:=
BC={xA[1],xB[1],yA[1](/.soli),yB[1](/.soli),q[1],q[Nx+1](1-theta)};
x
fA
x
fB
y
PAi
y
PBi
q
f
q
f
Forming the system of equations in total we have 106 equations
In[]:=
sys=Flatten@Join[Table[{eq1[i],eq2[i],eq3[i],eq4[i],eq5[i]},{i,2,Nx}],{eq4[1],eq5[Nx+1],eq1[Nx+1],eq2[Nx+1],eq3[Nx+1]},BC];
In[]:=
Length[sys]
Out[]=
106
This is the most challenging part - fetching a good initial guess for FindRoot to converge - here I adjust the initial guess using the solution obtained for the same problem using Chebyshev orthogonal collocation
In[]:=
initCHEBYSHEV={{xA[1],0.28146},{yA[1],0.651402},{xB[1],0.277548},{yB[1],0.25629},{q[1],0.627592},{xA[2],0.282556},{yA[2],0.6526},{xB[2],0.277484},{yB[2],0.255571},{q[2],0.629455},{xA[3],0.285819},{yA[3],0.656136},{xB[3],0.277282},{yB[3],0.253441},{q[3],0.635028},{xA[4],0.291172},{yA[4],0.661842},{xB[4],0.276915},{yB[4],0.249983},{q[4],0.644269},{xA[5],0.298483},{yA[5],0.669453},{xB[5],0.276344},{yB[5],0.245332},{q[5],0.657096},{xA[6],0.307565},{yA[6],0.678626},{xB[6],0.275525},{yB[6],0.239668},{q[6],0.673379},{xA[7],0.318176},{yA[7],0.688969},{xB[7],0.274422},{yB[7],0.233206},{q[7],0.692918},{xA[8],0.330024},{yA[8],0.700065},{xB[8],0.273013},{yB[8],0.226185},{q[8],0.715438},{xA[9],0.342774},{yA[9],0.711507},{xB[9],0.271298},{yB[9],0.218851},{q[9],0.740561},{xA[10],0.356065},{yA[10],0.72292},{xB[10],0.269304},{yB[10],0.211442},{q[10],0.767809},{xA[11],0.369524},{yA[11],0.733975},{xB[11],0.26708},{yB[11],0.204175},{q[11],0.796594},{xA[12],0.382785},{yA[12],0.744409},{xB[12],0.264698},{yB[12],0.197237},{q[12],0.826226},{xA[13],0.395505},{yA[13],0.754015},{xB[13],0.262243},{yB[13],0.190781},{q[13],0.855933},{xA[14],0.407376},{yA[14],0.762647},{xB[14],0.259809},{yB[14],0.184923},{q[14],0.884882},{xA[15],0.418135},{yA[15],0.770208},{xB[15],0.257487},{yB[15],0.179749},{q[15],0.912212},{xA[16],0.427562},{yA[16],0.776639},{xB[16],0.255365},{yB[16],0.175317},{q[16],0.937073},{xA[17],0.435484},{yA[17],0.781909},{xB[17],0.253521},{yB[17],0.171663},{q[17],0.958665},{xA[18],0.441769},{yA[18],0.786005},{xB[18],0.252018},{yB[18],0.16881},{q[18],0.976273},{xA[19],0.446322},{yA[19],0.788927},{xB[19],0.250909},{yB[19],0.166767},{q[19],0.9893},{xA[20],0.449077},{yA[20],0.790678},{xB[20],0.250229},{yB[20],0.16554},{q[20],0.997302},{xA[21],0.45},{yA[21],0.791261},{xB[21],0.25},{yB[21],0.165131},{q[21],1.},{theta,0.392408}};
In[]:=
initCHEBYSHEVsorted=Sort@initCHEBYSHEV;
In[]:=
initFD=Join[{{theta,initCHEBYSHEVsorted[[1]][[2]]}},Table[{q[2-i+Nx],initCHEBYSHEVsorted[[1+i]][[2]]},{i,Nx+1,1,-1}],Table[{xA[2-i+Nx],initCHEBYSHEVsorted[[22+i]][[2]]},{i,Nx+1,1,-1}],Table[{xB[2-i+Nx],initCHEBYSHEVsorted[[43+i]][[2]]},{i,Nx+1,1,-1}],Table[{yA[2-i+Nx],initCHEBYSHEVsorted[[64+i]][[2]]},{i,Nx+1,1,-1}],Table[{yB[2-i+Nx],initCHEBYSHEVsorted[[85+i]][[2]]},{i,Nx+1,1,-1}]]
Out[]=
{{theta,0.392408},{q[1],1.},{q[2],0.997302},{q[3],0.9893},{q[4],0.976273},{q[5],0.958665},{q[6],0.937073},{q[7],0.912212},{q[8],0.884882},{q[9],0.855933},{q[10],0.826226},{q[11],0.796594},{q[12],0.767809},{q[13],0.740561},{q[14],0.715438},{q[15],0.692918},{q[16],0.673379},{q[17],0.657096},{q[18],0.644269},{q[19],0.635028},{q[20],0.629455},{q[21],0.627592},{xA[1],0.45},{xA[2],0.449077},{xA[3],0.446322},{xA[4],0.441769},{xA[5],0.435484},{xA[6],0.427562},{xA[7],0.418135},{xA[8],0.407376},{xA[9],0.395505},{xA[10],0.382785},{xA[11],0.369524},{xA[12],0.356065},{xA[13],0.342774},{xA[14],0.330024},{xA[15],0.318176},{xA[16],0.307565},{xA[17],0.298483},{xA[18],0.291172},{xA[19],0.285819},{xA[20],0.282556},{xA[21],0.28146},{xB[1],0.25},{xB[2],0.250229},{xB[3],0.250909},{xB[4],0.252018},{xB[5],0.253521},{xB[6],0.255365},{xB[7],0.257487},{xB[8],0.259809},{xB[9],0.262243},{xB[10],0.264698},{xB[11],0.26708},{xB[12],0.269304},{xB[13],0.271298},{xB[14],0.273013},{xB[15],0.274422},{xB[16],0.275525},{xB[17],0.276344},{xB[18],0.276915},{xB[19],0.277282},{xB[20],0.277484},{xB[21],0.277548},{yA[1],0.791261},{yA[2],0.790678},{yA[3],0.788927},{yA[4],0.786005},{yA[5],0.781909},{yA[6],0.776639},{yA[7],0.770208},{yA[8],0.762647},{yA[9],0.754015},{yA[10],0.744409},{yA[11],0.733975},{yA[12],0.72292},{yA[13],0.711507},{yA[14],0.700065},{yA[15],0.688969},{yA[16],0.678626},{yA[17],0.669453},{yA[18],0.661842},{yA[19],0.656136},{yA[20],0.6526},{yA[21],0.651402},{yB[1],0.165131},{yB[2],0.16554},{yB[3],0.166767},{yB[4],0.16881},{yB[5],0.171663},{yB[6],0.175317},{yB[7],0.179749},{yB[8],0.184923},{yB[9],0.190781},{yB[10],0.197237},{yB[11],0.204175},{yB[12],0.211442},{yB[13],0.218851},{yB[14],0.226185},{yB[15],0.233206},{yB[16],0.239668},{yB[17],0.245332},{yB[18],0.249983},{yB[19],0.253441},{yB[20],0.255571},{yB[21],0.25629}}
In[]:=
Length[initFD]
Out[]=
106
In[]:=
Length[sys]
Out[]=
106
First part :
we find the stage cut (given in Table 2 of Ref. 1) when the total dimensionless area is equal to unity.
we find the stage cut (given in Table 2 of Ref. 1) when the total dimensionless area is equal to unity.
The solution - we get very good agreement with SHINDO 1985
In[]:=
sol=FindRoot[sys/.st1,initFD,MaxIterations10000,AccuracyGoal->20]
Out[]=
{theta0.372408,q[1]1.,q[2]0.978217,q[3]0.956742,q[4]0.93558,q[5]0.914736,q[6]0.894214,q[7]0.874019,q[8]0.854156,q[9]0.834628,q[10]0.815439,q[11]0.796594,q[12]0.778095,q[13]0.759945,q[14]0.742148,q[15]0.724706,q[16]0.707621,q[17]0.690895,q[18]0.674528,q[19]0.658522,q[20]0.642877,q[21]0.627592,xA[1]0.45,xA[2]0.442454,xA[3]0.434788,xA[4]0.427005,xA[5]0.419107,xA[6]0.411098,xA[7]0.402981,xA[8]0.394759,xA[9]0.386439,xA[10]0.378026,xA[11]0.369524,xA[12]0.360942,xA[13]0.352286,xA[14]0.343564,xA[15]0.334785,xA[16]0.325957,xA[17]0.317091,xA[18]0.308198,xA[19]0.299287,xA[20]0.29037,xA[21]0.28146,xB[1]0.25,xB[2]0.251853,xB[3]0.253685,xB[4]0.255492,xB[5]0.257272,xB[6]0.259018,xB[7]0.260726,xB[8]0.262392,xB[9]0.26401,xB[10]0.265574,xB[11]0.26708,xB[12]0.268522,xB[13]0.269892,xB[14]0.271186,xB[15]0.272396,xB[16]0.273517,xB[17]0.274542,xB[18]0.275463,xB[19]0.276276,xB[20]0.276973,xB[21]0.277548,yA[1]0.791261,yA[2]0.786447,yA[3]0.781451,yA[4]0.776264,yA[5]0.77088,yA[6]0.76529,yA[7]0.759487,yA[8]0.753462,yA[9]0.747207,yA[10]0.740715,yA[11]0.733975,yA[12]0.726982,yA[13]0.719726,yA[14]0.7122,yA[15]0.704396,yA[16]0.696308,yA[17]0.68793,yA[18]0.679254,yA[19]0.670277,yA[20]0.660994,yA[21]0.651402,yB[1]0.165131,yB[2]0.168501,yB[3]0.171982,yB[4]0.175576,yB[5]0.179288,yB[6]0.183119,yB[7]0.187074,yB[8]0.191154,yB[9]0.195363,yB[10]0.199702,yB[11]0.204175,yB[12]0.208782,yB[13]0.213525,yB[14]0.218405,yB[15]0.223421,yB[16]0.228573,yB[17]0.233859,yB[18]0.239278,yB[19]0.244826,yB[20]0.250498,yB[21]0.25629}
Theta the stage cut is equal to 0.372408
Permeate mole fractions - we recover the permeate mole fractions listed in Table 2
Second part:
we find the dimensionless total area (given in Table 3 of Ref . 1) when the stage cut is equal to 0.5.
we find the dimensionless total area (given in Table 3 of Ref . 1) when the stage cut is equal to 0.5.
st the total dimensionless area is equal to 1.47516
Permeate mole fractions - we recover the permeate mole fractions listed in Table 3
References
References
◼
[1] Y. Shindo, T. Hakuta, H. Yoshitome and H. Inoue H (1985) Calculation methods for multicomponent gas separation by permeation. Sep Sci Technol 20:445–459. https://doi.org/10.1080/01496398508060692
CITE THIS NOTEBOOK
CITE THIS NOTEBOOK
Solution of a gas permeation problem using the finite differences technique
by Housam Binous and Ahmed Bellagi
Wolfram Community, STAFF PICKS, December 25, 2024
https://community.wolfram.com/groups/-/m/t/3344456
by Housam Binous and Ahmed Bellagi
Wolfram Community, STAFF PICKS, December 25, 2024
https://community.wolfram.com/groups/-/m/t/3344456