import math
# Variables
T = 0 + 273.15; #[K] - Temperature
P = 20*10**(5); #[Pa] - Pressure
R = 8.314; #[J/mol*K] - Universal gas constant
#componenet 1 : methane (1)
#componenet 2 : methanol (2)
H_constant = 1022; #[bar] - Henry's law constant
H_constant = H_constant*10**(5); #[Pa]
# The second virial coefficients are
B_11 = -53.9; #[cm**(3)/mol]
B_11 = B_11*10**(-6); #[m**(3)/mol]
B_12 = -166; #[cm**(3)/mol]
B_12 = B_12*10**(-6); #[m**(3)/mol]
B_22 = -4068; #[cm**(3)/mol]
B_22 = B_22*10**(-6); #[m**(3)/mol]
den_meth = 0.8102; #[g/cm**(3)] - Density of methanol at 0 C
Mol_wt_meth = 32.04; # Molecular weight of methanol
P_2_sat = 0.0401; #[bar] - Vapour pressure of methanol at 0 C
# Calculations
#The molar volume of methanol can be calculated as
V_2_liq = (1/(den_meth/Mol_wt_meth))*10**(-6); #[m**(3)/mol]
#The phase equilibrium equation of the components at high pressure
#y1*phi_1*P = x_1*H_1
#y2*phi_2*P = x_2*H_2
#Since methane follows Henry's law therefore methanol follows the lewis-Rnadall rule
#f_2 is the fugacity of the compressed liquid which is calculated using
#f_2 = f_2_sat*math.exp[V_2_liq*(P - P_sat_2)/(R*T)]
#where f_2_sat can be calculated using virial equation
# math.log(phi_2_sat) = math.log(f_2_sat/P_2_sat) = (B_22*P_2_sat)/(R*T)
f_2_sat = P_2_sat*math.exp((B_22*P_2_sat*10**(5))/(R*T)); #[bar]
#Putting the value of 'f_2_sat' in the math.expression of f_2 , we get
f_2 = f_2_sat*math.exp(V_2_liq*(P - P_2_sat*10**(5))/(R*T)); #[bar]
#Now let us calculate the fugacity coefficients of the species in the vapour mixture
del_12 = 2*B_12 - B_11 - B_22; #[m**(3)/mol]
#math.log(phi_1) = (P/(R*T))*(B_11 + y2**(2)*del_12)
#math.log(phi_2) = (P/(R*T))*(B_22 + y1**(2)*del_12)
#The calculation procedure is to assume a value of y1, calculate 'phi_1' and 'phi_2' and calculate 'x_1' and 'x_2' from the phase equilibrium equations and see whether x_1 + x_2 = 1,if not then another value of y1 is assumed.
y2 = 0.1;
error=10;
while(error>0.001):
y1=1-y2;
phi_1 = math.exp((P/(R*T))*((B_11 + y2**(2)*del_12)));
phi_2 = math.exp((P/(R*T))*((B_22 + y1**(2)*del_12)));
x_1 = (y1*phi_1*P)/H_constant;
x_2 = (y2*phi_2*P)/(f_2*10**(5));
x = x_1 + x_2;
error=abs(1-x);
y2=y2 - 0.000001;
# Results
print " The solubility of methane in methanol is given by x1 = %f"%(x_1);
import math
from scipy.integrate import quad
# Variables
x_C2H6_1 = 0.33*10**(-4); # Solubility of ethane in water at 25 C and 1 bar
#componenet 1 : ethane (1)
#componenet 2 : water (2)
# Z = 1 - 7.63*10**(3)*P - 7.22*10**(-5)*P**(2)
#The phase equilibrium equation of ethane is
#f_1_V = x_1*H_1
#since vapour is pure gas, f_1_V = x_1*H_1 or, phi_1*P = x_1*H_1, where 'phi_1' is fugacity coefficient of pure ethane
# math.log(phi) = integral('Z-1)/P) from limit '0' to 'P'
P1 = 0;
P2 = 1;
P3 = 35;
# Calculations
def f51(P):
return (1-7.63*10**(-3)*P-7.22*10**(-5)*P**(2)-1)/P
intgral = quad(f51,P1,P2)[0]
phi_1_1 = math.exp(intgral); # - Fugacity coefficient of ethane at 1 bar
f_1_1 = phi_1_1*P2; #[bar] - Fugacity of ethane at 1 bar
#Similarly
def f52(P):
return (1-7.63*10**(-3)*P-7.22*10**(-5)*P**(2)-1)/P
intgral_1 = quad(f52,P1,P3)[0]
phi_1_35 = math.exp(intgral_1); # Fugacity coefficient of ethane at 35 bar
f_1_35 = phi_1_35*P3; #[bar] - Fugacity of ethane at 35 bar
# At ethane pressure of 1 bar , x_C2H6_1*H_1 = phi_1_1
H_1 = phi_1_1/x_C2H6_1; #[bar] - Henry's constant
# At ethane pressure of 35 bar , x_C2H6_35*H_1 = phi_1_35
x_C2H6_35 = f_1_35/H_1; # Solubility of ethane at 35 bar pressure
# Results
print "The solubility of ethane at 35 bar is given by x_C2H6 = %e "%(x_C2H6_35);
from scipy.optimize import fsolve
import math
# Variables
T = 200; #[K]
R = 8.314; #[J/mol*K] - universal gas constant
# G_E = A*x_1*x_2
A = 4000; #[J/mol]
x_1 = 0.6; # Mle fraction of feed composition
# Calculations and Results
# Since A is given to be independent of temperature
UCST = A/(2*R); #[K] - Upper critical solution temperature
print " The UCST of the system is %f K"%(UCST);
# Since the given temperature is less than UCST therefore two phase can get formed at the given temperature.
# x1_alpha = 1 - x1_beta
# We know that, x1_alpha*Y_1_alpha = x2_alpha*Y_2_alpha
# x1_alpha*math.exp[(A/(R*T))*(x2_alpha)**(2)] = (1 - x1_alpha)*math.exp[(A/(R*T))*(x1_alpha)**(2)]
# where use has beeen made of the fact that x1_alpha = 1 - x1_beta and x2_beta = 1 - x1_beta = x1_alpha .Taking math.logarithm of both side we get
# math.log(x1_alpha) + (A/(R*T))*(1 - x1_alpha)**(2) = math.log(1 - x1_alpha) + (A/(R*T))*x1_alpha**(2)
# math.log(x1_alpha/(1-x1_alpha)) = (A/(R*T))*(2*x1_alpha - 1)
def f(x1_alpha):
return math.log(x1_alpha/(1-x1_alpha)) - (A/(R*T))*(2*x1_alpha - 1)
x1_alpha = fsolve(f,0.1)
x1_beta = fsolve(f,0.9)
# Because of symmetry 1 - x1_beta = x1_alpha
# It can be seen that the equation, math.log(x1/(1-x1)) = (A/(R*T))*(2*x1 - 1) has two roots.
print " The composition of two liquid phases in equilibrium is given by x1_alpha = %f and x1_beta = %f"%(x1_alpha,x1_beta);
z1 = 0.6;
# z1 = x1_alpha*alpha + x1_beta*beta
# beta = 1 - alpha
alpha = (z1-x1_beta)/(x1_alpha-x1_beta); #[mol]
Beta = 1 - alpha; #[mol]
print " The relative amount of phases is given by alpha = %f mol and beta = %f mol"%(alpha,Beta);
# the relative amounts of the phases changes with the feed composition
#math.log(x1/(1-x1)) = (A/(R*T))*(2*x1 - 1)
print " math.logx1/1-x1 = A/R*T*2*x1 - 1";
print " If the above equation has two real roots of x1 one for phase alpha and the other for phase beta then two liquid phases get formed";
print " and if it has no real roots then a homogeneous liquid mixture is obtained";
# Variables
T = 300; #[K]
R = 8.314; #[J/mol*K] - universal gas constant
A = 7000; #[J/mol]
# math.log(x_1/(1-x_1)) = (A/(R*T))*(2*x_1-1)
# Calculations
def f(x_1):
return math.log(x_1/(1-x_1))-((A/(R*T))*(2*x_1-1))
x1_alpha=fsolve(f,0.1)
x1_beta=1-x1_alpha;
# Results
print "The equilibrium compositin of the two liquid phase system is given by x1_alpha \t = %f x1_beta \t = %f"%(x1_alpha,x1_beta);
R = 8.314; #[J/mol*K] - Universal gas constant
M_wt_meth = 32; # Molecular weight of methanol
M_wt_water = 18; # Molecular weight of water
m_meth = 0.01; #[g] - Mass of methanol added per cm**(3) of solution
den_sol = 1; #[g/cm**(3)]
# Calculations
#The mole fraction of solute is given by
#x_2 = (moles of solute in cm**(3) of solution)/(moles of solute + moles of water) in 1 cm**(3) of solution
x_2 = (m_meth/M_wt_meth)/((m_meth/M_wt_meth)+((1-m_meth)/M_wt_water));
#We know that heat of fusion of water is
H_fus = -80; #[cal/g] - Enthalpy change of fusion at 0 C
H_fus = H_fus*4.186*M_wt_water; #[J/mol]
#Therefore freezing point depression is given by
# T - T_m = (R*(T**(2))*x_2)/H_fus
T_f = 273.15; #[K] - Freezing point of water
delta_T_f = (R*(T_f**(2))*x_2)/H_fus; #[K]
# Results
print "The depression in freezing point is given by delta_T = %f K"%(delta_T_f);
# Variables
R = 8.314; #[J/mol*K] - universal gas constant
T_f = 273.15; #[K] - Freezing point of water
m_water = 100; #[g] - Mass of water
m_NaCl = 3.5; #[g] - Mass of NaCl
M_wt_water = 18.015; # Molecular weight of water
M_wt_NaCl = 58.5; # Molecular weight of NaCl
mol_water = m_water/M_wt_water; #[mol] - Moles of water
mol_NaCl = m_NaCl/M_wt_NaCl; #[mol] - Moles of NaCl
H_fus = -80; #[cal/g] - Enthalpy change of fusion at 0 C
H_fus = H_fus*4.186*M_wt_water; #[J/mol]
# Calculations
#Mole fraction of the solute (NaCl) is given by
x_2 = mol_NaCl/(mol_NaCl+mol_water);
x_2_act = 2*x_2; # Actual mole fraction
#Now depression in freezing point is given by
# T - T_m = (R*(T**(2))*x_2_act)/H_fus
delta_T_f = (R*(T_f**(2))*x_2_act)/H_fus; #[C]
#Thus freezing point of seawater = depression in freezing point
# Results
print "The freezing point of seawater is %f C"%(delta_T_f);
# Variables
R = 8.314; #[J/mol*K] - universal gas constant
T_b = 373.15; #[K] - Boiling point of water
m_water = 100.; #[g] - Mass of water
m_C12H22 = 5.; #[g] - Mass of glucise (C12H22)
M_wt_water = 18.015; # Molecular weight of water
M_wt_C12H22 = 342.30; # Molecular weight of C12H22
mol_water = m_water/M_wt_water; #[mol] - Moles of water
mol_C12H22 = m_C12H22/M_wt_C12H22; #[mol] - Moles of C12H22
# Calculations
H_vap = 540.; #[cal/g] - Enthalpy change of vaporisation
H_vap = H_vap*4.186*M_wt_water; #[J/mol]
#Mole fraction of the solute (C12H22) is given by
x_2 = mol_C12H22/(mol_C12H22+mol_water);
#The boiling point elevation is given by
# T - T_b = (R*T_b**(2)*x_2**(2))/H_vap**(2)
delta_T_b = (R*T_b**(2)*x_2)/(H_vap);
# Results
print "The elevation in boiling point is given by delta_T = %f C"%(delta_T_b)
# Variables
R = 8.314; #[J/mol*K] - Universal gas constant
T = 25 + 273.15; #[K] - Surrounding temperature
den_water = 1000.; #[kg/m**(3)] - Density of water
m_water = 100.; #[g] - Mass of water
m_C12H22 = 5.; #[g] - Mass of glucise (C12H22)
M_wt_water = 18.015; # Molecular weight of water
M_wt_C12H22 = 342.30; # Molecular weight of C12H22
mol_water = m_water/M_wt_water; #[mol] - Moles of water
mol_C12H22 = m_C12H22/M_wt_C12H22; #[mol] - Moles of C12H22
# Calculations
#Mole fraction of the water is given by
x_1 = mol_water/(mol_C12H22+mol_water);
#Molar volume of water can be calculated as
V_l_water = (1./den_water)*M_wt_water*10**(-3); #[m**(3)/mol]
#The osmotic pressure is given by
pi = -(R*T*math.log(x_1))/V_l_water; #[N/m**(2)]
pi = pi*10**(-5); #[bar]
# Results
print "The osmotic pressure of the mixture is %f bar"%(pi);
# Variables
R = 8.314; #[J/mol*K] - universal gas constant
T = 25 + 273.15; #[K] - Surrounding temperature
den_water = 1000.; #[kg/m**(3)] - Density of water
m_water = 100.; #[g] - Mass of water
m_NaCl = 3.5; #[g] - Mass of NaCl
M_wt_water = 18.015; # Molecular weight of water
M_wt_NaCl = 58.5; # Molecular weight of NaCl
mol_water = m_water/M_wt_water; #[mol] - Moles of water
mol_NaCl = m_NaCl/M_wt_NaCl; #[mol] - Moles of NaCl
# Calculations
H_fus = -80.; #[cal/g] - Enthalpy change of fusion at 0 C
H_fus = H_fus*4.186*M_wt_water; #[J/mol]
#Mole fraction of the solute (NaCl) is given by
x_2 = mol_NaCl/(mol_NaCl+mol_water);
x_2_act = 2*x_2; # Actual mole fraction
x_1 = 1 - x_2_act;
#Molar volume of water can be calculated as
V_l_water = (1/den_water)*M_wt_water*10**(-3); #[m**(3)/mol]
#The osmotic pressure is given by
pi = -(R*T*math.log(x_1))/V_l_water; #[N/m**(2)]
pi = pi*10**(-5); #[bar]
# Results
print "The minimum pressure to desalinate sea water is %f bar"%(pi);
# Variables
R = 8.314; #[J/mol*K] - universal gas constant
T = 173.15; #[K] - Surrounding temperature
P = 60; #[bar]
P = P*10**(5); #[Pa]
#componenet 1 : CO2 (1)
#componenet 2 : H2 (2)
P_1_sat = 0.1392; #[bar] - Vapour pressre of pure solid CO2
P_1_sat = P_1_sat*10**(5); #[bar]
V_s_1 = 27.6; #[cm**(3)/mol] - Molar volume of solid CO2
V_s_1 = V_s_1*10**(-6); #[m**(3)/mol]
n_1 = 0.01; #[mol] - Initial number of moles of CO2
n_2 = 0.99; #[mol] - Initial number of moles of H2
#Let us determine the fugacity of solid CO2 (1) at 60 C and 173.15 K
# f_1 = f_1_sat*math.exp(V_s_1*(P-P_1_sat)/(R*T))
# Calculations
#Since vapour pressure of pure solid CO2 is very small, therefore
f_1_sat = P_1_sat;
f_1 = f_1_sat*math.exp(V_s_1*(P-P_1_sat)/(R*T));
#Since gas phase is ideal therefore
# y1*P = f_1
y1 = f_1/P;
#Number of moles of CO2 in vapour phase at equilibrium can be calculated as
#y1 = (n_1_eq/(n_1_eq + n_2)). Therefore
n_1_eq = n_2*y1/(1-y1);
#Therefore moles of CO2 precipitated is
n_ppt = n_1 - n_1_eq; #[mol]
# Results
print "The moles of CO2 precipitated is %f mol"%(n_ppt);
# Variables
R = 8.314; #[J/mol*K] - universal gas constant
T = 350.; #[K] - Surrounding temperature
#componenet 1 : organic solid (1)
#componenet 2 : CO2 (2)
P_1_sat = 133.3; #[N/m**(2)] - Vapour pressre of organic solid
V_s_1 = 200.; #[cm**(3)/mol] - Molar volume of organic solid
V_s_1 = V_s_1*10.**(-6); #[m**(3)/mol]
#/At 350 K, the values of the coefficients
B_11 = -500.; #[cm**(3)/mol]
B_22 = -85.; #[cm****(3)/mol]
B_12 = -430.; #[cm**(3)/mol]
# Calculations
#From phase equilibrium equation of component 1, we get
# y1*P*phi_1 = f_1
# f_1 = f_1_sat*math.exp(V_s_1*(P-P_1_sat)/(R*T))
#Since vapour pressure of organic solid is very small, therefore
f_1_sat = P_1_sat;
# Now let us determine the fugacity coefficient of organic solid in the vapour mixture.
# math.log(phi_1) = (P/(R*T))*(B_11 + y2**(2)*del_12)
del_12 = (2*B_12 - B_11 - B_22)*10**(-6); #[m**(3)/mol]
#It is given that the partial pressure of component 1 in the vapour mixture is 1333 N?m**(2)
# y1*P = 1333 N/m**(2) or, y1 = 1333/P
# y2 = 1- 1333/P
# math.log(phi_1) = (P/(R*T))*(B_11 + (1- 1333/P)**(2)*del_12)
#The phase equilibrium equation becomes
# y1*P*phi_1 = f_1_sat*math.exp(V_s_1*(P-P_1_sat)/(R*T))
#Taking math.log on both side we have
# math.log(y1*P) + math.log(phi_1) = math.log(f_1_sat) + (V_s_1*(P-P_1_sat)/(R*T))
# (V_s_1*(P-P_1_sat)/(R*T)) - math.log(phi_1) = math.log(1333/133.3) = math.log(10)
#substituting for math.log(phi_1) from previous into the above equation we get
# (V_s_1*(P-P_1_sat)/(R*T)) - (P/(R*T))*(B_11 + (1- 1333/P)**(2)*del_12) - math.log(10) = 0
# On simplification we get,
# 975*P**(2) - 6.7*10**(9)*P + 4.89*10**(8) = 0
# Solving the above qyadratic equation using shreedharcharya rule
P3 = (6.7*10**(9) + ((-6.7*10**(9))**(2)-4*975*4.98*10**(8))**(1./2))/(2*975); #[Pa]
P4 = (6.7*10**(9) - ((-6.7*10**(9))**(2)-4*975*4.98*10**(8))**(1./2))/(2*975); #[Pa]
# The second value is not possible, therefore pressure of the system is P3
P3 = P3*10**(-5); #[bar]
# Results
print " The total pressure of the system is %f bar"%(P3);