求解偏微分方程
在一个偏微分方程(PDE),被解的函数依赖于几个变量,微分方程可以包括对每个变量的偏导数。偏微分方程用于模拟波动、热流、流体弥散和其他具有随时间变化的空间行为的现象。
什么类型的偏微分方程可以解决MATLAB?
MATLAB®PDE解算器pdepe
在一个空间变量中求解偏微分方程系统的初边值问题x和时间t.你可以把这些看作是一个变量的ode,它也随着时间的变化而变化。
pdepe
对所解的一维方程使用非正式的分类:
有时间导数的方程抛物线.一个例子是热方程 .
没有时间导数的方程椭圆.一个例子是拉普拉斯方程 .
pdepe
要求方程组中至少有一个抛物方程。换句话说,系统中至少有一个方程必须包含时间导数。
pdepe
也解决了某些2-D和3-D问题,这些问题由于角对称而简化为1-D问题(参见参数描述)对称不变米
更多信息)。
偏微分方程工具箱™将此函数扩展到具有狄利克雷和诺伊曼边界条件的二维和三维的广义问题。
求解一维偏微分方程
1-D PDE包含一个函数u(x,t)这取决于时间t还有一个空间变量x.MATLAB的PDE求解器pdepe
求解该形式的一维抛物型和椭圆型偏微分方程的方程组
该方程具有如下性质:
偏微分方程为t0≤t≤tf而且一个≤x≤b.
空间间隔[一个,b]必须是有限的。
米
取值为0、1或2,对应于板,圆柱,或球形分别对称。如果米> 0,然后一个≥0还必须坚持。的系数 是通量项和吗 是一个源项。
通量项必须依赖于偏导数∂u/∂x.
偏导数对时间的耦合只限于与对角矩阵相乘 .这个矩阵的对角元素要么为零,要么为正。一个为零的元素对应一个椭圆方程,任何其他元素对应一个抛物方程。至少有一个抛物方程。的一个因素c对应的抛物方程可以在的孤立值处消失x如果它们是网格点(对解决方案进行评估的点)。的不连续性c而且年代由于材料接口是允许的,只要在每个接口上放置一个网格点。
解决方案的过程
来解偏微分方程pdepe
,你必须定义方程系数c,f,年代,初始条件,解在边界处的行为,以及用来评价解的网格点。函数调用Sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan)
使用此信息计算指定网格上的解:
在一起,xmesh
而且tspan
向量形成二维网格pdepe
上计算解决方案。
方程
您必须用期望的标准形式表示pdepdepe
.用这种形式,你可以读出系数的值c,f,年代.
在MATLAB中,你可以用这种形式的函数对方程进行编码
函数[c,f,s] = pdefun(x,t,u,dudx) c = 1;F = dudx;S = 0;结束
pdefun
定义方程
.如果有多个方程,那么c,f,年代是每个元素对应一个方程的向量。
初始条件
在初始时间t=t0,对所有人来说x,解分量满足形式的初始条件
在MATLAB中,你可以用函数的形式对初始条件进行编码
函数U0 = icfun(x) U0 = 1;结束
U0 = 1
的初始条件u0(x,t0) = 1.如果有多个方程,那么情况
是一个向量,其中每个元素定义一个方程的初始条件。
边界条件
在边界处x=一个或x=b,对所有人来说t,解分量满足形式的边界条件
问(x,t)是一个对角矩阵,其元素要么为零,要么从不为零。注意,边界条件是用通量表示的f,而不是的偏导数u关于x.这两个系数p(x,t,u)而且问(x,t),只有p可以依赖u.
在MATLAB中,可以用函数形式对边界条件进行编码
函数[pL,qL,pR,qR] = bcfun(xL,uL,xR,uR,t) pL = uL;qL = 0;pR = uR - 1;qR = 0;结束
pL
而且qL
左边界的系数是,而公关
而且qR
是右边边界的系数。在这种情况下bcfun
定义边界条件
如果有多个方程,那么输出pL
,qL
,公关
,qR
是由每个元素定义一个方程的边界条件的向量。
集成选项
选择MATLAB PDE求解器中的默认积分属性来处理常见问题。在某些情况下,可以通过覆盖这些默认值来提高求解器性能。要做到这一点,使用odeset
创建一个选项
结构。然后,将结构传递给pdepe
作为最后一个输入参数:
Sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options)
底层ODE求解器的选项ode15s
,只有下表所示的选项可供选择pdepe
.
类别 |
选项名称 |
---|---|
错误控制 |
|
步长 |
|
事件日志记录 |
评估解决方案
当你解一个方程pdepe
, MATLAB以3-D数组的形式返回解索尔
,在那里索尔(i, j, k)
包含了k
求值的解的分量t(我)
而且x (j)
.一般来说,您可以提取k
解决方案组件的命令U = sol(:,:,k)
.
您指定的时间网格仅用于输出目的,并不影响求解器所采取的内部时间步骤。但是,您指定的空间网格可能会影响解决方案的质量和速度。解出一个方程后,你可以用pdeval
来计算返回的解决方案结构pdepe
用不同的空间网格。
例子:热方程
抛物型偏微分方程的一个例子是一维的热方程:
这个方程描述了热的耗散 而且 .目标是解出温度 .温度最初是一个非零常数,所以初始条件是
同时,温度在左边边界为零,在右边边界非零,所以边界条件是
要在MATLAB中求解这个方程,需要对方程、初始条件和边界条件进行编码,然后在调用求解器之前选择一个合适的解网格pdepe
.您可以将所需的函数作为本地函数包含在文件的末尾(如本例所示),或者将它们作为单独的命名文件保存在MATLAB路径的目录中。
代码方程
在编写方程式之前,您需要确保它的形式是pdepe
解算器预计:
在这种形式下,热方程是
所以系数值为:
的价值
作为参数传递给pdepe
,而其他系数则编码在方程的函数中,即
函数[c,f,s] = heatpde(x,t,u,dudx) c = 1;F = dudx;S = 0;结束
(注意:在示例的最后,所有函数都作为局部函数包含在内。)
代码初始条件
热方程的初始条件函数为 .此函数必须接受for的输入 ,即使它没有使用。
函数U0 =热(x) U0 = 0.5;结束
代码边界条件
边界条件的标准形式pdepe
解算器是
写出这种形式,这个问题的边界条件是
所以 而且 是
.
对应的函数是
函数[pl,ql,pr,qr] = heatbc(xl,ul,xr,ur,t) pl = ul;Ql = 0;Pr = ur - 1;Qr = 0;结束
选择解决方案网格
使用20点的空间网格和30点的时间网格。由于解迅速达到稳态,时间点接近 更紧密地间隔在一起,以捕获输出中的此行为。
L = 1;x = linspace(0,L,20);T = [linspace(0,0.05,20), linspace(0.5,5,10)];
解决方程
最后,利用对称性解方程 , PDE方程,初始条件,边界条件,网格 而且 .
M = 0;Sol = pdepe(m,@heatpde,@heatic,@heatbc,x,t);
策划解决方案
使用显示亮度图像
使解矩阵形象化。
colormap热Imagesc (x,t,sol)“距离x”,“翻译”,“乳胶”) ylabel (“t”,“翻译”,“乳胶”)标题(“$0 \le x \le 1$和$0 \le t \le 5$的热方程”,“翻译”,“乳胶”)
本地函数
函数[c,f,s] = heatpde(x,t,u,dudx) c = 1;F = dudx;S = 0;结束函数U0 =热(x) U0 = 0.5;结束函数[pl,ql,pr,qr] = heatbc(xl,ul,xr,ur,t) pl = ul;Ql = 0;Pr = ur - 1;Qr = 0;结束
PDE示例和文件
对于大多数常见的1-D PDE问题,有几个可用的示例文件可以作为很好的起点。要探索和运行示例,使用微分方程示例应用程序。要运行该应用程序,输入
odeexamples
要打开单个文件进行编辑,请键入
编辑exampleFileName.m
要运行一个示例,输入
exampleFileName
该表包含可用的PDE示例文件的列表。
示例文件 |
描述 |
例子链接 |
---|---|---|
|
简单的PDE,说明了解的公式、计算和绘图。 |
|
|
涉及不连续的问题。 |
|
|
需要计算偏导数值的问题。 |
|
|
两个偏微分方程的解在区间两端都有边界层,且小时变化迅速的系统t. |
|
|
阶跃函数为初始条件的偏微分方程系统。 |
参考文献
[1] Skeel, R. D.和M. Berzins,“抛物方程在一个空间变量上的空间离散化方法”,科学与统计计算SIAM杂志, 1990年第11卷,第1-32页。