QQ登录

只需要一步,快速开始

 注册地址  找回密码
查看: 2556|回复: 0
打印 上一主题 下一主题

[建模教程] 偏微分方程的数值解(四): 化工应用————扩散系统之浓度分布

[复制链接]
字体大小: 正常 放大
浅夏110 实名认证       

542

主题

15

听众

1万

积分

  • TA的每日心情
    开心
    2020-11-14 17:15
  • 签到天数: 74 天

    [LV.6]常住居民II

    邮箱绑定达人

    群组2019美赛冲刺课程

    群组站长地区赛培训

    群组2019考研数学 桃子老师

    群组2018教师培训(呼伦贝

    群组2019考研数学 站长系列

    跳转到指定楼层
    1#
    发表于 2020-6-10 10:29 |只看该作者 |正序浏览
    |招呼Ta 关注Ta |邮箱已经成功绑定
    ) Y" t- n5 x) E# X( \  I! K, o
    题意解析:0 o  b  O) ~8 Y

    $ y: ?- w/ F2 p. u. T. l3 j(a) 因气体 A 与液体 B 不发生反应,故其扩散现象的质量平衡方程如下:
    - J: O; o6 _+ p6 [) R
    * u( @' [, `5 h1 Q/ b: J) f7 M- {
    / ^4 D# u0 J4 v9 [; i# y# Z
    (b) 在气体 A 与液体 B 会发生一次反应的情况下,其质量平衡方程需改写为
    ( m3 f1 F* {# P4 f5 B% `/ T8 u8 N( m! U! g# E# Z' ]

      j* T) y# @) s4 H
    8 U' r) u$ T- f& q$ k1 i$ e而起始及边界条件同上。
    + W! v2 Y, y7 A% ^9 O8 R' r% _* U
    ( S" ]; A0 a4 y1 x在获得浓度分布后,即可以 Fick’s law1 w2 j. V9 n5 R1 R) y2 J
    : _* g" }( Q9 ]5 [: I+ G% n

    3 w0 `* Q4 W* C% D* u& Q1 w2 H5 C) x9 H5 s: i3 i
    计算流通量。
    : N- @1 q* M5 f' Q. i8 E! c5 Z" u( t/ z# a9 }- E  ~
    MATLAB 程序设计: 此问题依旧可以利用 pdepe 迅速求解。现就各状况的处理过程简述如下5 h2 [$ J) c  O- A4 n

    & Y5 j/ Q, u; }; O: i- ]
    1 Z! m) A- r% r# o6 c- G6 J
    + M$ l. L8 g+ r  H6 m" p8 y; a2 y利用以上的处理结果,可编写 MATLAB 参考程序如下:+ ~1 K, n% M0 S+ r  Q
    & Z/ z3 g" o& K# b8 Q" b
    function ex20_3_2
    ) M) S5 N8 f' W: [# N%*****************************1 `% `4 w; _) g) S$ G# i
    % 扩散系统之浓度分布
    ; J; r! v+ d+ n. y%*****************************
    ' ]0 I1 u) Q7 }clear( l/ e! [, n! A" h, w
    clc$ B" {% z4 e( C: Y  f) ^9 n: M
    global DAB k CA0
    8 V" b2 h! H8 L: y%******************************: Q. s( j: `7 G! `: G
    % 给定数据# A, V' j# ?" w# U
    %******************************' l$ I9 N: H) w( b; J) o
    CA0=0.01;. t: @- s. M0 g/ l7 W7 }/ J
    L=0.1;
    4 u: Y. j+ o* e3 |3 ^- n- U: zDAB=2e-9;0 u; F! ~; ~3 o- K# h# F1 r! ^# y  I: z
    k=2e-7;. @9 ?% P/ m- l/ R. ~5 g
    h=10*24*3600;
    - b2 H& \  K3 z  t9 m%*******************************
      y% r+ T4 O8 a# G0 v  {! l% 取点  B6 S" s% k. I4 Z2 L1 O4 A
    %*******************************  A$ i9 g- j, a
    t=linspace(0,h,100);& n* ~+ w3 }7 s
    z=linspace(0,L,10);& o" U% _, ~9 i) x6 I
    %*******************************4 b1 k3 l8 q3 N6 S0 N& {
    % case (a)
    & _( N% d/ Y$ j/ Z0 O& V! z%*******************************
    $ S8 e; [8 Z' Q* xm=0;5 L$ G( g8 V" T
    sol=pdepe(m,@ex20_3_2pdefuna,@ex20_3_2ic,@ex20_3_2bc,z,t);
    ; B* A! l$ J! C# a. z: g0 \CA=sol(:,:,1);: ~2 X* j4 O( E  S; S( n
    for i=1:length(t)
    6 B, Z' e/ Z6 }; j [CA_i,dCAdz_i]=pdeval(m,z,CA(i,,0);
    ! X- u) A6 y4 E NAz(i)=-dCAdz_i*DAB;
    9 `) l5 @3 b5 Z  b( \0 h- V. I. K* @end1 J% _' t) H% K5 j, C0 z
    figure(1)" k3 ~) |( z5 ^7 P# V
    subplot(211). w$ v- ^. F5 S/ p% _
    surf(z,t/(24*3600),CA)( l# [& p5 z9 B8 B6 c1 y% T. T
    title('case (a)') + S: K' [2 j' d, x& H; }8 J
    xlabel('length (m)')
    6 \, p! ]( w9 E4 l1 Dylabel('time (day)')+ h. @9 k3 e( p
    zlabel('conc. (mol/m^3)'): p. U2 P3 E+ b
    subplot(212)  F5 t( S/ \' F
    plot(t/(24*3600),NAz'*24*3600)
    7 S. o$ T' {' W. ?xlabel('time (day)')
    * O, ^, d8 I9 g. k5 t& C8 ?0 w; dylabel('flux (mol/m^2.day)')3 s: a; B, E6 ?* g5 [- B! j! S
    %************************************
    + v' F( U9 b# V6 Y- M" s4 c% case (b)( K% b8 O  A9 d5 U4 S, }8 m' \6 U( T9 T
    %************************************
    * ~4 }  U# X; U" J, l! e2 Gm=0;
    & l9 U) L1 ~! ~0 O) |3 Y' Msol=pdepe(m,@ex20_3_2pdefunb,@ex20_3_2ic,@ex20_3_2bc,z,t);( v  Q* z( l* _, w& d
    CA=sol(:,:,1);2 z. l7 U( U7 @/ A- r; Y( e& c8 W
    for i=1:length(t)4 l; Y$ Q7 P' X7 ]
    [CA_i,dCAdz_i]=pdeval(m,z,CA(i,,0);! H' b# |4 x5 f6 [# g
    NAz(i)=-dCAdz_i*DAB;. j) t! K3 h! l* W
    end
    . M; N5 `# f' Q* Y& q' U! H%8 N, p" }4 a8 |/ a" O: J; |
    figure(2)
    ; m2 q8 L! U  R/ j. v+ H  q. nsubplot(211)# a  {% O$ [/ l
    surf(z,t/(24*3600),CA)% i5 W4 l0 |4 _
    title('case (b)')
    * A. f8 }5 i1 Z$ }+ U3 F8 exlabel('length (m)')2 X( q) m& z" g8 u0 G) x/ y
    ylabel('time (day)')
    : L  J  f6 \& E2 {' ~7 R* gzlabel('conc. (mol/m^3)')
    0 }) `1 g5 A  C. Y: Usubplot(212)
    ( [& G# u9 C% V6 splot(t/(24*3600),NAz'*24*3600)
    3 a, }0 y: P# x# x$ ]0 u. nxlabel('time (day)')
    + n& K! {5 r: J: `2 h! Z' D7 D3 oylabel('flux (mol/m^2.day)')- R/ h- ]" t& V6 i* t2 }+ a
    %********************************************, h+ `5 S3 l* A4 n
    % PDE 函数0 r5 X5 Y5 L* _3 k4 z5 I4 g- @+ ^; A& ?
    %********************************************& S  A& ~+ a0 f/ G
    % case (a)1 l( O* i3 j0 _( L; A+ K* [7 K
    %********************************************7 L+ X# D$ m6 W* L# T3 E
    function [c,f,s]=ex20_3_2pdefuna(z,t,CA,dCAdz)/ I7 m# D! T, X! |7 K. J  e
    global DAB k CA0
    7 V3 n, S! R2 z" K7 Uc=1;
    # S# E6 {4 r* kf=DAB*dCAdz;
    1 `- O3 l/ S1 b0 h' r$ j' Ps=0;% V7 \1 U5 h. {4 Q! _2 t
    %*********************************************3 C6 O; H* ?2 b+ U
    % case (a)
    ! U6 A: T. R: `3 h%*********************************************" V$ O1 N- I8 E2 F$ K
    function [c,f,s]=ex20_3_2pdefunb(z,t,CA,dCAdz): V: A/ {% \3 q' I& I
    global DAB k CA0
    , _( E: [* J7 P& p0 {' ?: d7 Ec=1;5 @# M' l# K- L! S  ^0 m) U3 _
    f=DAB*dCAdz;7 Q- i, @9 c1 |/ k
    s=k*CA;1 s* I2 e0 w1 L
    %**********************************************. X6 M+ Y! B* N% c: J" x
    % 初始条件函数
    / q. y& q2 ~8 p; w$ P: i- Z%**********************************************
    & m+ ~& {+ ?* m6 sfunction CA_i=ex20_3_2ic(z)
    2 [" ?8 ?* c. N5 k" Q* s, w& zCA_i=0;
    + @7 x. M; g- f7 Z8 `0 k6 W5 h%************************************************ + _5 p( ^, c+ T+ T+ `$ a
    % 边界条件函数) S# Z; f: I1 Z; E( b
    %************************************************
    8 q  K) E0 F( S) s: H% `# p4 {function [pl,ql,pr,qr]=ex20_3_2bc(zl,CAl,zr,CAr,t)
    5 \7 f' n9 T0 N1 }: |global DAB k CA0) P; H0 T! t1 C* V
    pl=CAl-CA0;2 X6 N' [6 K. ?
    ql=0;$ R: y* d% e5 {$ M
    pr=0;
    ' r+ w/ m- L  L( b! q" v( N5 V' Vqr=1/DAB;
    7 q$ G# N8 G, C. D  Q
    , y+ J  v7 [7 `8 B# ^9 I————————————————
    ; d. I& i% z7 s6 j6 Z版权声明:本文为CSDN博主「wamg潇潇」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。
    + w  u# Q2 c9 g原文链接:https://blog.csdn.net/qq_29831163/article/details/89711694
    ; e  t' u' K) g+ S! N
    8 m& b8 }: g/ q: K  F& }) E9 i9 z* E. L
    zan
    转播转播0 分享淘帖0 分享分享0 收藏收藏0 支持支持0 反对反对0 微信微信
    您需要登录后才可以回帖 登录 | 注册地址

    qq
    收缩
    • 电话咨询

    • 04714969085
    fastpost

    关于我们| 联系我们| 诚征英才| 对外合作| 产品服务| QQ

    手机版|Archiver| |繁體中文 手机客户端  

    蒙公网安备 15010502000194号

    Powered by Discuz! X2.5   © 2001-2013 数学建模网-数学中国 ( 蒙ICP备14002410号-3 蒙BBS备-0002号 )     论坛法律顾问:王兆丰

    GMT+8, 2026-7-30 14:35 , Processed in 0.350540 second(s), 51 queries .

    回顶部