function [v,err,count] = Clapeyrons(R,T,a,b,P);
R = 0.082057;
T = 293;
Tc = 416.90;
Pc = 78.72918;
a = ((R^2)*(Tc^(5/2)))/(9*Pc*(2^(1/3)-1));
b = (R*Tc*(2^(1/3)-1))/(3*Pc);
P = 1;
if (~isscalar(P)) || (~isreal(P)) || P <= 0
error('Input argument P must be positive real scalar.')
end
Iteration_limit = 20;
Tolerance=10^-7;
A = (a*P)/((R^2)*(T^(5/2)));
B = (b*P)/(R*T);
v = R * T / P;
Z = (P*v)/(R*T);
for count = 1:Iteration_limit + 1
if count == Iteration_limit + 1
error('Iteration limit reached. Iteration did not converge.')
break
end
f = (Z^3) - (Z^2) + (A-B-(B^2))*Z - (A*B);
if abs(f)<Tolerance
break
end
f_dash = 3*Z^2 - 2*Z + (A - B - (B^2));
Z = Z - (f/f_dash);
v = Z*R*T/P;
end
err=abs(f);
end
I need to change P=1 so...
P = [1,1.5,2,2.5,3,5,10,15,25,50,100];
R = 0.082057;
T = 293;
Tc = 416.90;
Pc = 78.72918;
V_real =zeros(11,1);
V_ideal =zeros(11,1);
for count = 1:11
V_real(count) = Clapeyrons(P(count));
V_ideal(count) = (R*T)/P(count);
end
plot(P,V_real)
hold on
plot(P,V_ideal,'x')
title('Molar Volume vs Pressure for Cl2')
xlabel('Pressure P (atm)')
ylabel('Molar Volume v (L/mol)')
legend({'Redlich-Kwong Equation','Ideal Gas Law'},'Location','northeast')
...this graph for the Redlich Kwon Equation line will go through the same P values as in the Ideal Gas Law