QQ登录

只需要一步,快速开始

 注册地址  找回密码
查看: 2559|回复: 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 |邮箱已经成功绑定
    8 _4 g4 C8 |, C
    题意解析:
    7 e: c; u% A, n7 Z8 }( P+ }
    ' R+ N( y; d0 Q' c0 p7 T(a) 因气体 A 与液体 B 不发生反应,故其扩散现象的质量平衡方程如下:/ T* k9 K( Y# f
    - k' L: k0 p' G1 x8 j" {

    , F9 `) T' B9 ]
    ' L# g& N4 }: _0 W( B$ {(b) 在气体 A 与液体 B 会发生一次反应的情况下,其质量平衡方程需改写为2 m' Z, L: v9 b( d0 f) r

    9 B7 X) V4 A0 x0 w! K* H8 o- E1 l7 I* t: P

    : e- L4 k. P8 {2 c2 v+ g  d而起始及边界条件同上。
    6 m$ R" N7 g% N/ Q2 S3 L
    ' |" Z8 d3 m! H, {5 w: J在获得浓度分布后,即可以 Fick’s law
    ; V3 V6 t$ R( D$ R( u6 _5 M7 z, o7 h* T2 n+ O& P0 o% w) b# |* P

    0 B8 S0 ~2 n. ~' y7 b! N
    ( a) U0 G" F* z& m0 l  Q$ m计算流通量。6 L0 [4 Q8 F% o. |# |' k( W

    9 U% @( h7 ]% L! V. B% [9 cMATLAB 程序设计: 此问题依旧可以利用 pdepe 迅速求解。现就各状况的处理过程简述如下
    $ F& s5 D0 X5 B( T+ W% n3 k
    $ N8 V6 ?9 E/ N( k: ^$ u
    2 @  ?: {) y3 u; {" t1 W
    3 V8 _1 f: U- H' h9 l  A, t9 `利用以上的处理结果,可编写 MATLAB 参考程序如下:' l7 o- Y+ @7 m& B( c$ ?2 D
    ; l( {5 T8 O8 W2 ]* R
    function ex20_3_2
    6 f% u3 R8 [+ j/ k- [' ~. _$ f4 s%*****************************; z% \9 `, o$ D: }5 m( i
    % 扩散系统之浓度分布7 C4 ?) Q( q$ b% s0 X
    %*****************************
      p, v, V( I; W* Fclear3 X; ?! }  v. c6 M' t' r
    clc
    " ?7 j8 v$ s( s5 ~; w# iglobal DAB k CA0
    6 [2 T4 r6 h9 n- D%******************************
    9 ?* ]* G/ w; l! o& J9 j% 给定数据+ I: l5 U  B; h
    %******************************, ~, s; m5 _7 s
    CA0=0.01;4 z) X  V9 I( t6 Q
    L=0.1;
    ( ^& a3 \- ?& l3 k; P& eDAB=2e-9;
    . ?& w3 G, Q4 g6 R6 P4 bk=2e-7;
    ; [) _# k' M0 Y+ X2 c, rh=10*24*3600;2 D4 v- I, r, }: w) n
    %*******************************
    3 X& Z' |' I  E% 取点
    - S. i) L/ b6 H0 B2 t%*******************************6 |  V- m! s$ }. q
    t=linspace(0,h,100);
    2 D; u6 m/ v& K3 g% cz=linspace(0,L,10);
    . D* y1 |- {" E  d) B%******************************** g2 W' v/ u0 O
    % case (a)0 {1 L5 `1 P9 x$ z/ T" }
    %*******************************
    ! F( i$ V: M- ?# Fm=0;
    $ B; U; @) G) F# u9 k/ L2 v, jsol=pdepe(m,@ex20_3_2pdefuna,@ex20_3_2ic,@ex20_3_2bc,z,t);
    : P% f9 c/ G( U* rCA=sol(:,:,1);9 e2 X# Y( k! x  z6 v( E
    for i=1:length(t)
    $ k9 ^. M7 i+ r+ S, s, g [CA_i,dCAdz_i]=pdeval(m,z,CA(i,,0);
    5 R4 U6 X' C; c& l% w NAz(i)=-dCAdz_i*DAB;
    9 H! `' Y# ~. X) P$ G% E+ Gend
    / z% j8 ^$ O5 l. S& |figure(1)
    . f9 p# {6 d8 U" Lsubplot(211)
    - f# F8 p& l% x1 v! ]8 J+ lsurf(z,t/(24*3600),CA). c1 i9 H) A9 w( Z  b$ V  m
    title('case (a)')
    & U2 M7 M2 v; o8 J% z% q. o! o% K# H5 qxlabel('length (m)')
    6 D( c& c& E7 R# x! l- jylabel('time (day)')2 M7 m& f' `% d3 h. a
    zlabel('conc. (mol/m^3)')) D/ e% l- z+ L$ k
    subplot(212)
    + T) a& ]1 m# C/ v- _" nplot(t/(24*3600),NAz'*24*3600)) W9 x: h; o, _6 ?8 o
    xlabel('time (day)')% m$ f6 b+ E1 B) G1 `
    ylabel('flux (mol/m^2.day)')
    & Q  j: g5 l" D/ s( ?6 n%************************************" z" s7 \; {4 W1 I
    % case (b)
    # s9 C" Z; a1 B, P%************************************6 G* ?' g# t3 y: p) {' f) X4 q
    m=0;' k3 ]+ `. \0 P% p; x2 _" y' g( e2 P6 w
    sol=pdepe(m,@ex20_3_2pdefunb,@ex20_3_2ic,@ex20_3_2bc,z,t);, g; B2 @* j" G0 Q  b' O
    CA=sol(:,:,1);3 s% e0 F6 S: B2 [1 G0 M) g
    for i=1:length(t)( h, e3 |, f* b/ g
    [CA_i,dCAdz_i]=pdeval(m,z,CA(i,,0);& y( N9 H- |. e, P& A+ G6 w
    NAz(i)=-dCAdz_i*DAB;
    % G3 k& J& N# S6 vend
    8 }! r" p& @% Y' c%
    , [+ k! u( g; C0 ]figure(2)
    ' B8 Q- `( n) Rsubplot(211)
    8 u! w8 H% c* jsurf(z,t/(24*3600),CA)
    5 {& p1 }4 A4 k) @" F* }; utitle('case (b)'), o- L+ H& Y9 t
    xlabel('length (m)')
    , O& x. p6 r% e' U8 vylabel('time (day)')
    / u0 y' ^5 ?! m0 z1 W  tzlabel('conc. (mol/m^3)')0 h* ?$ S: I) }$ |0 i" \- v+ R: Z
    subplot(212)
    " n* I7 l' m# i  Oplot(t/(24*3600),NAz'*24*3600)( m# c" p4 V: ?6 u# J0 Z' n
    xlabel('time (day)')
    # [' W' n2 R+ U3 S' J! hylabel('flux (mol/m^2.day)')* \' Q, x4 d7 R# k( |9 M
    %********************************************: V: d# X3 V- K4 R. ]: N
    % PDE 函数
    # y! F) d' K/ W# \" X8 N- h/ ^+ r%********************************************
    : C7 U! ~* p% f* m# b% case (a)
    : v# n% f1 @& v  r%********************************************
    ! f) @% [% z2 O- }% qfunction [c,f,s]=ex20_3_2pdefuna(z,t,CA,dCAdz)' d2 D- f& W  a/ ^
    global DAB k CA0
    * p; x* n/ F" Y1 u7 }3 W$ `c=1;/ z/ i+ \+ G* c" {  s/ g3 N) c& L
    f=DAB*dCAdz;& O& _' T, c, p1 N
    s=0;) E; C2 W! N1 k3 i- k1 x3 i2 _
    %*********************************************( N- R; l6 U$ j5 x/ T9 D
    % case (a)) f3 c0 g, ]1 t7 l2 k
    %*********************************************
    ( V& ~7 `: {: F6 v' W7 @. O" Lfunction [c,f,s]=ex20_3_2pdefunb(z,t,CA,dCAdz)+ Y! f# t! ^- k! D; q
    global DAB k CA01 X: |3 Y, Q+ A3 O! F
    c=1;% L( Y$ ]  C* D. I$ v5 G3 T) W
    f=DAB*dCAdz;8 Q0 b5 @2 o6 b- Y
    s=k*CA;
    / B' R8 D) K9 m9 z- l& V: j$ f0 @%**********************************************- a9 K* l1 d; K
    % 初始条件函数
    ) n- ^8 I6 j* {! k/ w8 s%**********************************************: d. {9 G+ n& Z6 ^
    function CA_i=ex20_3_2ic(z)
    ; x) f! v8 m, r$ r% S& U8 bCA_i=0;
    * a2 t: [" j& X( t2 K%************************************************ 3 N1 x, T/ X, M
    % 边界条件函数
    8 ~1 B1 a0 v8 V5 {8 l% @+ _- T%************************************************% r6 ]7 E3 D, M& }/ L$ q: r
    function [pl,ql,pr,qr]=ex20_3_2bc(zl,CAl,zr,CAr,t)" \' T, v; ^- m1 \# Q8 E7 B
    global DAB k CA0
    : V  Y0 _/ R* `3 L8 b8 `pl=CAl-CA0;
    # U2 h2 t: ?% P/ e) Hql=0;* _( d* n- ]' T/ o, z4 u
    pr=0;5 X# k# L) `) @: G. T
    qr=1/DAB;
    2 M" {% r9 W7 p9 M5 {! Q* k1 C# q7 l/ l5 Z3 a* B0 `
    ————————————————
    : f$ N' }/ ^, x2 ?* L3 J版权声明:本文为CSDN博主「wamg潇潇」的原创文章,遵循CC 4.0 BY-SA版权协议,转载请附上原文出处链接及本声明。- X  T2 j/ u: T6 T$ F/ u
    原文链接:https://blog.csdn.net/qq_29831163/article/details/89711694, A, T4 \; m) a5 \
    $ r- w* r7 C' k# k3 o7 N9 K
    5 e: ?( u! r/ v5 R9 L  l* q
    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-31 07:50 , Processed in 0.321997 second(s), 50 queries .

    回顶部