MATLAB图解精馏塔理论塔板数程序代码.docx

MATLAB图解精馏塔理论塔板数程序代码.docx

《MATLAB图解精馏塔理论塔板数程序代码.docx》由会员分享，可在线阅读，更多相关《MATLAB图解精馏塔理论塔板数程序代码.docx（6页珍藏版）》请在冰点文库上搜索。

MATLAB图解精馏塔理论塔板数

MATLAB图解精馏塔理论塔板数程序代码

functiondistillation%文件名“distillation”可以更改

%输入计算所需参数

q=1;%输入进料热状况参数

R=1.5;%输入回流比

xD=0.95;%输入塔顶轻组分摩尔分数

xW=0.04;%输入塔底轻组分摩尔分数

xF=0.52;%输入进料轻组分摩尔分数

%以下输入相平衡数据

x0=[0

0.0196078

0.0392156

0.0588235

0.0784313

0.0980392

0.1176471

0.1372549

0.1568627

0.1764706

0.1960784

0.2156863

0.2352941

0.254902

0.2745098

0.2941176

0.3137255

0.3333333

0.3529412

0.372549

0.3921569

0.4117647

0.4313725

0.4509804

0.4705882

0.4901961

0.5098039

0.5294118

0.5490196

0.5686275

0.5882353

0.6078431

0.627451

0.6470588

0.6666667

0.6862745

0.7058824

0.7254902

0.745098

0.7647059

0.7843137

0.8039216

0.8235294

0.8431373

0.8627451

0.8823529

0.9019608

0.9215686

0.9411765

0.9607843

0.9803922

1];

y0=[0

0.0437029

0.0856224

0.1258286

0.1643911

0.2013788

0.2368595

0.2708994

0.303563

0.3349129

0.3650094

0.3939109

0.4216732

0.4483501

0.4739928

0.4986506

0.5223702

0.5451963

0.5671715

0.5883362

0.6087289

0.6283862

0.6473428

0.6656317

0.6832842

0.70033

0.7167974

0.7327131

0.7481026

0.76299

0.7773982

0.791349

0.8048631

0.8179601

0.8306587

0.8429768

0.8549313

0.8665382

0.877813

0.8887702

0.899424

0.9097874

0.9198734

0.9296939

0.9392607

0.9485847

0.9576768

0.966547

0.9752052

0.9836608

0.9919228

1];

Yr=@（x）R/（R+1）.*x+xD/（R+1）;%精馏段操作线

fun=@（x）（q-1）*（R/（R+1）.*x+xD/（R+1））-（q*（x-xF）+（q-1）*xF）;

xQ=fzero（fun,0.5）;%求操作点

yQ=Yr（xQ）;

xOP=[xW,xQ,xD];

yOP=[xW,yQ,xD];

yfit=linspace（0,1,1001）;

xfit=interp1（y0,x0,yfit,'pchip'）;

%%绘制图形

holdon

boxon

plot（[01],[01],'k'）;

xlabel（'x'）

ylabel（'y'）

plot（x0,y0,'r'）

plot（xfit,yfit,'r-'）

plot（xF,xF,'b*'）

plot（xQ,yQ,'bo'）

plot（xOP,yOP,'b-'）

k=1;

yn

（1）=xD;

xn

（1）=interp1（y0,x0,yn

（1）,'pchip'）;

plot（[xD,xn

（1）],[yn

（1）,yn

（1）],'b-'）

text（xn

（1）,yn

（1）,num2str

（1）,...

'HorizontalAlignment','center','VerticalAlignment','bottom'）

whilexn（k）>xW

yn（k+1）=interp1（xOP,yOP,xn（k））;

k=k+1;

xn（k）=interp1（y0,x0,yn（k）,'pchip'）;

plot（[xn（k-1）,xn（k-1）],[yn（k-1）,yn（k）],'b-'）

plot（[xn（k-1）,xn（k）],[yn（k）,yn（k）],'b-'）

text（xn（k）,yn（k）,num2str（k）,...

'HorizontalAlignment','center','VerticalAlignment','bottom'）

end

N=k;

plot（[xn（N）,xn（N）],[yn（N）,xn（N）],'b-'）

text（xn（N）,yn（N）,num2str（N）,...

'HorizontalAlignment','center','VerticalAlignment','bottom'）

N_Feed=find（xn

N_Feed=N_Feed

（1）;

text（0.5,0.5,{strcat（'所需理论板：

',num2str（N））,...

strcat（'进料板位置：

',num2str（N_Feed））},...

'HorizontalAlignment','left','VerticalAlignment','top'）

%以下代码是为了去掉顶端和右边坐标轴的刻度

boxoff

ax2=axes（'Position',get（gca,'Position'）,...

'XAxisLocation','top',...

'YAxisLocation','right',...

'Color','none',...

'XColor','k','YColor','k'）;

set（ax2,'YTick',[]）;

set（ax2,'XTick',[]）;

boxon

运行结果如图1所示：

图1图解苯-甲苯精馏塔理论塔板数和进料位置

VI