数学建模社区-数学中国
标题:
偏微分方程的数值解(三): 化工应用实例 ----------触煤反应装置内温度及转换率的分布
[打印本页]
作者:
浅夏110
时间:
2020-6-10 10:27
标题:
偏微分方程的数值解(三): 化工应用实例 ----------触煤反应装置内温度及转换率的分布
例 4 触煤反应装置内温度及转换率的分布
$ v2 I+ Q2 y( A4 I4 Q7 m
# l( U7 E" z) t: _
以外部热交换式的管形固定层触煤反应装置，进行苯加氢反应产生环己烷。此反应 系统之质量平衡及热平衡方程式如下：
: |. C: f& D6 E' |; F; h) ]
4 Y) T; s9 y/ x- K
5 \. ~* q4 I: b: q, ~0 j- ], q
% C. K% S6 T4 v
其中T 为温度(℃)， f 为反应率，L 为轴向距离，r 为径向距离。此系统的边界条件为
3 G- u$ G, d y* d4 q' m9 f
" W3 R2 T$ U" F" Q- M0 V* Z* v
3 _0 P1 G- a% N8 T+ }# ]5 r2 ]
+ h/ S; q- I) L
此外，式中之相关数据及操作条件如下：
$ ` o' l3 P* @ R8 N
! S0 M+ Q5 [' n1 {$ ?) k
(i)反应速率式
& T# r9 j* V# Z: X) R
4 s+ U+ e; O5 m/ R
- Y) a8 L9 p2 Q* Q9 j
# M S( O, f1 c& }+ e
其中 P 表示分压(atm)，而速率参数为
, O/ D U# }2 Z9 J x1 W% j1 q
* u0 r \& _9 z2 J: u* H! m
/ g M9 ~) A7 Q
# U9 z: `1 |, V& V
上式中，下标 B，H 及 C 分别代表苯，氢及环己烷。R 为理想气体常数(1.987cal/mol·K)。
9 G" y. l# M) J1 T
) \& S ~% L) H* w3 A- Y
(ii)操作条件及物性数据
3 Z( V4 W2 ?# W' m* s
+ I& ^9 S+ u# F" v& E' u
8 R! P- @* |8 G+ g6 ]$ v
2 ]* V! c% q% ]2 K6 f
; h4 x% ]/ q2 x& C+ _* \# X
0 a/ C- i) i# Z% t ~1 T' d
题意解析：
$ T& o9 r: W: F3 S
" R9 H# u7 Q, \/ z+ T
3 j9 y& Z$ X5 `* @1 c" C
( G5 e' w, D2 D+ T$ [
" \2 x/ Y" ~. H- Q( L
将上式，连同反应速率式，带入平衡方程式中，配合边界条件，可利用 pdepe 求解。
- U. K+ l. e) _; ^; y4 p1 u
Q4 K4 a8 F h. s; x6 C6 ]
MATLAB 程序设计 将原方程改写成如式（35）的标准式
. t* K# ^- Z9 C3 V/ V/ J9 G& y
7 ?! t$ e% l2 c2 d: s$ _
# J- P" {) P/ U1 Q4 i# ^
7 G0 h9 V. D! {) q: h+ ~; i7 P
因此
: U+ @! l( _. U7 {- w W
& m; m2 b7 c1 I: s4 n
$ C3 A" e* I8 z& _7 b" F* ~- N' ]
2 C6 ~- \. ?1 _0 Q
根据以上的分析，可编写 MATLAB 程序求解此 PDE 问题，其参考程序如下：
. t: s- X2 J) P' Z8 X
4 A2 m3 j- Z: w* v
function ex60_3_1
6 Y d' k5 G9 a( e
%******************************
/ T1 b7 r2 j$ ^, R5 l( U& v5 z$ y6 ?
% 触媒反应器内温度及转化率的分布
) X1 J: Y% @# x. ^5 b* U
%******************************
6 U; o" k0 H5 M1 O2 b# s4 Q
global Pt rw Tw G M y0 Mav rho_B Cp dHr h0 u R ke hw De
5 i) K* A- L6 p5 X% k$ Q& E4 J
%******************************
# [ t, s( `( W7 b# T5 k
% 给定数据
, x" o" L& i! j3 r) M. A; `' O
%******************************
; K) B- V: e. I. z
Pt=1.25; %总压(atm)
# X N$ p0 D" U' X# o" w5 }
rw=0.025; %管径(m)
# }$ @4 i9 }( X8 ?, W) G, V H/ N
Tw=100+273; %壁温(℃)
7 d8 l6 O2 ^- E1 L
G=631; %质量流率(kg/m2hr)
$ i m9 u+ n7 E' Q
M=30;
9 r0 ?# |8 Q) H- K/ S6 V. N
y0=0.0323;
4 h/ L, j; t/ Q" b; m1 {
Mav=4.47;
% Y5 L) _2 g6 N' v
rho_B=1200;
, ?8 o- i: O# y7 G
Cp=1.74;
, v9 m1 v8 d& {/ Q
dHr=-49250;
2 j; V# c2 M; e, Y
h0=65.8;
1 h- }0 U+ x3 Q7 ^9 e. O# I
T0=125+273;
: `/ q; a* @4 \# F
Lw=1;
1 _4 _0 }: V8 V6 u0 q
u=8.03;
- J- a. J$ N" h2 a7 C7 Z8 G0 F
R=1.987;
6 x) P. N* f* l7 d
ke=0.65;
$ E6 f! Y6 N; i8 P8 R$ z2 X7 v
hw=112;
5 _) X' p; q5 J9 P- Z# x
De=0.755;
3 W. s6 |9 F9 x: m3 {( L
%********************
& B0 j' Y1 @$ w! }0 i6 k) S
m=1;
8 c+ M) r( {( [; k) |
%********************
V$ c) |# h4 O- \5 k6 r4 J
% 取点
. Q0 N% v% z3 y+ ^0 H3 {" G3 J
%********************
! H n4 V/ u+ Z
r=linspace(0,rw,10);
0 G) q$ c# [/ n/ I9 u8 D4 o. m
L=linspace(0,Lw,10);
7 F1 n7 e! O' {) U1 b3 a4 J+ W
%***********************
, q# A1 `% _& g% _, i
% 利用 pdepe 求解
" I' R- x, J( ~3 n' e+ z3 R
%***********************
% c; h; [& Q) P4 @1 M
sol=pdepe(m,@ex20_3_1pdefun,@ex20_3_1ic,@ex20_3_1bc,r,L);
5 L! \7 a. A& a! b$ v2 c
T=sol(:,:,1); %温度
! Q6 Z; A& X0 u
f=sol(:,:,2); %反应率
( c0 V% i7 U4 ]' V# Y
%***********************
- E5 R9 K2 g# W: }( A, P$ s+ r
% 绘图输出
6 T4 _: G6 s! }; c7 x9 O- S
%***********************
' o& G: t( f- ]$ g S
figure(1)
! W* W6 g# b& g
surf(L,r,T'-273)
8 P" }8 j, n1 P
title('temp')
" {" v9 v) Q7 S, @. M
xlabel('L')
- v" [4 x- `1 S3 l' {/ Z& a
ylabel('r')
' s* ?6 D& ?+ I5 O6 D
zlabel('temp (0C)')
5 W8 {. w8 E/ x5 i$ ]+ X
%
8 \4 Z3 Q1 R; ~1 E6 m0 U% ~1 n7 T1 v
figure(2)
* W9 Y" K- l# m. [8 Q$ h1 ]: D- M
surf(L,r,f')
& ]# N) h: q+ J' ?/ J$ ], C8 G
title('reaction rate')
9 ?% }/ I; o9 _
xlabel('L')
* ~0 ]" S, |; n# K! `
%初始条件函数
& N7 z5 X, G( ~; B0 H6 b3 F
%**********************************
7 R2 z$ m- k6 p
function u0=ex20_3_1ic(x)
4 O2 H, y7 d% i( q
u0=[125+273 0]';
! I0 k8 p" x) @ G% W+ d1 \
%**********************************
# j/ }, E9 Z8 G1 V9 C; l
% 边界条件档
/ T2 X- e! s- x M
%**********************************
! c9 S* |1 o( u. Y
function [pl,ql,pr,qr]=ex20_3_1bc(rl,ul,rr,ur,L)
$ |: X0 r, q2 ?: ^ X
global Pt rw Tw G M y0 Mav rho_B Cp dHr h0 u R ke hw De
& F- a. V# E$ \
pl=[0 0]';
2 A f* L. C; i" L8 \; P) e( R
ql=[1 1]';
% {0 [$ {0 H t" a
pr=[hw*(ur(1)-Tw) 0]';
6 x0 {- W$ W' [. X
qr=[G*Cp 1]';
! @! a3 s$ y+ H8 G: W3 o R
ylabel('r')
9 v9 N9 d( g# w' Z2 m1 a; x' U
zlabel('reaction rate')
$ ?# u, u1 i" s0 `$ {4 K
%*************************************************
! I) e: S0 X; _
% PDE 函数
9 ~" b- a- M! k( x1 e
%*************************************************
2 n- Y) Y# ~$ [3 {% B/ n
function [c1,f1,s1]=ex20_3_1pdefun(r,L,u1,DuDr)
0 O |6 h$ L1 G1 W4 n
global Pt rw Tw G M y0 Mav rho_B Cp dHr h0 u R ke hw De
6 ^6 V `+ w7 P+ R9 G
T=u1(1);
% S, N3 T7 M# I Z& k* d2 S' Q
f=u1(2);
# d8 {& |0 A+ q
%
) q% C5 w. u8 G" x6 w1 s. ^* K }
k=exp(-12100/(R*T)+32.3/R);
( A! b( a4 z, y/ ?7 v2 O
Kh=exp(15500/(R*T)-31.9/R);
5 B" ], F3 V. D
Kb=exp(11200/(R*T)-23.1/R);
R0 y) @! y; z. c" I: j; A
Kc=exp(8900/(R*T)-19.4/R);
" Y3 A+ y" a- A0 t
%
& w: @- A- [" p0 I
a=1+M-3*f;
0 K: \3 U7 j. N$ F3 d' d
ph=Pt*(M-3*f)/a;
" }# w+ j) z* I+ @2 m0 g
pb=Pt*(1-f)/a;
0 B3 m7 m. O, H% y% F
pc=Pt*f/a;
: K, e! q" A8 C8 W' N# m
%
; D+ k# [* y$ f1 Q8 d
rA=k*Kh^3*Kb*ph^3*pb/(1+Kh*ph+Kb*pb+Kc*pc)^4;
0 L3 @" Y1 S1 G5 G
%
3 c' I, u& G+ R- s5 M; d9 O1 `
c1=[1 1]';
4 g% b8 v# V O
f1=[ke/(G*Cp) De/u]'.*DuDr;
, ]! G: I2 k2 G5 T0 Y: G# Q
%s1=[ke/(G*Cp*r)*DuDr(1)-rA*rho_B*dHr/(G*Cp)-2*h0*(T-Tw)/(rw)
- M( I0 G; B2 r. P
s1=[-rA*rho_B*dHr/(G*Cp);rA*rho_B*Mav/(G*y0)];
N4 }0 U0 q+ W1 Y9 ~
%**********************************
! V* x3 ?8 {: U0 k6 n* ~+ L
& Q" ^ L& V2 v% z5 z
————————————————
$ O! Z8 k- L; o7 @7 `' k
版权声明：本文为CSDN博主「wamg潇潇」的原创文章，遵循CC 4.0 BY-SA版权协议，转载请附上原文出处链接及本声明。
8 ]( ?3 l# l) w2 N6 l
原文链接：https://blog.csdn.net/qq_29831163/article/details/89711536
2 K9 B% d, O$ {/ |
( v) j6 m% M: [5 F" |- r$ T
) o1 f6 P4 Q3 `8 o( ^% o& y
欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)
Powered by Discuz! X2.5
5 \. ~* q4 I: b: q, ~0 j- ], q
8 R! P- @* |8 G+ g6 ]$ v
; h4 x% ]/ q2 x& C+ _* \# X
3 j9 y& Z$ X5 `* @1 c" C
# J- P" {) P/ U1 Q4 i# ^