Iteration using for loop and if statement

3 次查看（过去 30 天）
GAGANDEEP KAUR
GAGANDEEP KAUR 2021-1-18
评论： GAGANDEEP KAUR 2021-1-27
I am writing a code where I need to calculate P( as P_Initial) with all the variables calculated in other functions, then using this P_initial calculte PhiN, assign PhiN to Phi and once again calculate P(as P_total) with new value of Phi. Now if condition needs to be applied if true then print output ,if false recalculte PhiN using recent P_total and then check for condition again. This is what I need to do but while writing it in loop and to assign new value iteratively I am getting some issues. I am tring as follows. My concern is from "for loop o=1:N". Please i could get some suggestion to write it properly
%%To calculate system pressure and vapor phase composition
%yi=vapor composition
%P=Pressure(From Antoine's Equation)
%e=tolerance
e=0.4
[f_sat]=FugacityS(T,P,Tc,Pc)
[Gamma_t]=Gamma_sum(taua,x1,x2,x3,rho,zc,za,xc,xa,xm,T)
[Phi]=Redlich_kwong(T,Tc,Pc)
N=length(T);
M=length(Pc);
xm=[0.5096 0.5092 0.5087
0.0963 0.0964 0.0965
0.3941 0.3944 0.3948];
a0=zeros(1,M);aHI=zeros(1,N);aH2O=zeros(1,N);aI2=zeros(1,N);a=zeros(N,1);
for i=1:M
a0(i)=(0.42748*(((R^2)*(Tc(i)^2.5))/Pc(i)));
end
for j=1:N
for i=1
aHI(j)=(((a0(i)*a0(i+1))^0.5)*xm(i,j)*xm(i+1,j))+(((a0(i)*a0(i+2))^0.5)*xm(i,j)*xm(i+2,j)); %applying mixing rules
aH2O(j)=(((a0(i+1)*a0(i))^0.5)*xm(i+1,j)*xm(i,j))+(((a0(i+1)*a0(i+2))^0.5)*xm(i+1,j)*xm(i+2,j));
aI2(j)=(((a0(i+2)*a0(i))^0.5)*xm(i+2,j)*xm(i,j))+(((a0(i+2)*a0(i+1))^0.5)*xm(i+2,j)*xm(i+1,j));
end
end
format long g
for j=1:N
a(j)=(aHI(j)+aH2O(j)+aI2(j));
a(j)=(a(j))';
end
b0=zeros(N,M);b=zeros(1,N);A0=zeros(N,M);A=zeros(1,N);B0=zeros(N,M);
B=zeros(1,N);
for j=1:N
for i=1:M
b0(i,j)=(0.08664*((R*Tc(i))/Pc(i)))*(xm(i,j));
A0(i,j)=(((0.4278*Tc(i)^2.5)/(Pc(i)*T(j)^2.5))^(1/2))*(xm(i,j));
B0(i,j)=((0.08664*Tc(i))/(Pc(i)*T(j))*(xm(i,j)));
end
b(j)=sum(b0(:,j));
A(j)=sum(A0(:,j));
B(j)=sum(B0(:,j));
end
Asq=((A).^2);
P_calI=zeros(N,M); y=zeros(N,M); P_Initial=zeros(N,1);y_total=zeros(N,1);
for o=1:N
for p=1:M
P_calI(p,o)= (xm(p,o)*Gamma_t(p,o)*f_sat(p,o))/Phi(p,o)
end
P_Initial(p,o)=sum(P_calI(:,o))
Z=zeros(1,M);h=zeros(1,M);PhiN=zeros(1,M);
h(p,o)=b(o)./Vm(p,o);
Z(p,o)=(1./1-h(p,o))-((Asq(o)/B(o)).*(h(p,o)./(1+h(p,o))));
PhiN(p,o)=exp(Z(p,o)-1-log(Z(p,o)-(B(o).*P_cal(p,o)))-(Asq(o)/B(o)).*(log(1+((B(o).*P_calI(p,o))./Z(p,o)))))
Phi=PhiN; P_total=zeros(N,M);P_cal=zeros(N,M)
P_cal(p,o)= (xm(p,o)*Gamma_t(p,o)*f_sat(p,o))/Phi(p,o)
P_total(p,o)=sum(P_cal(:,o))
if ((ptotal-p_Initial/p_Initial)<e)
y(p,o)= (xm(p,o)*Gamma_t(p,o)*f_sat(p,o))/(Phi(p,o)*P_Initial(p,o))
disp(p_total)
disp(y)
else
Z=zeros(1,M);h=zeros(1,M);PhiN=zeros(1,M);
h(p,o)=b(o)./Vm(p,o);
Z(p,o)=(1./1-h(p,o))-((Asq(o)/B(o)).*(h(p,o)./(1+h(p,o))));
PhiN(p,o)=exp(Z(p,o)-1-log(Z(p,o)-(B(o).*P_cal(p,o)))-(Asq(o)/B(o)).*(log(1+((B(o).*P_cal(p,o))./Z(p,o)))))
Phi=PhiN; P_total1=zeros(N,M)
P_cal(p,o)= (xm(p,o)*Gamma_t(p,o)*f_sat(p,o))/Phi(p,o)
P_total1(p,o)=sum(P_cal(:,o))
end
end
  3 个评论
显示 1更早的评论
GAGANDEEP KAUR
GAGANDEEP KAUR 2021-1-19
I understand your concern and I was also considering it quite odd but all other variables are calculated from 6-7 separate function files and needs to be called here which will be quite complicated to share here.
Still if any general commands can be suggested for writing, as we generally perform these steps in numerical methods,like to get initial guess( here P_initial) then from that PhiN and for second value of P(P_totoal here) recently calculated Phi(PhiN) then check error if lesser than accepted end by displying calculated variables if not calculate Phi again using recent P. Mthematically I have idea but in writing I am not able .
GAGANDEEP KAUR
GAGANDEEP KAUR 2021-1-27
Here, I am sharing quite edited version of aforementioned code with values of variables needed. Actually I am getting problem in re-entering the for loop need to enter p loop again where ym is to calculted again with new parameters calculated until condition is not satisfied.After the condition get satisfied, then it should enter O loop(of ym) and same calculations to be done.
%%To calculate system pressure and vapor phase composition
e=0.3
f_sat= [ 32.2587857596876 34.0501961493468 35.8746854026426
2.00528670528962 2.33558884889119 2.70783015514365
0.157445079883318 0.187083254240374 0.221158510598425];
Gamma_t=[3.19702443124609 2.54070700562388 2.08314683871606;
3.31850686110491 2.66546741032482 2.20838701999553;
3.44415347633114 2.79709939502664 2.34012176373733];
Vm= [ 68.0045214121783 70.0736822139801 72.5465660407398
18.1158800129555 18.2324752448415 18.3515283669796
51.5543766322469 51.7464896478157 51.9409677300663];
Pc= [831,220.58,116.54];
Tc=[394.15 399.15 404.15];
xm=[0.09630,0.0964,0.0965;
0.5096,0.5092,0.5087;
0.3941,0.3944,0.3948];
N=length(xm);
M=length(Pc);
Phi(1:3,1:3)=1;
P_calI=zeros(N,M);ym=zeros(N,M);P_0=zeros(N,1);
for o=1:N
for p=1:M
P_calI(p,o)= (xm(p,o)*Gamma_t(p,o)*f_sat(p,o))/Phi(p,o)
end
P_0(o)=sum(P_calI(:,o))
end
for o=1:N
for p=1:M
ym(p,o)= (xm(p,o)*Gamma_t(p,o)*f_sat(p,o))/(Phi(p,o)*P_0(o))
end
a0=zeros(1,M);aHI=zeros(1,N);aH2O=zeros(1,N);aI2=zeros(1,N);a=zeros(N,1);
R=0.08314;
for i=1:M
a0(i)=(0.42748*(((R^2)*(Tc(i)^2.5))/Pc(i)));
end
for i=1
aHI(o)=(((a0(i)*a0(i+1))^0.5)*ym(i,o)*ym(i+1,o))+(((a0(i)*a0(i+2))^0.5)*ym(i,o)*ym(i+2,o)); %applying mixing rules
aH2O(o)=(((a0(i+1)*a0(i))^0.5)*ym(i+1,o)*ym(i,o))+(((a0(i+1)*a0(i+2))^0.5)*ym(i+1,o)*ym(i+2,o));
aI2(o)=(((a0(i+2)*a0(i))^0.5)*ym(i+2,o)*ym(i,o))+(((a0(i+2)*a0(i+1))^0.5)*ym(i+2,o)*ym(i+1,o));
end
format long g
a(o)=(aHI(o)+aH2O(o)+aI2(o));
a(j)=(a(j))';
b0=zeros(N,M);b=zeros(1,N);A0=zeros(N,M);A=zeros(1,N);B0=zeros(N,M);
B=zeros(1,N);
for i=1:M
b0(i,o)=(0.08664*((R*Tc(i))/Pc(i)))*(ym(i,o));
A0(i,o)=(((0.4278*Tc(i)^2.5)/(Pc(i)*T(o)^2.5))^(1/2))*(ym(i,o));
B0(i,o)=((0.08664*Tc(i))/(Pc(i)*T(o))*(ym(i,o)));
end
b(o)=sum(b0(:,o));
A(o)=sum(A0(:,o));
B(o)=sum(B0(:,o));
Asq=((A).^2);
Z=zeros(N,M);h=zeros(N,M);Phi=zeros(N,M); P_1=zeros(N,1);P_cal=zeros(N,M)
for p=1:M
h(p,o)=b(o)./Vm(p,o)
Z(p,o)=(1./1-h(p,o))-((Asq(o)/B(o)).*(h(p,o)./(1+h(p,o))))
Phi(p,o)=exp(Z(p,o)-1-log(Z(p,o)-(B(o).*P_calI(p,o)))-(Asq(o)/B(o)).*(log(1+((B(o).*P_calI(p,o))./Z(p,o)))))
P_cal(p,o)= (xm(p,o)*Gamma_t(p,o)*f_sat(p,o))/Phi(p,o)
end
P_1(o)=sum(P_cal(:,o))
if ((P_1(o)-P_0(o)/P_0(o))<e)
break
end
P_0=P_1
P_calI=P_cal
end

请先登录，再进行评论。

请先登录，再回答此问题。

回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Loops and Conditional Statements 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by