function PltFigOfEquil
%PltFigOfEquil.m n     
%用圖解法求二元物連續(xù)精餾所需的理論板數(shù)
%Authou:LI ZiGuang
%Date:2005/3/31

clear all;clc
global y               
global Locat           %Locat:進(jìn)料板的位置
global BanOfRect       %BanOfRect:精餾段板數(shù)
global BanOfRas        %BanOfRas:提餾段板數(shù)
global BanOfTotal      %BanOfTotal:總的塔板數(shù)


%已知的參數(shù)
xd = 0.95                %塔頂餾出物中輕組分的摩爾分?jǐn)?shù)
xw = 0.05              %塔底餾出物中輕組分的摩爾分?jǐn)?shù)
xf = 0.55                %進(jìn)料中輕組分的摩爾分?jǐn)?shù)
R = 1.25          %全塔回流比
q = 1               %進(jìn)料的熱狀況參數(shù)

%用equilib函數(shù)（相平衡方程）求出一系列汽液平衡數(shù)據(jù)
for i = 1:11
    y = 0.1*(i-1)
    ye(i) = 0.1*(i-1)
    xe(i) = fzero(@equilib,0.5)
end

%求出精餾段操作線與進(jìn)料線的交點(diǎn)
yi=(xf+xd*q/R)/(1+q/R)
xi=(-(q-1)*(1-R/(R+1))*xd-xf)/((q-1)*R/(R+1)-q)

%-------------------------------------------------------------------
%求出所有折線與相平衡線、精餾操作線、提餾操作線的交點(diǎn)
i = 1
xp(1) = xd
yp(1) = xd
y = xd
slope = (yi-xw)/(xi-xw)
constant = xw-slope*xw

while xp(i)>xi
    xp(i+1) = fzero(@equilib,0.5)
    yp(i+1) = (R*xp(i+1))/(1+R)+xd/(1+R)
    y = yp(i+1)
    i = i+1
end
Locat = i-1
BanOfRect = i-1

yp(i) = slope*xp(i)+constant
y = yp(i)
while xp(i)>xw
    xp(i+1) = fzero(@equilib,0.5)
    yp(i+1) = slope*xp(i+1)+constant
    y = yp(i+1)
    i = i+1
end
BanOfRas = i-2-BanOfRect
BanOfTotal = i-2

%-----------------------------------------------------------------------

%---------------------------------------------------------------——————
%繪圖
figure(1); hold on; axis([0 1 0 1]);axis('square')
xlabel('x'); ylabel('y'); title('正戊烷-正己烷二元精餾塔的理論塔板數(shù)')
set(line([0 1],[0 1]),'Color',[0 1 0])                    %畫(huà)對(duì)角線
plot(xe,ye,'b')                                           %畫(huà)相平衡線
set(line([xi xd],[yi xd]),'Color',[1 0 1])                %畫(huà)精餾操作線
set(line([xf,xi],[xf,yi]),'Color',[1 0 1])                %畫(huà)進(jìn)料線
set(line([xw xi],[xw yi]),'Color',[1 0 1])                %畫(huà)提餾操作線
for ii = 1:1:(i-1)                                                    %%%%%%%%%%%%%%%
    set(line([xp(ii) xp(ii+1)],[yp(ii) yp(ii)]),'Color',[0 0 1])      %    繪制折線  %  
    if xp(ii+1)<xw
       break  
    end  
    set(line([xp(ii+1) xp(ii+1)],[yp(ii) yp(ii+1)]),'Color',[0 0 1])  % %%%%%%%%%%%%%
 end  
%------------------------------------------------------------------------

%輸出結(jié)果
disp('計(jì)算結(jié)果：')
fprintf('進(jìn)料口的位置:%s%d%s\n','在第',Locat,'塊塔板處')
fprintf('精餾段所需塔板數(shù):%d%s\n',BanOfRect,'(塊)')
fprintf('提餾段所需塔板數(shù):%d%s\n',BanOfRas,'(塊)')
fprintf('所需塔板總數(shù):%d%s\n',BanOfTotal,'(塊)')