MATLAB: How to write a code for iteration

chemical engineeringchemistrythermodynamics

Hello. I'm new using Matlab. I have to programming the below diagram in Matlab. I think I'm doing well, so far. But I've got a doubt, that is how to write a code for iteration (the step in the third block)? I guess I have to use a 'while' but I didn't know how. Hope you can help me. I let what I've done so far, in case you need it in order to answer my question. I'm writing this code with the purpose of plotting a x-y diagram and a T-x-y diagram. Thanks.
    clear all; clc; close all
% Bubble temperature for P = 101.33 kPa
P = 101.33 %kPa
X = linspace (0,1,100)
phi = 1.0
x1 = X(1);        x2 = X(2);
% Wilson
    %PARAMETERS - System dependent
    a12 = 437.98; %cal/mol

    a21 = 1238.00; %cal/mol
    R = 1.9872; %cal/mol?K 
    V1 = 76.92; %cm3/mol 

    V2 = 18.07; %cm3/mol 
    comp = 2; 
    % Constants for Antoine equation for 2-propane(1), water(2)
    A1 = 16.6780;   B1 = 3640.20;   C1 = 53.54;
    A2 = 16.2887;   B2 = 3816.44;   C2 = 46.13;
    % Saturate temperature
    T1sat = (B1/(A1-(log(P)))-C1);
    T2sat = (B2/(A2-(log(P)))-C2);
    %Temperature
    T1 = sum(X*T1sat)
    T2 = sum(X*T2sat)
    % Activities
    A12 = (V2/V1)*exp(-a12/(R*T1)); 
    A21 = (V1/V2)*exp(-a21/(R*T2));
    % Vapour pressure
    P1sat = exp((A1-B1)/((T1)+C1));
    P2sat = exp((A2-B2)/((T2)+C2));
    for i=1:size(X,2)-1;
    x1 = X(i);        x2 = X(i+1);
    % Calc for Gammas
    ft1 = -log(x1 + 1-*A12);
    st1 = x2*((A12/(x1 + x2*A12))-(A21/(x2 + x1*A21)));
    gamma_1 = exp(ft1 + st1);
    ft2 = -log(x2 + x1*A21);
    st2 = x1*((A12/(x1 + x2*A12))-(A21/(x2 + x1*A21)));
    gamma_2 = exp(ft2 - st2);
    % j = Water
    % Vapour pressure (Pj), kPa

       Pjsat (i)= P/sum(((x2*gamma_2(i))/phi)*(P2sat/P2sat))
    % Temperature, K 
       Tj(i) = (A2/(B2 -(log(Pjsat (i))))-C2);
    % Vapour pressure, Pi_sat
       Pi_sat = exp((A2-B2)/((Tj(i))+C2));
    % Total pressure, P = 101.33 kPa
       Ptotal = P1sat*gamma_1(i)*x1 + P2sat*gamma_2(i)*x2;
       error = abs((101.33-Ptotal)/101.33);   
       counter = 0;  
       % Calc for y
          y (i) = ((x2*gamma_2(i)*Pi_sat)/phi*Ptotal);
       % Vapour pressure (Pj), kPa
          Pj_sat = P/sum(((x2*gamma_2(i))/phi)*(Pjsat/Pjsat))
       % Temperature, K
          Tj=(A2/(B2 -(log(Pj_sat)))-C2);
       % Iteration 
          if i>1
          while abs(Tj(i)-Tj(i-1))<1e-10
          break
          end
          else 
     end;

Best Answer

After the first box, you write
deltaT = epsilon + 1; % so that the while loop is entered
while deltaT >= epsilon
    % put code of the second box here
    deltaT = % add equation for new updated value of deltaT
 end
Related Question