porousspherecat_1ord.m

matlab在化学工程中的应用实例程序,对大家研究学习有指导作用。

function PorousSphereCat_1Ord
% 球形催化剂中的一级不可逆等温反应的模拟计算
% Diffusion with Simultaneous Reaction in a porous sphere catalyst pellet
%
%   Author: HUANG Huajiang
%   Copyright 2003 UNILAB Research Center, 
%   East China University of Science and Technology, Shanghai, PRC
%   $Revision: 1.0 $  $Date: 2003/05/28 $

clear all
clc

global fai Cs

Cs = 1;
rp = 5e-3/2;        % m
k = 4;              % 1/s
D = 0.05e-4;        % m^2/s
fai = rp*sqrt(k/D)  % fai: Thiele模数
a = 0;
b = 1;

% 求解ODE-BVP问题
% initialize of solution with a guess of y1(r)=0, y2(r)=0 
solinit = bvpinit(linspace(a,b,100),[0 0]);
sol = bvp4c(@ODEs,@BCfun,solinit);

% 绘制浓度分布图
plot(sol.x,sol.y(1,:),'b-')
xlabel('Dimensionless Radius')
ylabel('Dimensionless Concentration')

% 计算催化剂的有效因子
I1 = quadl(@IntFunc1,a,b,[],[],sol.x,sol.y(1,:))
I2 = quadl(@IntFunc2,a,b)
eta = I1/I2
fprintf('\t催化剂的有效因子为: η= %.4f',eta)

% ------------------------------------------------------------------
function dydr = ODEs(r,y)
global fai Cs
dy1dr = y(2);
if r == 0           % 避免除数为零
    dy2dr = 0;
else
    dy2dr = fai^2*y(1)/Cs - 2/r*y(2);
end
dydr = [dy1dr; dy2dr];

% ------------------------------------------------------------------
function bc = BCfun(ya,yb)
bc = [yb(1)-1; ya(2)];

% ------------------------------------------------------------------
function f = IntFunc1(r,ri,yi)
f = spline(ri,yi,r) .* r.^2;

% ------------------------------------------------------------------
function f = IntFunc2(r)
global Cs
f = Cs * r.^2;