سلام
من این برنامه رو نوشتم:
clcclear all
syms Z Zv Zl
R=8.314;
n=0;
n=n+1;
Tc=304.2;
Pc=73.8;
vc=93.9;
w=0.225;
ac=0.457235*(R*Tc)^2/Pc ;
b=0.077796*R*Tc/Pc;
k=0.37464+1.54226*w-0.26993*w^2;
for tn=-50:20:200;
T=tn+273.15;
Tr=T/Tc;
p=10^((7.5322-((835.06)*(tn+268.223)^(-1)))*(1.013/760));
p
x=(1+k*(1-sqrt(Tr)))^2;
a=ac*x;
A=a*p*(R*T).^(-2);
B=b*p*(R*T)^(-1);
Z=[1 -(1-B) (A-2*B-3*B^2) (A*B-B^2-B^3)];
z=roots(Z)
if z>0
Zl=min(zz)
Zv=max(zz)
end
phiv=(Zv-1)-log(Zv-B)-A*log((Zv+B(1+sqrt(2)))*(Zv+B(1-sqrt(2))))*(sqrt(8)*B)^(-1);
phil=(Zl-1)-log(Zl-B)-A*log((Zl+B(1+sqrt(2)))*(Zl+B(1-sqrt(2))))*(sqrt(8)*B)^(-1);
P(n+1) =Pn * (phil/phiv);
error=(P(n+1)-Pn);
if (error < 10^-4); % this statement has to be verify if not
print Pn tn
else % |
P(n+1) =Pn * (phil/phiv); % V
P(n+1) =Pn; %( this pressure should update the previous assume pressure )
end
end
باید به این شکل باشه که در بخشی که نوشتم z=roots(z) در نتیجه به من ۳ تا Z بده. دوتا مثبت و یکی منفی.
باید عدد اول که مثبت هست رو به Zv بده٬ عدد دوم مثبت که کوچکتر از Zv هست رو به Zl بده. اما اتفاقی که میوفته اینه که عدد اول بزرگ رو به درستی نشون میده اما به جای عدد مثبت دوم در Zl عدد منفی رو جاش نشون میده.
در حالی که من نوشتم: if z>0 اما باز هم اشتباه رخ میده.
ممنون میشم بررسی کنید ایراد منو بگید.
تشکر