Help with lsqnonline VLE (vapor liquid equilibrium) data.
2 vues (au cours des 30 derniers jours)
Afficher commentaires plus anciens
I am a beginner in the use of MATlLAB, and I need to find parameters A and B of the Redlich Kister equation and the new composition in equilibrium. My objective function is (Pe-Pc) ^ 2. Pe: Experimental pressure, Pc: Calculated Pressure.
my code is as follows
clc
T=input('Temperature=');%Temperatura [K]
R=8314; %[cm3*kPa/mol*K]
filename = 'DataVLE.xlsx';
sheet = 1;
PeRange = 'B2:B7';
Ps1Range= 'F2:F7';
m1Range = 'C2:C7';
m2Range = 'D2:D7';
rho1Range= 'G2:G7';
rho2Range= 'I2:I7';
B11Range= 'L2:L7';
%Cargar los datos experimentales de P, masa1 y masa 2
Pe=xlsread(filename,sheet,PeRange);
Ps1=xlsread(filename,sheet,Ps1Range);
m1=xlsread(filename,sheet,m1Range); %[gramos]
m2=xlsread(filename,sheet,m2Range); %[gramos]
M1=17.031; %masa molar del componente 1
M2=108.11; %maaa molar del componente 2
rho1=xlsread(filename,sheet,rho1Range); %densidad del líquido iónico [g/cm3]
rho2=xlsread(filename,sheet,rho2Range); %densidad del líquido iónico [g/cm3]
B11=xlsread(filename,sheet,B11Range); %segundo coeficiente del virial [cm3/mol]
VT=10.73566178; %Volumen total del equilibrio [cm3]
%Calculo de las moles y composición molar
n1=(m1./M1);
n2=(m2./M2);
z1=n1./(n1+n2);
z2=1-z1;
%suposiciones
y1=1;
rhom1=rho1/M1;
rhom2=rho2/M2; %[mol/cm3]
vL=n2.*rhom2; %[cm3/mol]
VL=m2./rho2; %[cm3]
%Volumen del vapor
VV=VT-VL;
%Parametros de Redlich - Kister (A' B')
A=1;B=1;
%Valores semilla
x1=z1;Pc=Ps1;A=1;B=1;
x0=[x1 Pc A B];
%función nonlinear solver
x = lsqnonlin(optim(x1,A,B,R,T,B11,VV,Ps1,vL,n1,n2,Pe),x0);
the first function is
function Pc= BARKER(x1,A,B,R,T,B11,VV,Ps1,vL,n1,n2,Pe)
gamma1=exp((A+3.*B).*(1-x1).^2-4.*B.*(1-x1).^3);
gamma2=exp((A-3.*B).*x1.^2+4.*B.*x1.^3);
Ps1c = Ps1.*exp((vL.^2-B11).*(Pe-Ps1)./(R*T));
Pc=gamma1.*x1.*Ps1c;
n1V=(Pc.*VV)/(R*T+B11.*Pc);
n1L=n1-n1V;
x1c=n1L./(n1L+n2);
end
The optimization function es
function min = optim(x1,A,B,R,T,B11,VV,Ps1,vL,n1,n2,Pe)
Pc = BARKER(x1,A,B,R,T,B11,VV,Ps1,vL,n1,n2,Pe);
min = sum(Pe - Pc).^2;
Voir également
Catégories
En savoir plus sur Particle & Nuclear Physics dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!