MATLAB Function Reference

pdepe

在一空間變量和時間中，解拋物線與橢圓部份微分等式(parabolic and elliptic partial differential equations，PDEs)系統的初始-臨界值問題(initial-boundary value problems)

Syntax

sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan)
sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options)
sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options,p1,p2...)

Arguments

    

m	與系統對稱性有關的參數。m 可以是 slab = 0，cylindrical = 1 或是 spherical = 2。
pdefun	定義 PDE 元素的函數。
icfun	定義初始狀態的函數。
bcfun	定義臨界狀態的函數。
xmesh	定義點在哪一個數值解時被 tspan 的值所要求的向量 [x0, x1, ..., xn]。xmesh 的元素必須滿足 x0 < x1 < ... < xn。xmesh 的長度必須  3.
tspan	定義點在哪一個數值解時被 xmesh 的值所要求的向量 [t0, t1, ..., tf]。tspan 的元素必須滿足 t0 < t1 < ... < tf。tspan 的長度必須 3。
options	一些下面 ODE 解題器的選項，可在 pdepe 獲得：RelTol、AbsTol、NormControl、InitialStep 和 MaxStep。在大部分的情形，這些選項的預設值提供符合要求的解。詳細情形請看 odeset。
p1,p2,...	傳到 pdefun、icfun 和 bcfun 的選擇性參數。

Description

sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan) 在一空間變量 和時間 中，解拋物線與橢圓 PDEs 系統的初始-臨界值問題。常微分等式(ODEs)經過積分，在 tspan 時間得到近似解。pdepe 函數在 xmesh 格線上傳回解答的值。

pdepe 可以解下列形式的 PDEs：

(1-1)

PDEs 須保持 和 。在 間必須是有限的。 可以是 0、1 或 2，分別對應到 slab、cylindrical 或 spherical symmetry。如果 ，則 必須 0。

在 等式 1-1 中， 是變量而 是i來源項。對時間部份微分的連結被限制要乘上對角線矩陣 。矩陣的對角線元素不是零就是正數。是零的元素對應到橢圓等式，其餘的對應到拋物線等式至少必須有一個拋物線等式。如果 值是格線點， 對應到拋物線的的元素可能在孤立值的 消失。不連續的 和/或 會導致容許的物質界面提供在每個界面的格線點。

對 和所有 ，解答元素需滿足下面形式的初始狀態

(1-2)

對所有 ，不論是 或是 ，解答元素滿足下面形式的初始狀態

(1-3)

的元素不是零就是永不為零。注意，臨界狀態是由變量 表示，而不是 。除此之外，兩個係數裡只有 可以與 相依。

在呼叫 sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan) 時：

• m 對應到 。
• xmesh(1) 和 xmesh(end) 分別對應到 和 。
• tspan(1) 和 tspan(end) 分別對應到 和 。
• pdefun 計算 、 和 (等式 1-1)。它有下列的形式

[c,f,s] = pdefun(x,t,u,dudx)

輸入參數有純量 x 和 t 及向量 u 和 dudx，它們分別求取解答 和解答對 部份積分的近似值。c、f 和 s 都是直行向量。c 儲存矩陣 的對角線元素(等式 1-1)。

• icfun 求取初始狀態。它可以是下面的形式

u = icfun(x)

當以參數 x、icfun 呼叫時，計算然後傳回在直行向量 u 中 x 的解答元素初始值。

• bcfun 求取臨界狀態的 和 (等式 1-3)。它有下面的形式

[pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t)

ul 是在左臨界 xl = 的近似值而 ur 是在右臨界 xr = 的近似值。pl 和 ql 是直行向量，它們對應到在 xl 求得的 和 ，同樣的，pr 和 qr 對應到 xr。當 和 時，靠近 的解答界線要求變量 在 消失。pdepe 自動地利用這個臨界狀態且它忽略 pl 和 ql傳回的值。

pdepe 以多維陣列 sol 的形式傳回解答。 = ui = sol(:,:,i) 是解答向量 第 i 個元素的近似值元素 ui(j,k) = sol(j,k,i) 在  = (tspan(j),xmesh(k)) 近似 。

ui = sol(j,:,i) 在時間 tspan(j)、格線點 xmesh(:) 近似於解答的 i 元素 用 pdeval 來計算近似值和它不在 xmesh 點的部份微分。詳細情形請看 pdeval。

sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options) 像之前一樣解答，不過用 odeset 函數產生的參數 options 值取代預設值。只有基本 ODE 解題器的部份選項可用於 pdepe：RelTol、AbsTol、NormControl、InitialStep 和 MaxStep。移走輸入參數 options 的預設情形通常可以得到令人滿意的結果。詳細情形請看 odeset。

sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options,p1,p2...) 傳另外的參數 p1、p2、... 到函數 pdefun、icfun 和 bcfun。如果沒有選項被設定，用 options = []。

Remarks

• 陣列 xmesh 和 tspan 在 pdepe 中扮演不同的角色。

tspan - pdepe 函數用動態選擇時間步距和公式的 ODE 解題器做時間的積分。tspan 的元素只定義在你想要有解答的地方，且計算花費和 tspan 只有微弱的相關。

xmesh - 在 xmesh 格線解答的二次近似值。一般來說，最好在解答變化快速的地方用近的空間格線。pdepe 沒有 自動選擇 中的格線。你必須在 xmesh 提供近似的固定格線。計算花費和 xmesh 有很大的關係。當 時，不需要用靠近 的格線評估協調單一性。

• 時間積分由 ode15s 完成。 pdepe 利用ode15s 的能力來解當 等式 1-1 含有橢圓等式時出現的微分-代數等式，還有解含有稀疏模式的 Jacobians。
• 在判斷之後，橢圓等式變成代數等式。如果對應到橢圓等式的初始狀態向量元素判斷後不是"前後一致的"，pdepe 會試著在時間積分之前調整它們。因為這個理由，傳回的初始時間解和其他時間相比可能會有判斷錯誤。如果格線很好，pdepe 可以找到接近給定值的前後一致初始狀態。如果 pdepe 顯示一個訊息表示找到前後一致初始狀態有困難，試著改進格線。

對應到拋物線等式的初始狀態向量元素沒必要調整。

Examples

Example 1. 這個例子示範明確的公式化、計算和圖示單一 ODE 的解。

這個等式在 和時間 時成立。

PDE 滿足初始狀態

和臨界狀態

我們可以很方便的使用副函數來安置所有需要被 pdepe 用到的函數在一個 M-檔。

function pdex1

m = 0;
x = linspace(0,1,20);
t = linspace(0,2,5);

sol = pdepe(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t);
% Extract the first solution component as u.
u = sol(:,:,1);

% A surface plot is often a good way to study a solution.
surf(x,t,u)    
title('Numerical solution computed with 20 mesh points.')
xlabel('Distance x')
ylabel('Time t')

% A solution profile can also be illuminating.
figure
plot(x,u(end,:))
title('Solution at t = 2')
xlabel('Distance x')
ylabel('u(x,2)')
% --------------------------------------------------------------
function [c,f,s] = pdex1pde(x,t,u,DuDx)
c = pi^2;
f = DuDx;
s = 0;
% --------------------------------------------------------------
function u0 = pdex1ic(x)
u0 = sin(pi*x);
% --------------------------------------------------------------
function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t)
pl = ul;
ql = 0;
pr = pi * exp(-t);
qr = 1;

在這個例子中，PDE、初始狀態和臨界狀態都以副函數 pdex1pde、pdex1ic 和 pdex1bc 描述。

表面圖顯示解答的行為。

以下的圖顯示解答在 t 最終值的側面(i.e., t = 2)。

Example 2. 這個例子示範 PDEs 系統的解法。問題在端點和之間都有臨界層。對小的 ，解答變化迅速。

PDEs 如下

其中

.

等式在 和時間 時成立。

PDE 滿足初始狀態

和臨界狀態

在 pdepe 期待的形式下，等式變成

的部份微分臨界狀態必須以變量的形式描述。在 pdepe 期待的形式下，左臨界狀態變成

而右臨界狀態變成

對小的 ，解答變化迅速。程式選擇時間步距來解決迅速的改變，但是如果要看這個行為的圖，範例必須選擇輸出時間。因為有在 [0,1] 兩端點的解答的臨界層，所以範例在靠近 0 和 1 的地方設置格線點，來解決迅速變化的問題。通常，需要一些實驗來選擇一個可以顯露出解答行為的格線。

function pdex4
m = 0;
x = [0 0.005 0.01 0.05 0.1 0.2 0.5 0.7 0.9 0.95 0.99 0.995 1];
t = [0 0.005 0.01 0.05 0.1 0.5 1 1.5 2];

sol = pdepe(m,@pdex4pde,@pdex4ic,@pdex4bc,x,t);
u1 = sol(:,:,1);
u2 = sol(:,:,2);

figure
surf(x,t,u1)
title('u1(x,t)')
xlabel('Distance x')
ylabel('Time t')

figure
surf(x,t,u2)
title('u2(x,t)')
xlabel('Distance x')
ylabel('Time t')
% --------------------------------------------------------------
function [c,f,s] = pdex4pde(x,t,u,DuDx)
c = [1; 1];                                  
f = [0.024; 0.17] .* DuDx;                   
y = u(1) - u(2);
F = exp(5.73*y)-exp(-11.47*y);
s = [-F; F];                                                          
% --------------------------------------------------------------
function u0 = pdex4ic(x);
u0 = [1; 0];                                 
% --------------------------------------------------------------
function [pl,ql,pr,qr] = pdex4bc(xl,ul,xr,ur,t)
pl = [0; ul(2)];                               
ql = [1; 0];                                  
pr = [ur(1)-1; 0];                            
qr = [0; 1];                                  

在這個例子中，PDEs、初始狀態和臨界狀態被副函數 pdex4pde、 pdex4ic 和 pdex4bc 描述。

表面圖顯示解答元素的行為。

See Also

function_handle, pdeval, ode15s, odeset, odeget

References

[1]  Skeel, R. D. and M. Berzins, "A Method for the Spatial Discretization of Parabolic Equations in One Space Variable," SIAM Journal on Scientific and Statistical Computing, Vol. 11, 1990, pp.1-32.

pcolor

pdeval