Chapter 8 : Phase Equilibria III Phase Diagrams¶
Example 8.1 Page No : 470¶
In [5]:
from scipy.optimize import fsolve
import math
# Variables
A_C5H12 = 9.2131 ; #From table E8.2A
B_C5H12 = 2477.07 ; #From table E8.2A
C_C5H12 = -39.94 ; #From table E8.2A
A_C6H12 = 9.1325 ; #From table E8.2A
B_C6H12 = 2766.63 ; #From table E8.2A
C_C6H12 = -50.50 ; #From table E8.2A
A_C6H14 = 9.2164 ; #From table E8.2A
B_C6H14 = 2697.55 ; #From table E8.2A
C_C6H14 = -48.78 ; #From table E8.2A
A_C7H16 = 9.2535 ; #From table E8.2A
B_C7H16 = 2911.32 ; #From table E8.2A
C_C7H16 = -56.51 ; #From table E8.2A
x_C5H12 = 0.3 ;
x_C6H12 = 0.3 ;
x_C6H14 = 0.2 ;
x_C7H16 = 0.2 ;
# Calculations
def f82(T):
return -1 + (x_C5H12 * math.exp(A_C5H12 - B_C5H12 / (T + C_C5H12)) + x_C6H12 *\
math.exp(A_C6H12 - B_C6H12 / (T + C_C6H12)) + x_C6H14 * \
math.exp(A_C6H14 - B_C6H14 / (T + C_C6H14)) + x_C5H12 * math.exp(A_C5H12 - B_C5H12 / (T \
+ C_C5H12)) + x_C7H16 * math.exp(A_C7H16 - B_C7H16 / (T + C_C7H16)));
y =fsolve(f82,300)
# Results
print " The temperature at which the liquid develops the first bubble of vapour = %d K"%(y);
# Note : answer may vary because of fsolve function of python.
Example 8.2 Page No : 471¶
In [6]:
from scipy.optimize import fsolve
import math
# Variables
A_C5H12 = 9.2131 ; #From table E8.2A
B_C5H12 = 2477.07 ; #From table E8.2A
C_C5H12 = -39.94 ; #From table E8.2A
A_C6H12 = 9.1325 ; #From table E8.2A
B_C6H12 = 2766.63 ; #From table E8.2A
C_C6H12 = -50.50 ; #From table E8.2A
A_C6H14 = 9.2164 ; #From table E8.2A
B_C6H14 = 2697.55 ; #From table E8.2A
C_C6H14 = -48.78 ; #From table E8.2A
A_C7H16 = 9.2535 ; #From table E8.2A
B_C7H16 = 2911.32 ; #From table E8.2A
C_C7H16 = -56.51 ; #From table E8.2A
y_C5H12 = 0.3 ;
y_C6H12 = 0.3 ;
y_C6H14 = 0.2 ;
y_C7H16 = 0.2 ;
P = 1 ; #[bar]
# Calculations and Results
def f83(T):
return -1 + P * ( y_C5H12 / math.exp(A_C5H12 - B_C5H12 / (T +
C_C5H12)) + y_C6H12 / math.exp(A_C6H12 - B_C6H12 / (T + C_C6H12)) + \
y_C6H14 / math.exp(A_C6H14 - B_C6H14 / (T + C_C6H14)) + y_C7H16 / math.exp(A_C7H16 - B_C7H16 / (T + C_C7H16)));
y =fsolve(f83,300)
print "The temperature at which vapour develops the first drop of liquid = %.f K"%(y) ;
T = y ;
P_sat_C5H12 = math.exp(A_C5H12 - B_C5H12 / (T + C_C5H12)) ;
p_sat_C6H12 = math.exp(A_C6H12 - B_C6H12 / (T + C_C6H12)) ;
P_sat_C6H14 = math.exp(A_C6H14 - B_C6H14 / (T + C_C6H14)) ;
P_sat_C7H16 = math.exp(A_C7H16 - B_C7H16 / (T + C_C7H16)) ;
x_C5H12 = y_C5H12 * P / P_sat_C5H12 ;
x_C6H12 = y_C6H12 * P / p_sat_C6H12 ;
x_C6H14 = y_C6H14 * P / P_sat_C6H14 ;
x_C7H16 = y_C7H16 * P / P_sat_C7H16 ;
print "x_C5H12 = %.3f x_C6H12 = %.3f x_C6H14 = %.3f x_C7H16 = %.3f"%(x_C5H12,x_C6H12,x_C6H14,x_C7H16) ;
Example 8.3 Page No : 471¶
In [8]:
from scipy.integrate import *
from numpy import *
import math
# Variables
P_a_sat = 0.53 ; #[bar]
P_b_sat = 0.16 ; #[bar]
X = 1./3 ;
Y = 1- X ;
x_a_feed = 0.5 ;
x_b_feed = 0.5 ;
# Calculations
a = Y * -(x_a_feed + x_b_feed) + Y**2 ;
b = X * Y *(P_a_sat + P_b_sat) - (x_a_feed * P_b_sat + x_b_feed * P_a_sat)*X ;
c = P_a_sat * P_b_sat * X**2;
k=poly1d([a,b,c])#'k',0);
#P = c + b*k**1 + a*k**2 ;
M = roots(k);
# Results
for i in range(2):
ans = sign(M[0]) ;
if ans == 1:
print " Pressure = %.2f bar"%(M[i]) ;
Xa = x_a_feed / (P_a_sat / M[0] * X + Y) ; #....E8.4D
Ya = Xa * P_a_sat / M[0] #.....E8.4B
print " Xa = %.2f Ya = %.1f"%(Xa,Ya);
break
Example 8.5 Page No : 476¶
In [5]:
from scipy.optimize import fsolve
import math
# Variables
P1_sat = 0.72 ; #[bar]
P2_sat = 0.31 ; #[bar]
A = 3590. ;
B = -1180. ;
R = 8.314 ;
T = 70. + 273 ; #[K]
# Calculations
def f85(x1):
return -.48 + ( x1 * math.exp((A + 3*B) * (1 - x1)**2 / (R * T) -
4 * B * (1 - x1)**3 / (R * T)) * P1_sat) / ( x1 * math.exp((A + 3*B) *
(1 - x1)**2 / (R * T) -4 * B * (1 - x1)**3 / (R * T)) * P1_sat +(1 - x1) *
math.exp((A - 3*B) * x1**2 / (R * T) -4 * B * x1**3 / (R * T)) * P2_sat ) ;
y = fsolve(f85,0.1)[0]
x1 = y ;
P = ( x1 * math.exp((A + 3*B) * (1 - x1)**2 / (R * T) - 4 * B *
(1 - x1)**3 / (R * T)) * P1_sat) + (1 - x1) * math.exp((A - 3*B) * x1**2 /
(R * T) -4 * B * x1**3 / (R * T)) * P2_sat ;
# Results
print "The value of x1 = %.3f"%(y) ;
print "Pressure = %.2f bar"%(P) ;
Example 8.7 Page No : 488¶
In [6]:
# Variables
import math
P = 0.223 ; #[bar]
P_a_sat = 0.156 ; # [bar]
P_b_sat = 0.124 ; #[bar]
R = 8.314 ;
T = 50 + 273 ;
Xa = 0.554 ;
Xb = 1 - Xa ;
# Calculations
gama_a = P / P_a_sat ;
A1 = R * T * math.log(gama_a) / (Xb**2) * 10**-3 ;
gama_b = P / P_b_sat ;
A2 = R * T * math.log(gama_b) / (Xa**2) * 10**-3 ;
A = math.ceil((A1 + A2) / 2) ;
# Results
print "Value of two suffix Marguels parameter = %.1f kJ/mol"%(A);
Example 8.9 Page No : 491¶
In [10]:
# Variables
import math
R = 8.314 ;
T = 10 + 273 ; #[K]
A_B = 9.2806 ; # From Appendix A , Table A1.1
B_B = 2788.5 ; # From Appendix A , Table A1.1
C_B = -52.36 ; # From Appendix A , Table A1.1
A_C = 9.1325 ; # From Appendix A , Table A1.1
B_C = 2766.63 ; # From Appendix A , Table A1.1
C_C = -50.50 ; # From Appendix A , Table A1.1
x1 = [0 ,0.0610 ,0.2149 ,0.3187 ,0.4320 ,0.5246 ,0.6117 ,0.7265 ,0.8040 ,0.8830 ,0.8999 ,1] ; #From table E8.9A
P_exp = [6344 ,6590 ,6980 ,7140 ,7171 ,7216 ,7140 ,6974 ,6845 ,6617 ,6557 ,6073] ; #From table E8.9A
P_1_sat = 6072.15 ; #[Pa]
P_2_sat = 6344 ; #[Pa]
A = [1390 ,1391 ,1392 ,1393 ,1394 ,1395 ,1396 ,1397 ,1398 ,1399 ,1400 ,1401 ,1402 ,1403 ,1404 ,1405 ,1406 ,1407 ,1408 ,1409 ,1410 ] ;
C = zeros([21,12])
for k in range(21):
y = A[k] ;
P = zeros(12)
for i in range(12):
P[i] = x1[i] * math.exp( y / (R * T ) * (1 - x1[i])**2) * P_1_sat+(1 - x1[i]) * math.exp(y / (R * T ) * x1[i]**2) * P_2_sat ;
C[k,i] = (P[i] - P_exp[i])**2 ;
R = zeros(21)
for k in range(21):
y = 0 ;
for i in range(12):
y = y + C[k,i] ;
R[k] = y ;
k = 100000. ;
for i in range(21):
K = R[i] ;
if K < k:
k = K ;
for i in range(21):
if R[i] == k :
print "The two suffix Margules co-efficient is = %g J/mol" %(A[i]) ;
#Note : answer may vary because of rounding off error.
Example 8.10 Page No : 494¶
In [3]:
# Variables
import math
R = 8.314 ;
T = 10 + 273.15 #[K]
A_B = 9.2806 ; # From Appendix A , Table A1.1
B_B = 2788.5 ; # From Appendix A , Table A1.1
C_B = -52.36 ; # From Appendix A , Table A1.1
A_C = 9.1325 ; # From Appendix A , Table A1.1
B_C = 2766.63 ; # From Appendix A , Table A1.1
C_C = -50.50 ; # From Appendix A , Table A1.1
x1 = [0,0.0610 ,0.2149 ,0.3187 ,0.4320 ,0.5246 ,0.6117 ,0.7265 ,0.8040 ,0.8830 ,0.8999 ,1] ; #From table E8.9A
P_exp = [6344 ,6590 ,6980 ,7140 ,7171 ,7216 ,7140 ,6974 ,6845 ,6617 ,6557 ,6073] ; #From table E8.9A
P_1_sat = 6073. ; #[Pa]
P_2_sat = 6344. ; #[Pa]
A = linspace(1390,1410,21) ;
B = linspace(60,80,21)
w = 1 / (R * T) ;
S = zeros([21,21])
for k in range(21):
y = A[k] ;
for i in range(21):
z = B[i] ;
P = zeros(12)
R = zeros(12)
for j in range(12):
P[j] = x1[j] * math.exp((y + 3 * z) * (1 - (x1[j]))**2 *w-4*z*(1-x1[j])**3* w )* P_1_sat + (1-x1[j])*exp((y -3*z)*(x1[j])**2 * w + 4 * z * (x1[j]**3) * w )*P_2_sat ;
R[j] =(P[j] - P_exp[j])**2 ;
m = 0 ;
for l in range(12):
m = m + R[l]
S[k,i] = m ;
D = zeros(21)
for i in range(21):
k = S[i,0]
for l in range(1,21):
if S[i,l] < k:
k = S[i,l]
D[i] = k ;
a = D[0] ;
for i in range(2,21):
if D[i] < a:
a = D[i] ;
# Results
for i in range(21):
if D[i] == a :
for l in range(21):
if S[i,l] == a:
print " A = %g J/mol"%(A[i]);
print " B = %g J/mol"%(B[l]) ;
Example 8.11 Page No : 494¶
In [1]:
%matplotlib inline
from matplotlib.pyplot import *
# Given
R = 8.314 ;
T = 10 + 273.15 ; #[K]
x1 = [0 ,0.0610 ,0.2149 ,0.3187 ,0.4320 ,0.5246 ,0.6117 ,0.7265 ,0.8040 ,0.8830 ,0.8999 ,1] ; #From table E8.9A
P_exp = [6344 ,6590 ,6980 ,7140 ,7171 ,7216 ,7140 ,6974 ,6845 ,6617 ,6557 ,6073 ,6073] ; #From table E8.9A
y1 = [ 1 ,0.0953 ,0.2710 ,0.3600 ,0.4453,0.5106 ,0.5735 ,0.6626 ,0.7312 ,0.8200 ,0.8382, 0 ] ; #From table E8.9A
P_1_sat = 6073. ; #[Pa]
P_2_sat = 6344. ; #[Pa]
n = 0 ;
x2 = zeros(11)
y2 = zeros(11)
g_E = zeros(11)
ydata = zeros(11)
xdata = zeros(11)
for i in range(1,11):
x2[i] = 1 - x1[i] ;
y2[i] = 1 - y1[i] ;
g_E[i] = R * T *( x1[i] * math.log (( y1[i] * P_exp[i]) / (x1[i]* P_1_sat)) + x2[i] * math.log((y2[i] * P_exp[i]) / (x2[i] * P_2_sat)) ) ;
n = n + g_E[i] / ((x1[i] * x2[i]) * 10) ;
ydata[i-1] = (g_E[i]/(x1[i]*x2[i]));
xdata[i-1] = x1[i] - x2[i] ;
m= 0
n=0
o = 0
p= 0
N = 10
for i in range(1,11):
m = m + g_E[i] * (2 * x1[i] - 1) / ( x1[i] * x2[i]) ;
n = n + g_E[i] / ( x1[i] * x2[i]) ;
o = o + (2 * x1[i] - 1) ;
p = p + (2 * x1[i] - 1)**2 ;
x_bar = o / N ;
y_bar = n / N ;
a1 = (N * m - n * o)/(N * p - o**2) ;
a0 = y_bar - a1 * x_bar ;
ydata2 = zeros(11)
for i in range(11):
ydata2[i] = a0 + a1*xdata[i] ;
# Results
plot(xdata,ydata,"+") ;
plot(xdata,ydata2) ;
xlabel("x1-x2")
ylabel("g**E/x1x2 [J/mol]")
#suptitle("Fit of g E /x 1 x 2 vs. x 1 2 x 2 for benzene (1)–cyclohexane (2) at 10°C to straight line")
ylim(1320,1500)
print " From average, the value of A = %d J/mol"%(n/10) ;
print " From linear regression best fit line the values of A and B are \
%.1f J/mol & %.1f J/mol respectively ."%(a0, a1) ;
#Readers can refer figure E8.11 .
show()
Example 8.12 Page No : 500¶
In [15]:
# Variables
H_O2 = 44253.9 #[bar] , From table 8.1
p_O2 = 0.21 #[bar]
# Calculations
x_O2 = p_O2 / H_O2 ;
v_H2O = 1/(1/0.001 * 1/0.018 * 0.001 );
_O2_ = x_O2 / v_H2O ; #[M]
# Results
print "Mole fraction of O2 = %.2e"%(x_O2 ) ;
print "Concentration of O2 = %.2e M "%(_O2_) ;
Example 8.13 Page No : 500¶
In [13]:
# Variables
import math
P = 300. ; #[bar]
V_bar_inf_N2 = 3.3 * 10**-5 ;
R = 8.314 ;
T = 298. ; #[K]
y_N2 = 1. ; # At 25*C vapour pressure of water is small
H_N2_1 = 87365. ; #[bar]
P_c = 33.8 ; #[bar]
T_c = 126.2 ; # [K]
w = 0.039 ; # From Appendix A.1
log_w_0 = 0.013 ;
log_w_1 = 0.210 ;
# Calculations
H_N2_300 = H_N2_1 * math.exp((V_bar_inf_N2 * (P -1) * 10**5 )/ (R * T)) ;
k = log_w_0 + w * log_w_1 ;
sai_N2 = 10**k ;
x_N2 = y_N2 * sai_N2 * P / H_N2_300 ;
# Results
print "Solubility of N2 in water = %.5f"%(x_N2) ;
Example 8.20 Page No : 518¶
In [2]:
%matplotlib inline
import math
from matplotlib.pyplot import *
from numpy import *
# Variables
R = 8.314 ;
T = 20. + 273 ; #[K]
A = 6000. ; #[J/mol]
B = -384. ; #[J/mol]
x_a = [0.001 ,0.03 ,0.05 ,0.06 ,0.075 ,0.1 ,0.12 , 0.13 ,0.15 ,0.2 ,0.25 ,0.3 ,0.35 ,0.4 ,0.45,0.475 ,0.5 ,0.55 ,0.6 ,0.65 ,0.7 ,0.75 ,0.8 ,0.8475 ,0.85 ,0.9 ,0.925 ,0.95 ,0.975 ,0.999] ;
y_data = zeros(30)
y_data2 = zeros(30)
# Calculations
for i in range(30):
y_data[i] = R * T * ( x_a[i] * math.log(x_a[i]) + (1 - x_a[i]) * math.log(1- x_a[i])) + x_a[i] * (1 - x_a[i]) * (A + B * (2*x_a[i] - 1 )) ;
y_data2[i] = - 82 * x_a[i]- 185.6 ;
m = min(y_data) ;
for i in range(30):
if y_data[i] == m:
a = x_a[i] ;
for i in range(30):
y_data2[i] = -(R * T *( math.log(a) - math.log(1 - a)) + A * (1 - 2*a) + B * (6 * a - 1 - 6 * a**2)) * (x_a[i] - a) + m ;
y_data3 = zeros(20)
for i in range(20):
y_data3[i] = y_data[i] - y_data2[i] ;
n = min(y_data3) ;
for i in range(20):
if y_data3[i] == n:
b = x_a[i] ;
# Results
plot(x_a ,y_data) ;
plot(x_a ,y_data2) ;
xlabel("Xa")
ylabel("g-XaGa-XbGb")
suptitle("the binary system described ")
print "The equilibrium composition can be found by drawing a line tangent to\
the minima . In this case the answer is %.2f and %.1f ." %(a ,b)
show()
Example 8.21 Page No : 519¶
In [15]:
from scipy.optimize import fsolve
import math
# Variables
A = 6349 ; #[J/mol]
B = -384. ; #[J/mol]
R = 8.314 ;
T = 20 + 273 ; #[K]
k = 0.000001 ;
# Calculations
def f816(x_a):
return R * T * (1/x_a + 1/(1 - x_a)) - 2 * A +6 * B * (1 - 2 * x_a) + k
ans1 = fsolve(f816,0.1)[0]
ans2 = fsolve(f816,0.5)[0]
print "%.3f < x_a < %.3f "%(ans1,ans2)
Example 8.22 Page No : 522¶
In [1]:
from scipy.optimize import fsolve
# Variables
import math
T = 300. #[K]
A = 6235. #[J/mol]
P_a_sat = 100. * 10**3 ; #[Pa]
P_b_sat = 50. * 10**3 ; #{Pa}
R = 8.314 ;
w = 1./(R * T)
def f817(R):
x_a_a = R[0]
x_a_b = R[1]
Z817 = [0,0]
Z817[0] = x_a_b * math.exp(A * (1 - x_a_b) ** 2 * w) - x_a_a * math.exp(A*(1-x_a_a)**2*w)
Z817[1] = (1 - x_a_b) * math.exp(A * ( x_a_b) ** 2 * w) - (1 - x_a_a) * math.exp(A * (x_a_a) ** 2 * w )
return Z817
x0 = [0.75 , 0.1] ;
z = fsolve(f817,x0)
# Results
print "The compositions are : x_a_a = %.3f and x_a_b = %.3f"%(z[0],z[1]) ;
P = z[0] * math.exp(A * z[1] ** 2 * w) * P_a_sat + z[1] * math.exp(A * z[0] ** 2 * w) * P_b_sat ;
print "Total pressure = %d kPa"%(P * 10**-3) ;
y_a = z[0] * math.exp(A * z[1] ** 2 * w) * P_a_sat / P ;
print "y_a = %.3f" %( y_a ) ;
# Note : incorrect answer in textbook
Example 8.23 Page No : 535¶
In [17]:
# Variables
T_b = 373.15 ; #[K]
del_h_vap = 2257. ; #[J/g]
MW_salt = 58.5 ; #[g/mol]
MW_water = 18. ; #[g/mol]
w_salt = 3.5 ;
w_water = 100 - w_salt ;
R = 8.314 ;
# Calculations
x_salt = (w_salt / MW_salt) / (w_salt / MW_salt + w_water / MW_water) ;
x_b = 2 * x_salt
del_T = R * T_b**2 / (del_h_vap * MW_water)* x_b ;
# Results
print "The temperature that sea water boils is = %.2f degreeC"%(100 + del_T);
Example 8.24 Page No : 538¶
In [18]:
# Variables
rho_w = 1000. # [kg/m**3]
g = 9.8 ; # [m/s**2]
h = 0.0071 ; #[m]
m_b = 1.93 * 10**-3 ; # [kg]
V = 520. * 10**-6 ; #[m**3]
R = 8.314
T = 298.
# Calculations
PI = rho_w * g * h ;
C_b = m_b / V ;
MW_b = R * T * C_b / PI ;
# Results
print "The molecular weight of the protein = %d kg/mol"%( MW_b );