
Calculating Mixture Critical Temperature And Pressure Using Equation O
#1
Posted 07 August 2018  07:58 PM
#2
Posted 08 August 2018  01:47 PM
You state that you "trying to calculate critical Temperature and Pressure of binary and ternary mixtures using an equation of state". However, you have a problem: you are "way off experimental data".
What is your question or query? It can't be that your calculations are wrong or not the correct ones because you haven't furnished your calculations. If your question is: "why are my calculations not giving correct results?", then you must furnish your calculations  assuming that you have already had them proofed and checked for accuracy and applicability. Our Forum members can't do much without a copy of your work in order to determine if it is the correct method and if it is applied correctly.
#3
Posted 09 August 2018  07:02 AM
http://www.scielo.br...322000000400040
https://www3.nd.edu/...kst/bas2000.pdf
Hi ,
you may consider the links above .
Breizh
Edited by breizh, 09 August 2018  07:05 AM.
#4
Posted 09 August 2018  08:13 PM
Hi,
Thank you so much for the response.
Breizh: The links you have given use the Gibb's method of finding the critical T and P. I am not following that method.
The goal is to try and use the inflection point property and find the critical points. Which states dP/dv = 0 and d^{2}P/dv^{2 }= 0 [these are partial derivatives of pressure with respect to molar volume]
I am giving my method here: Will be grateful to anyone who can point put a bug and show me the correct way
1. I start with Peng Robinson EOS:
P = { RT/(v }  { a / (v^{2} + 2vb  b^{2}) }
2. I find expressions for dP/dv = 0 and d^{2}P / dv^{2}^{}= 0.
3. I equate them to end up with a cubic equation: v^{3 } (3b)v^{2 } (3b^{2})v  (3b^{3}) = 0
This has 3 roots but one is real and rest two are complex conjugates.
The real root gives me the critical volume
4. I then put back this value of 'v' in the equation of dP/dv = 0. T is the only unknown in here now. I solve it to get T.
5. Now I put this T back in equation of state to get pressure.
I am giving my code below.
% Getting Pseudo Critical properties % Using PengRobinson Equation of State % Binary mixture: nmethane + propane % Author: Suman Chakraborty clear, clc; %% User Inputs % Universal gas constant % R = 8.31446 J/mol K; 83.144 cm^3 bar/mol K R = 83.144; % mole fraction liquid phase x_C1 = 0:0.001:1; x_C3 = 1  x_C1; % size of mole fraction array N = size(x_C1); N = N(1,2); %% methane Data % Critical Temperature (K) Tc_C1 = 190.6; % Critical Pressure (bar) Pc_C1 = 46.1; % Molecular weight (gm/mole) MW_C1 = 16.0425; % Acentric Factor (w) w_C1 = 0.011; %0.010 %% propane Data % Critical Temperature (K) Tc_C3 = 369.9; % Critical Pressure (bar) Pc_C3 = 42.5; % Molecular weight (gm/mole) MW_C3 = 44.0956; % Acentric Factor (w) w_C3 = 0.153; %0.152 %% methane Pure component calculation % a(Tc) [cm^6 * bar/mol^2] a_Tc_C1 = (0.45724 * R * R * Tc_C1 * Tc_C1)/Pc_C1; % b(Tc) [cm^3 / mole] b_C1 = (0.0778 * R * Tc_C1) / Pc_C1; % kappa kappa_C1 = 0.37464 + (1.54226*w_C1)  (0.26992*w_C1*w_C1); %% propane Pure component calculation % a(Tc) [cm^6 * bar/mol^2] a_Tc_C3 = (0.45724 * R * R * Tc_C3 * Tc_C3)/Pc_C3; % b(Tc) [cm^3 / mole] b_C3 = (0.0778 * R * Tc_C3) / Pc_C3; % kappa kappa_C3 = 0.37464 + (1.54226*w_C3)  (0.26992*w_C3*w_C3); %% Main Loop Starts for i=1:N % Application of mixing rules % Repulsion parameter of mixture  [process 1] %bm = (x_C1(i)*b_C1) + (x_C3(i)*b_C3); % Repulsion parameter of mixture  [process 2] term1 = (x_C1(i)*(b_C1^(1/3))) + (x_C3(i)*(b_C3^(1/3))); bm = term1^3; % Solving cubic equation a3 = 1; % coeff of v^3 a2 = (3*bm); % coeff of v^2 a1 = (3*bm*bm); % coeff of v a0 = 3*(bm^3); % constant coeff = [a3 a2 a1 a0]; % Z is compressibility factor soln_v = roots(coeff); % Roots real and >0 v_real = []; for j=1:3 %>B to get real fugacity, since fhi has a term log(ZB) if isreal(soln_v(j)) && soln_v(j) > 0 temp = real(soln_v(j)); v_real = [v_real temp]; end end index = length(v_real); if index >= 1 % Vapor Phase v_final = max(v_real); end % Guess Temperature [K] T_guess = 2000; % Convergence variable err = 1.0; eps = 0.01; while true % Reduced Temperature Tr_C1 = T_guess / Tc_C1; % Reduced Temperature Tr_C3 = T_guess / Tc_C3; % a(T) of C1 [cm^6 * bar/mol^2] term1 = 1  (Tr_C1^0.5); term2 = 1 + (kappa_C1*term1); a_T_C1 = a_Tc_C1 * term2 * term2; % a(T) of C3 [cm^6 * bar/mol^2] term1 = 1  (Tr_C3^0.5); term2 = 1 + (kappa_C3*term1); a_T_C3 = a_Tc_C3 * term2 * term2; % Attraction parameter of mixture % Assuming K_ij=0 % 1> C1 2> C3 a11 = (x_C1(i)*x_C1(i)) * a_T_C1; a12 = (x_C1(i)*x_C3(i)) * ((a_T_C1*a_T_C3)^0.5); a21 = (x_C3(i)*x_C1(i)) * ((a_T_C3*a_T_C1)^0.5); a22 = (x_C3(i)*x_C3(i)) * a_T_C3; am = a11 + a12 + a21 + a22; %LHS of condition term1 = (v_finalbm) ^ 2; LHS = (R*T_guess)/term1; %RHS of condition term1 = 2*am*(v_final+bm); term2 = (v_final*v_final) + (2*v_final*bm)  (bm*bm); RHS = term1/(term2*term2); err = (LHSRHS)/LHS; err = abs(err); if err<eps break; end T_guess = T_guess  0.5; end Tc(i) = T_guess; %calculating pressure term1 = (R*Tc(i))/(v_finalbm); term2 = (v_final*v_final) + (2*v_final*bm)  (bm*bm); Pc(i) = term1  (am/term2); end
#5
Posted 12 August 2018  02:26 AM
http://pillars.che.p...ns_of_State.pdf
Hi.
Consider the document above in addition to my previous email.
Breizh
Similar Topics
Liquid PressureStarted by Guest_dhindsa_mexican_* , 17 Aug 2018 



Keeping A Constant PressureStarted by Guest_Kangari_* , 06 Aug 2018 



Max Differential Pressure Across An On/off Valve Prior To Open ItStarted by Guest_Roark_* , 14 Aug 2018 



Vessel Vacuum Due To Pressure Control FailureStarted by Guest_Butterfly_* , 12 Aug 2018 



Understanding Liquid PressureStarted by Guest_dhindsa_mexican_* , 08 Aug 2018 

