【正文】 A+ B→C ， 其反應(yīng)速率常數(shù)。氣體 A 在液體 B 中之?dāng)U散系數(shù)為。管中儲放靜止液體 B，高度為 L=10 ㎝，放置于充滿 A 氣體的環(huán)境中。 qr=[G*Cp 1]39。 ql=[1 1]39。 %********************************** %初始條件函數(shù) %********************************** function u0=ex20_3_1ic(x) u0=[125+273 0]39。 %s1=[ke/(G*Cp*r)*DuDr(1)rA*rho_B*dHr/(G*Cp)2*h0*(TTw)/(rw) s1=[rA*rho_B*dHr/(G*Cp)。 f1=[ke/(G*Cp) De/u]39。 % rA=k*Kh^3*Kb*ph^3*pb/(1+Kh*ph+Kb*pb+Kc*pc)^4。 pb=Pt*(1f)/a。 % a=1+M3*f。 Kb=exp(11200/(R*T))。 % k=exp(12100/(R*T)+)。) %************************************************* % PDE 函數(shù) %************************************************* function [c1,f1,s1]=ex20_3_1pdefun(r,L,u1,DuDr) global Pt rw Tw G M y0 Mav rho_B Cp dHr h0 u R ke hw De T=u1(1)。) zlabel(39。) ylabel(39。) xlabel(39。) title(39。temp (0C)39。r39。L39。temp39。 %反應(yīng)率 %*********************** % 繪圖輸出 %*********************** figure(1) surf(L,r,T39。 T=sol(:,:,1)。 L=linspace(0,Lw,10)。 %******************** m=1。 hw=112。 R=。 Lw=1。 h0=。 Cp=。 Mav=。 %質(zhì)量流率(kg/m2hr) M=30。 %管徑(m) Tw=100+273。另外，左邊界條件( r = 0處)寫成 即 同理右邊界條件()可寫成 即 根據(jù)以上的分析，可編寫 MATLAB 程序求解此 PDE 問題，其參考程序如下： function ex60_3_1 %****************************** % 觸媒反應(yīng)器內(nèi)溫度及轉(zhuǎn)化率的分布 %****************************** global Pt rw Tw G M y0 Mav rho_B Cp dHr h0 u R ke hw De %****************************** % 給定數(shù)據(jù) %****************************** Pt=?；谝韵碌姆磻?yīng)方程 則各分壓與總壓之關(guān)系為 將上式，連同反應(yīng)速率式，帶入平衡方程式中，配合邊界條件，可利用 pdepe 求解。(ii)操作條件及物性數(shù)據(jù) 題意解析： 因反應(yīng)速率式 A r 與分壓有關(guān)，而分壓又與反應(yīng)率 f 有關(guān)。R 為理想氣體常數(shù)(此反應(yīng) 系統(tǒng)之質(zhì)量平衡及熱平衡方程式如下： 其中T 為溫度(℃)， f 為反應(yīng)率，L 為軸向距離，r 為徑向距離。 qr=[0 1]39。 ql=[1 0]39。 %**************************************** %初始條件函數(shù) %**************************************** function u0=ex20_2ic(x) u0=[1 0]39。 F=exp(*y)exp(*y)。.*dudx。) %*************************************** %pde 函數(shù) %*************************************** function [c,f,s]=ex20_2pdefun(x,t,u,dudx) c=[1 1]39。) ylabel(39。) xlabel(39。) % figure(2) surf(x,t,u2) title(39。) ylabel(39。) xlabel(39。 %第二個狀態(tài)之?dāng)?shù)值解輸出 %************************************* %繪圖輸出 %************************************* figure(1) surf(x,t,u1) title(39。 u1=sol(:,:,1)。 t=[0 1 2]。此問題的參考程序如下： function ex20_2 %***************************************%求解一維偏微分方程組的一個綜合函數(shù)程序 %*************************************** m=0。 t=[0 1 2]。同時對于 t 小處亦可取密一些。 步驟 5： 取點。 pr=[ur(1)1 0]39。 步驟 4：編寫邊界條件函數(shù) function [pl,ql,pr,qr]=ex20_2bc(xl,ul,xr,ur,t) pl=[0 ul(2)]39。 s=[F F]39。 y=u(1)u(2)。 f=[ ]39。寫成 即同理，右邊界條件( x =1處)為 即步驟 2：編寫偏微分方程的系數(shù)向量函數(shù)。 (i)初值條件 (ii)邊值條件 解 步驟 1：改寫偏微分方程為標(biāo)準(zhǔn)式 因此 和m = 0。 例 3 試解以下聯(lián)立的偏微分方程系統(tǒng) 其中且0 ≤ x ≤ 1和t ≥ 0。 pr=pi*exp(t)。 %****************** %邊界條件函數(shù) %****************** function [pl,ql,pr,qr]=ex20_1bc(xl,ul,xr,ur,t) pl=ul。 s=0。) %****************** %pde 函數(shù) %****************** function [c,f,s]=ex20_1pdefun(x,t,u,dudx) c=pi^2。) ylabel(39。) xlabel(39。 %繪圖 title(39。 %輸出點位置 [uout,dudx]=pdeval(m,x,u(M,:),xout)。) %***************** %t=tf=2 時各位置之解 %***************** figure(3) M=length(t)。 ) zlabel(39。) ylabel(39。) xlabel(39。 title(39。) %************* %與解析解做比較 %************* figure(2) surf(x,t,exp(t)39。 ) zlabel(39。) ylabel(39。) xlabel(39。 %取出答案 %************ %繪圖輸出 %************ figure(1) surf(x,t,u) title(39。 %tspan %************ %以 pde 求解 %************ sol=pdepe(m,ex20_1pdefun,ex20_1ic,ex20_1bc,x,t)。 x=linspace(0,1,20)。) 綜合以上各步驟，可寫成一個程序求解例 2。) ylabel(39。) xlabel(39。 %繪圖 title(39。 %輸出點位置 [uout,dudx]=pdeval(m,x,u(M,:),xout)。 %繪成圖 2 M=length(t)。) 若要顯示特定點上的解，可進(jìn)一步指定 x 或 t 的位置，以便繪圖。 ) zlabel(39。) ylabel(39。) xlabel(39。 surf(x,t,u) title(39。 步驟 7 顯示結(jié)果。 m=0。 %x 取 20 點 t=linspace(0,2,5)。 步驟 5 取點。 pr=pi*exp(t)。在編寫之前，先將邊界條件改寫成標(biāo)準(zhǔn)形式，如式（37）， 找出相對應(yīng)的 p(.)和q(.)函數(shù)，然后寫出 MATLAB 的邊界條件函數(shù)，例如，原邊界條 件可寫成 即 pl = u(0,t), ql = 0,和 因而，邊界條件函數(shù)可編寫成 function [pl,ql,pr,qr]=ex20_1bc(xl,ul,xr,ur,t) pl=ul。 function u0=ex20_1ic(x) u0=sin(pi*x)。 s=0。 function [c,f,s]=ex20_1pdefun(x,t,u,dudx) c=pi^2。此即 和m = 0。當(dāng)完成以下各步驟后，可進(jìn)一步將其匯總為一主程序 ，然后求解。 ref. Keel,. and M. Berzins,“A Method for the Spatial Discritization of Parabolic Equations in One Space Variable”，SIAM J. Sci. and Sat. Comput.,1990. 以下將以數(shù)個例子，詳細(xì)說明 pdepe 的用法。 uout 為基于所指定位置 xout ，固定