why why why why why ..........​..........​..........​..........​..

3 views (last 30 days)
Nasir Qazi
Nasir Qazi on 27 Feb 2012
I have an code and when I try to run it, the status say busy , and it wont show me the answer to my solution , what exactly is the problem ?
[EDITDED, Jan Simon]: Copied from the comment section:
This program does not converge:
% Calculate Bubble Point and Dew Point SRK-EOS
% P = RT / V - b - alpha*a / V(V-b)
%----------------------------------------------
% liquid phase composition
y = [0.86 0.05 0.05 0.02 0.01 0.0095 0.0005]';
% vapor phase composition
x = [0.45 0.05 0.05 0.03 0.01 0.01 0.4]';
% critical pressure
Pc =[4.5947 4.8711 4.2472 3.6397 3.3688 3.1233 1.965]';
% critical temperature
Tc =[190.7389 305.5111 370.0333 425.3444 469.8889 512.7778 644.4444]';
% Accentric Factor
w = xlsread('myfile.xlsx','Sheet1','E2:E8');
%---------------------------------------------
L = 1;
R = 8.314e-3; % MPa m3 / mole.K
T = 200; % K
Tr = T./Tc; % K
P(L) = 5.0; %MPa
m = 0.480+1.574.*w -0.176.*w.^2;
al=1+ m.*(1-sqrt(Tr)).^2;
a = 0.42748.*(((R.*Tc).^2)./Pc);
b = 0.08664.*(R.*Tc./Pc);
q=0;
while q==0
% Z3 - Z2 + (A-B-B^2)Z-AB=0
% for gas phase
SA = 0;
RR = zeros(7,1);
for i = 1:7
for j = 1:7
SA = SA + y(i).*y(j).*sqrt(a(i).*a(j).*al(i).*al(j));
RR(i) = RR(i) + y(j).*sqrt(a(i).*a(j).*al(i).*al(j));
end
end
SB = sum(y.*b);
% For Liquid Phase
SAA = 0;
MM = zeros(7,1);
for i = 1:7
for j = 1:7
SAA = SAA + x(i).*x(j)*sqrt(a(i).*a(j)*al(i).*al(j));
MM(i) = MM(i) +x(j).*sqrt(a(i).*a(j).*al(i).*al(j));
end
end
SBB = sum(x.*b);
%for Gas Phase
A = SA.*P(L)./(R*T).^2;
B = SB.*P(L)./(R*T);
ZV=roots([1 -1 A-B-B^2 -A*B]);
Zv=max(ZV);
Zvc=isreal(Zv);
if Zvc == false
Zv=abs(Zv);
disp('error 1')
end
%for liquid Phase
AA = (SAA.*P(L))./(R*T).^2;
BB = (SBB.*P(L)./(R*T));
Z=roots([1 -1 AA-BB-BB^2 -AA*BB]);
Zl=min(Z);
Zvl=isreal(Zl);
if Zvl == false
Zl=abs(Zl);
disp('error 2')
end
phiv =exp((b.*(Zv-1)/SB) - log(Zv-B)-(A/B).*((2.*RR/SA)-
(b/SB)).*log(1 + B/Zv));
phil =exp((b.*(Zl-1)./SBB) - log(Zl-BB)-(AA./BB).*((2.*MM/SAA)-
(b./SBB)).*log(1 + BB./Zl));
K = phil./phiv;
RA1= sum(y./K)-1;
if abs(RA1) < 1e-5 && abs((y./K) - x) < 1e-5
break
else
if L==1
L = L+1;
P(L) = 1.01.*P(L-1);
RA1;
else
if L==1
P(L)=P(L-1)-(RA1./(RA1-RA2))*(P(L-1)-P(L-2));
RA2 = RA1;
end
%adjustment
x2 = (1-0.5)*(y./K) + 0.5*x;
x = x2;
end
end
end
disp(P(L));
  2 Comments
Jiro Doke
Jiro Doke on 27 Feb 2012
Can you please rephrase the subject line so that it is specific to your question? Please read here:
http://www.mathworks.com/matlabcentral/about/answers/
and read the "Tips" under "Getting Started" -> "Ask a Question"
Walter Roberson
Walter Roberson on 27 Feb 2012
Duplicate is at http://www.mathworks.com/matlabcentral/answers/30397-trial-and-error-help-is-not-giving-me-solution-to-that-still-waiting

Sign in to comment.

Answers (3)

Grzegorz Knor
Grzegorz Knor on 27 Feb 2012
It means, that program is still running. Run it in debug mode, and you will see step by step how it works.

Jan
Jan on 27 Feb 2012
You can insert a waitbar or some output in the code to be informed about the status of process. Either the computations take a lot of time or you have created an infinite loop by accident.
  4 Comments
Nasir Qazi
Nasir Qazi on 27 Feb 2012
% Calculate Bubble Point and Dew Point SRK-EOS
% P = RT / V - b - alpha*a / V(V-b)
%----------------------------------------------
% liquid phase composition
y = [0.86 0.05 0.05 0.02 0.01 0.0095 0.0005]';
% vapor phase composition
x = [0.45 0.05 0.05 0.03 0.01 0.01 0.4]';
% critical pressure
Pc =[4.5947 4.8711 4.2472 3.6397 3.3688 3.1233 1.965]';
% critical temperature
Tc =[190.7389 305.5111 370.0333 425.3444 469.8889 512.7778 644.4444]';
% Accentric Factor
w = xlsread('myfile.xlsx','Sheet1','E2:E8');
%---------------------------------------------
L = 1;
R = 8.314e-3; % MPa m3 / mole.K
T = 200; % K
Tr = T./Tc; % K
P(L) = 5.0; %MPa
m = 0.480+1.574.*w -0.176.*w.^2;
al=1+ m.*(1-sqrt(Tr)).^2;
a = 0.42748.*(((R.*Tc).^2)./Pc);
b = 0.08664.*(R.*Tc./Pc);
q=0;
while q==0
% Z3 - Z2 + (A-B-B^2)Z-AB=0
% for gas phase
SA = 0;
RR = zeros(7,1);
for i = 1:7
for j = 1:7
SA = SA + y(i).*y(j).*sqrt(a(i).*a(j).*al(i).*al(j));
RR(i) = RR(i) + y(j).*sqrt(a(i).*a(j).*al(i).*al(j));
end
end
SB = sum(y.*b);
% For Liquid Phase
SAA = 0;
MM = zeros(7,1);
for i = 1:7
for j = 1:7
SAA = SAA + x(i).*x(j)*sqrt(a(i).*a(j)*al(i).*al(j));
MM(i) = MM(i) +x(j).*sqrt(a(i).*a(j).*al(i).*al(j));
end
end
SBB = sum(x.*b);
%for Gas Phase
A = SA.*P(L)./(R*T).^2;
B = SB.*P(L)./(R*T);
ZV=roots([1 -1 A-B-B^2 -A*B]);
Zv=max(ZV);
Zvc=isreal(Zv);
if Zvc == false
Zv=abs(Zv);
disp('error 1')
end
%for liquid Phase
AA = (SAA.*P(L))./(R*T).^2;
BB = (SBB.*P(L)./(R*T));
Z=roots([1 -1 AA-BB-BB^2 -AA*BB]);
Zl=min(Z);
Zvl=isreal(Zl);
if Zvl == false
Zl=abs(Zl);
disp('error 2')
end
phiv =exp((b.*(Zv-1)/SB) - log(Zv-B)-(A/B).*((2.*RR/SA)-
(b/SB)).*log(1 + B/Zv));
phil =exp((b.*(Zl-1)./SBB) - log(Zl-BB)-(AA./BB).*((2.*MM/SAA)-
(b./SBB)).*log(1 + BB./Zl));
K = phil./phiv;
RA1= sum(y./K)-1;
if abs(RA1) < 1e-5 && abs((y./K) - x) < 1e-5
break
else
if L==1
L = L+1;
P(L) = 1.01.*P(L-1);
RA1;
else
if L==1
P(L)=P(L-1)-(RA1./(RA1-RA2))*(P(L-1)-P(L-2));
RA2 = RA1;
end
%adjustment
x2 = (1-0.5)*(y./K) + 0.5*x;
x = x2;
end
end
end
disp(P(L));
% I am trying to using trial and error to get the results of P(L) , but my solution is not converging and keep going on. I need to know what can I do to resolve the issue .
Help Please
Jan
Jan on 27 Feb 2012
Well, this is an important detail: As far as I understand your problem is the lack of convergence. This is another story.

Sign in to comment.


Sean de Wolski
Sean de Wolski on 27 Feb 2012
The bald and not excessively bald and not excessively smart hamster obeyed a terrified and not excessively terrified hamster.
  1 Comment
Jan
Jan on 27 Feb 2012
@Nasir: If you are wondering, what the deeper meaning of this answer could be, type "why" in the command line - repeatedly.

Sign in to comment.

Categories

Find more on Programming in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!