MatlabPDE工具箱有限元法求解偏微分方程教学提纲
Matlab求解微分方程(组)及偏微分方程(组)
第四讲 Matlab 求解微分方程（组）理论介绍：Matlab 求解微分方程（组）命令 求解实例：Matlab 求解微分方程（组）实例实际应用问题通过数学建模所归纳得到的方程，绝大多数都是微分方程，真正能得到代数方程的机会很少.另一方面，能够求解的微分方程也是十分有限的，特别是高阶方程和偏微分方程（组）.这就要求我们必须研究微分方程（组）的解法：解析解法和数值解法. 一．相关函数、命令及简介1.在Matlab 中，用大写字母D 表示导数，Dy 表示y 关于自变量的一阶导数，D2y 表示y 关于自变量的二阶导数，依此类推.函数dsolve 用来解决常微分方程（组）的求解问题，调用格式为：X=dsolve(‘eqn1’,’eqn2’,…)函数dsolve 用来解符号常微分方程、方程组，如果没有初始条件，则求出通解，如果有初始条件，则求出特解.注意，系统缺省的自变量为t2.函数dsolve 求解的是常微分方程的精确解法，也称为常微分方程的符号解.但是，有大量的常微分方程虽然从理论上讲，其解是存在的，但我们却无法求出其解析解，此时，我们需要寻求方程的数值解，在求常微分方程数值解方面，MATLAB 具有丰富的函数，我们将其统称为solver ，其一般格式为：[T,Y]=solver(odefun,tspan,y0)说明：(1)solver 为命令ode45、ode23、ode113、ode15s 、ode23s 、ode23t 、ode23tb 、ode15i 之一.(2)odefun 是显示微分方程'(,)y f t y =在积分区间tspan 0[,]f t t =上从0t 到ft 用初始条件0y 求解.(3)如果要获得微分方程问题在其他指定时间点012,,,,f t t t t 上的解，则令tspan 012[,,,]f t t t t =(要求是单调的).(4)因为没有一种算法可以有效的解决所有的ODE 问题，为此，Matlab 提供了多种求解器solver ，对于不同的ODE 问题，采用不同的solver.表1 Matlab中文本文件读写函数说明：ode23、ode45是极其常用的用来求解非刚性的标准形式的一阶微分方程（组）的初值问题的解的Matlab常用程序，其中：ode23采用龙格-库塔2阶算法，用3阶公式作误差估计来调节步长，具有低等的精度.ode45则采用龙格-库塔4阶算法，用5阶公式作误差估计来调节步长，具有中等的精度.3．在matlab命令窗口、程序或函数中创建局部函数时，可用内联函数inline，inline函数形式相当于编写M函数文件，但不需编写M-文件就可以描述出某种数学关系.调用inline函数，只能由一个matlab表达式组成，并且只能返回一个变量，不允许[u,v]这种向量形式.因而，任何要求逻辑运算或乘法运算以求得最终结果的场合，都不能应用inline函数，inline函数的一般形式为：FunctionName=inline(‘函数内容’, ‘所有自变量列表’)例如：（求解F(x)=x^2*cos(a*x)-b ,a,b是标量；x是向量）在命令窗口输入:Fofx=inline(‘x .^2*cos(a*x)-b ’ , ‘x ’,’a ’,’b ’); g= Fofx([pi/3 pi/3.5],4,1) 系统输出为：g=-1.5483 -1.7259注意：由于使用内联对象函数inline 不需要另外建立m 文件，所有使用比较方便，另外在使用ode45函数的时候，定义函数往往需要编辑一个m 文件来单独定义，这样不便于管理文件，这里可以使用inline 来定义函数. 二．实例介绍1.几个可以直接用Matlab 求微分方程精确解的实例 例1 求解微分方程2'2x y xy xe -+=程序：syms x y; y=dsolve(‘Dy+2*x*y=x*exp(-x^2)’,’x ’)例 2 求微分方程'0x xy y e +-=在初始条件(1)2y e =下的特解并画出解函数的图形.程序：syms x y; y=dsolve(‘x*Dy+y-exp(1)=0’,’y(1)=2*exp(1)’,’x ’);ezplot(y)例 3 求解微分方程组530tdx x y e dtdy x y dt⎧++=⎪⎪⎨⎪--=⎪⎩在初始条件00|1,|0t t x y ====下的特解并画出解函数的图形.程序：syms x y t[x,y]=dsolve('Dx+5*x+y=exp(t)','Dy-x-3*y=0','x(0)=1','y(0)=0','t')simple(x); simple(y)ezplot(x,y,[0,1.3]);axis auto2.用ode23、ode45等求解非刚性标准形式的一阶微分方程（组）的初值问题的数值解（近似解）例 4 求解微分方程初值问题2222(0)1dy y x xdx y ⎧=-++⎪⎨⎪=⎩的数值解，求解范围为区间[0,0.5].程序：fun=inline('-2*y+2*x^2+2*x','x','y'); [x,y]=ode23(fun,[0,0.5],1); plot(x,y,'o-')例 5 求解微分方程22'2(1)0,(0)1,(0)0d y dyy y y y dt dtμ--+===的解，并画出解的图形.分析：这是一个二阶非线性方程，我们可以通过变换，将二阶方程化为一阶方程组求解.令12,,7dyx y x dtμ===，则 121221212,(0)17(1),(0)0dx x x dtdx x x x x dt⎧==⎪⎪⎨⎪=--=⎪⎩ 编写M-文件vdp.m function fy=vdp(t,x)fy=[x(2);7*(1-x(1)^2)*x(2)-x(1)]; end在Matlab 命令窗口编写程序 y0=[1;0][t,x]=ode45(@vdp,[0,40],y0);或[t,x]=ode45('vdp',[0,40],y0); y=x(:,1);dy=x(:,2); plot(t,y,t,dy)练习与思考：M-文件vdp.m 改写成inline 函数程序？ 3.用Euler 折线法求解Euler 折线法求解的基本思想是将微分方程初值问题00(,)()dyf x y dxy x y ⎧=⎪⎨⎪=⎩ 化成一个代数(差分)方程，主要步骤是用差商()()y x h y x h +-替代微商dydx，于是00()()(,())()k k k k y x h y x f x y x h y y x +-⎧=⎪⎨⎪=⎩记1,(),k k k k x x h y y x +=+=从而1(),k k y y x h +=+于是0011(),,0,1,2,,1(,).k k k k k k y y x x x h k n y y hf x y ++=⎧⎪=+=-⎨⎪=+⎩例 6 用Euler 折线法求解微分方程初值问题22(0)1dyx y dxy y ⎧=+⎪⎨⎪=⎩的数值解（步长h 取0.4），求解范围为区间[0,2].分析：本问题的差分方程为00110,1,0.4,0,1,2,,1(,).k k k k k k x y h x x h k n y y hf x y ++===⎧⎪=+=-⎨⎪=+⎩程序：>> clear >> f=sym('y+2*x/y^2'); >> a=0; >> b=2; >> h=0.4; >> n=(b-a)/h+1; >> x=0; >> y=1;>> szj=[x,y];%数值解 >> for i=1:n-1y=y+h*subs(f,{'x','y'},{x,y});%subs ，替换函数 x=x+h;szj=[szj;x,y]; end>>szj>> plot(szj(:,1),szj(:,2))说明：替换函数subs 例如：输入subs(a+b,a,4) 意思就是把a 用4替换掉，返回 4+b ，也可以替换多个变量，例如：subs(cos(a)+sin(b),{a,b},[sym('alpha'),2])分别用字符alpha 替换a 和2替换b ，返回 cos(alpha)+sin(2)特别说明：本问题可进一步利用四阶Runge-Kutta 法求解，Euler 折线法实际上就是一阶Runge-Kutta 法，Runge-Kutta 法的迭代公式为001112341213243(),,(22),6(,),0,1,2,,1(,),22(,),22(,).k k k k k k k k k k k k y y x x x h h y y L L L L L f x y k n h h L f x y L h h L f x y L L f x h y hL ++=⎧⎪=+⎪⎪=++++⎪⎪=⎪=-⎨⎪=++⎪⎪⎪=++⎪⎪=++⎩相应的Matlab 程序为：>> clear >> f=sym('y+2*x/y^2'); >> a=0; >> b=2; >> h=0.4; >> n=(b-a)/h+1; >> x=0; >> y=1;>> szj=[x,y];%数值解 >> for i=1:n-1l1=subs(f, {'x','y'},{x,y});替换函数 l2=subs(f, {'x','y'},{x+h/2,y+l1*h/2}); l3=subs(f, {'x','y'},{x+h/2,y+l2*h/2}); l4=subs(f, {'x','y'},{x+h,y+l3*h});y=y+h*(l1+2*l2+2*l3+l4)/6; x=x+h;szj=[szj;x,y]; end >>szj>> plot(szj(:,1),szj(:,2))练习与思考：(1)ode45求解问题并比较差异. (2)利用Matlab 求微分方程(4)(3)''20y y y -+=的解.(3)求解微分方程''2',2(1)0,030,(0)1,(0)0y y y y x y y --+=≤≤==的特解. (4)利用Matlab 求微分方程初值问题2''''00(1)2,|1,|3x x x y xy y y ==+===的解. 提醒：尽可能多的考虑解法 三．微分方程转换为一阶显式微分方程组Matlab 微分方程解算器只能求解标准形式的一阶显式微分方程（组）问题，因此在使用ODE 解算器之前，我们需要做的第一步，也是最重要的一步就是借助状态变量将微分方程（组）化成Matlab 可接受的标准形式.当然，如果ODEs 由一个或多个高阶微分方程给出，则我们应先将它变换成一阶显式常微分方程组.下面我们以两个高阶微分方程组构成的ODEs 为例介绍如何将它变换成一个一阶显式微分方程组.Step 1 将微分方程的最高阶变量移到等式左边，其它移到右边，并按阶次从低到高排列.形式为：()'''(1)'''(1)()'''(1)'''(1)(,,,,,,,,,,)(,,,,,,,,,,)m m n n m n x f t x x x x y y y y y g t x x x x y y y y ----⎧=⎨=⎩Step 2 为每一阶微分式选择状态变量，最高阶除外'''(1)123'''(1)123,,,,,,,,,m m n m m m m n x x x x x x x x x y x y x y x y--++++========注意：ODEs 中所有是因变量的最高阶次之和就是需要的状态变量的个数，最高阶的微分式不需要给它状态变量.Step 3 根据选用的状态变量，写出所有状态变量的一阶微分表达式''''122334123''12123,,,,(,,,,,),,(,,,,,)m m n m m m nm n x x x x x x x f t x x x x xx xg t x x x x +++++======练习与思考：(1)求解微分方程组**'''3312*'''3312()()22x x x y x r r y y y x y r r μμμμμμ⎧+-=+--⎪⎪⎨⎪=+--⎪⎩其中2r =1r =*1,μμ=-1/82.45,μ=(0) 1.2,x =(0)0,y ='(0)0,x ='(0) 1.049355751y =-(2)求解隐式微分方程组''''''''''''2235x y x y x y x y xy y ⎧+=⎨++-=⎩ 提示：使用符号计算函数solve 求'''',x y ，然后利用求解微分方程的方法 四．偏微分方程解法Matlab 提供了两种方法解决PDE 问题，一是使用pdepe 函数，它可以求解一般的PDEs,具有较大的通用性，但只支持命令形式调用；二是使用PDE 工具箱，可以求解特殊PDE 问题，PDEtoll 有较大的局限性，比如只能求解二阶PDE 问题，并且不能解决片微分方程组，但是它提供了GUI 界面，从复杂的编程中解脱出来，同时还可以通过File —>Save As 直接生成M 代码.1.一般偏微分方程（组）的求解(1)Matlab 提供的pdepe 函数，可以直接求解一般偏微分方程（组），它的调用格式为：sol=pdepe(m,@pdefun,@pdeic,@pdebc,x,t)@pdefun 是PDE 的问题描述函数，它必须换成标准形式：(,,)[(,,,)](,,,)m m u u u uc x t x x f x t u s x t u x t x x x-∂∂∂∂∂=+∂∂∂∂∂ 这样，PDE 就可以编写入口函数：[c,f,s]=pdefun(x,t,u,du)，m,x,t 对应于式中相关参数，du 是u 的一阶导数，由给定的输入变量可表示出c,f,s 这三个函数.@pdebc 是PDE 的边界条件描述函数，它必须化为形式：(,,)(,,).*(,,,)0up x t u q x t u f x t u x∂==∂ 于是边值条件可以编写函数描述为：[pa,qa,pb,qb]=pdebc(x,t,u,du),其中a 表示下边界，b 表示上边界.@pdeic 是PDE 的初值条件，必须化为形式：00(,)u x t u =，故可以使用函数描述为：u0=pdeic(x)sol 是一个三维数组，sol(:,:,i)表示i u 的解，换句话说，k u 对应x(i)和t(j)时的解为sol(i,j,k)，通过sol ，我们可以使用pdeval 函数直接计算某个点的函数值.(2)实例说明 求解偏微分2111222221220.024()0.17()u u F u u t xu u F u u tx ⎧∂∂=--⎪⎪∂∂⎨∂∂⎪=+-⎪∂∂⎩ 其中， 5.7311.46()xx F x e e -=-且满足初始条件12(,0)1,(,0)0u x u x ==及边界条件1(0,)0,u t x ∂=∂221(0,)0,(1,)1,(1,)0uu t u t t x∂===∂ 解：(1)对照给出的偏微分方程和pdepe 函数求解的标准形式，原方程改写为111221220.024()1.*()10.17u u F u u x u F u u u t x x ∂⎡⎤⎢⎥--⎡⎤⎡⎤⎡⎤∂∂∂=+⎢⎥⎢⎥⎢⎥⎢⎥-∂∂∂⎣⎦⎣⎦⎣⎦⎢⎥⎢⎥∂⎣⎦可见1121220.024()10,,,()10.17u F u u x m c f s F u u u x ∂⎡⎤⎢⎥--⎡⎤⎡⎤∂====⎢⎥⎢⎥⎢⎥-∂⎣⎦⎣⎦⎢⎥⎢⎥∂⎣⎦%目标PDE 函数function [c,f,s]=pdefun(x,t,u,du) c=[1;1];f=[0.024*du(1);0.17*du(2)];temp=u(1)-u(2);s=[-1;1].*(exp(5.73*temp)-exp(-11.46*temp)) end(2)边界条件改写为：下边界2010.*00f u ⎡⎤⎡⎤⎡⎤+=⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦上边界1110.*000u f -⎡⎤⎡⎤⎡⎤+=⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦%边界条件函数function [pa,qa,pb,qb]=pdebc(xa,ua,xb,ub,t) pa=[0;ua(2)]; qa=[1;0]; pb=[ub(1)-1;0]; qb=[0;1]; end(3)初值条件改写为：1210u u ⎡⎤⎡⎤=⎢⎥⎢⎥⎣⎦⎣⎦%初值条件函数 function u0=pdeic(x) u0=[1;0]; end(4)编写主调函数 clcx=0:0.05:1; t=0:0.05:2; m=0;sol=pdepe(m,@pdefun,@pdeic,@pdebc,x,t); subplot(2,1,1) surf(x,t,sol(:,:,1)) subplot(2,1,2) surf(x,t,sol(:,:,2))练习与思考： This example illustrates the straightforward formulation, computation, and plotting of the solution of a single PDE.2()u u t x xπ∂∂∂=∂∂∂ This equation holds on an interval 01x ≤≤ for times 0t ≥. The PDE satisfies the initial condition (,0)sin u x x π= and boundary conditions(0,)0;(1,)0t u u t e t xπ-∂=+=∂ 2.PDEtool 求解偏微分方程 (1)PDEtool （GUI ）求解偏微分方程的一般步骤在Matlab 命令窗口输入pdetool ，回车，PDE 工具箱的图形用户界面(GUI)系统就启动了.从定义一个偏微分方程问题到完成解偏微分方程的定解，整个过程大致可以分为六个阶段Step 1 “Draw 模式”绘制平面有界区域Ω，通过公式把Matlab 系统提供的实体模型：矩形、圆、椭圆和多边形，组合起来，生成需要的平面区域.Step 2 “Boundary 模式”定义边界，声明不同边界段的边界条件.Step 3 “PDE 模式”定义偏微分方程，确定方程类型和方程系数c,a,f,d ，根据具体情况，还可以在不同子区域声明不同系数.Step 4 “Mesh 模式”网格化区域Ω，可以控制自动生成网格的参数，对生成的网格进行多次细化，使网格分割更细更合理.Step 5 “Solve 模式”解偏微分方程，对于椭圆型方程可以激活并控制非线性自适应解题器来处理非线性方程；对于抛物线型方程和双曲型方程，设置初始边界条件后可以求出给定时刻t 的解；对于特征值问题，可以求出给定区间上的特征值.求解完成后，可以返回到Step 4，对网格进一步细化，进行再次求解.Step 6 “View 模式”计算结果的可视化，可以通过设置系统提供的对话框，显示所求的解的表面图、网格图、等高线图和箭头梯形图.对于抛物线型和双曲线型问题的解还可以进行动画演示.(2)实例说明用法求解一个正方形区域上的特征值问题：12|0u u u u λ∂Ω⎧-∆-=⎪⎨⎪=⎩ 正方形区域为：11,1 1.x x -≤≤-≤≤(1)使用PDE工具箱打开GUI求解方程(2)进入Draw模式，绘制一个矩形，然后双击矩形，在弹出的对话框中设置Left=-1,Bottom=-1,Width=2,Height=2,确认并关闭对话框(3)进入Boundary模式，边界条件采用Dirichlet条件的默认值(4)进入PDE模式，单击工具栏PDE按钮，在弹出的对话框中方程类型选择Eigenmodes,参数设置c=1,a=-1/2,d=1,确认后关闭对话框(5)单击工具栏的 按钮，对正方形区域进行初始网格剖分，然后再对网格进一步细化剖分一次(6)点开solve菜单，单击Parameters选项，在弹出的对话框中设置特征值区域为[-20,20](7)单击Plot菜单的Parameters项，在弹出的对话框中选中Color、Height(3-D plot)和show mesh项，然后单击Done确认(8)单击工具栏的“=”按钮，开始求解。
(完整版)偏微分方程的MATLAB解法
引言偏微分方程定解问题有着广泛的应用背景。
人们用偏微分方程来描述、解释或者预见各种自然现象，并用于科学和工程技术的各个领域fll。
然而，对于广大应用工作者来说，从偏微分方程模型出发，使用有限元法或有限差分法求解都要耗费很大的工作量，才能得到数值解。
现在，MATLAB PDEToolbox已实现对于空间二维问题高速、准确的求解过程。
偏微分方程如果一个微分方程中出现的未知函数只含一个自变量，这个方程叫做常微分方程，也简称微分方程；如果一个微分方程中出现多元函数的偏导数，或者说如果未知函数和几个变量有关，而且方程中出现未知函数对几个变量的导数，那么这种微分方程就是偏微分方程。
常用的方法有变分法和有限差分法。
变分法是把定解问题转化成变分问题，再求变分问题的近似解；有限差分法是把定解问题转化成代数方程，然后用计算机进行计算；还有一种更有意义的模拟法，它用另一个物理的问题实验研究来代替所研究某个物理问题的定解。
虽然物理现象本质不同，但是抽象地表示在数学上是同一个定解问题，如研究某个不规则形状的物体里的稳定温度分布问题，由于求解比较困难，可作相应的静电场或稳恒电流场实验研究，测定场中各处的电势，从而也解决了所研究的稳定温度场中的温度分布问题。
随着物理科学所研究的现象在广度和深度两方面的扩展，偏微分方程的应用范围更广泛。
从数学自身的角度看，偏微分方程的求解促使数学在函数论、变分法、级数展开、常微分方程、代数、微分几何等各方面进行发展。
从这个角度说，偏微分方程变成了数学的中心。
一、MATLAB方法简介及应用1.1 MATLAB简介MATLAB是美国MathWorks公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和Simulink两大部分。
1.2 Matlab主要功能数值分析数值和符号计算工程与科学绘图控制系统的设计与仿真数字图像处理数字信号处理通讯系统设计与仿真财务与金融工程1.3 优势特点1) 高效的数值计算及符号计算功能，能使用户从繁杂的数学运算分析中解脱出来；2) 具有完备的图形处理功能，实现计算结果和编程的可视化；3) 友好的用户界面及接近数学表达式的自然化语言，使学者易于学习和掌握；4) 功能丰富的应用工具箱(如信号处理工具箱、通信工具箱等) ，为用户提供了大量方便实用的处理工具。
(完整版)偏微分方程的MATLAB解法
引言偏微分方程定解问题有着广泛的应用背景。
人们用偏微分方程来描述、解释或者预见各种自然现象，并用于科学和工程技术的各个领域fll。
然而，对于广大应用工作者来说，从偏微分方程模型出发，使用有限元法或有限差分法求解都要耗费很大的工作量，才能得到数值解。
现在，MATLAB PDEToolbox已实现对于空间二维问题高速、准确的求解过程。
偏微分方程如果一个微分方程中出现的未知函数只含一个自变量，这个方程叫做常微分方程，也简称微分方程；如果一个微分方程中出现多元函数的偏导数，或者说如果未知函数和几个变量有关，而且方程中出现未知函数对几个变量的导数，那么这种微分方程就是偏微分方程。
常用的方法有变分法和有限差分法。
变分法是把定解问题转化成变分问题，再求变分问题的近似解；有限差分法是把定解问题转化成代数方程，然后用计算机进行计算；还有一种更有意义的模拟法，它用另一个物理的问题实验研究来代替所研究某个物理问题的定解。
虽然物理现象本质不同，但是抽象地表示在数学上是同一个定解问题，如研究某个不规则形状的物体里的稳定温度分布问题，由于求解比较困难，可作相应的静电场或稳恒电流场实验研究，测定场中各处的电势，从而也解决了所研究的稳定温度场中的温度分布问题。
随着物理科学所研究的现象在广度和深度两方面的扩展，偏微分方程的应用范围更广泛。
从数学自身的角度看，偏微分方程的求解促使数学在函数论、变分法、级数展开、常微分方程、代数、微分几何等各方面进行发展。
从这个角度说，偏微分方程变成了数学的中心。
一、MATLAB方法简介及应用1.1 MATLAB简介MATLAB是美国MathWorks公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和Simulink两大部分。
1.2 Matlab主要功能数值分析数值和符号计算工程与科学绘图控制系统的设计与仿真数字图像处理数字信号处理通讯系统设计与仿真财务与金融工程1.3 优势特点1) 高效的数值计算及符号计算功能，能使用户从繁杂的数学运算分析中解脱出来；2) 具有完备的图形处理功能，实现计算结果和编程的可视化；3) 友好的用户界面及接近数学表达式的自然化语言，使学者易于学习和掌握；4) 功能丰富的应用工具箱(如信号处理工具箱、通信工具箱等) ，为用户提供了大量方便实用的处理工具。
《应用matlab的pdetoolbox求解偏微分方程综合实验》
《应用matlab的pdetoolbox求解偏微分方程综合实验》编辑整理：尊敬的读者朋友们：这里是精品文档编辑中心，本文档内容是由我和我的同事精心编辑整理后发布的，发布之前我们对文中内容进行仔细校对，但是难免会有疏漏的地方，但是任然希望（《应用matlab的pdetoolbox求解偏微分方程综合实验》）的内容能够给您的工作和学习带来便利。
同时也真诚的希望收到您的建议和反馈，这将是我们进步的源泉，前进的动力。
本文可编辑可修改，如果觉得对您有帮助请收藏以便随时查阅，最后祝您生活愉快业绩进步，以下为《应用matlab的pdetoolbox求解偏微分方程综合实验》的全部内容。
2011届信计专业学生综合实验题目（要求按照所附开题报告表的格式填写提交开题报告。
一个小组选做一题，小组全体成员共同完成，每个小组只提交一份实验报告,按照出力多少排名。
提交时间在本学期18周以前.）3 应用Matlab 的PDE Toolbox 求解偏微分方程熟悉Matlab 的PDE 工具箱的功能，并用其求解具有工程背景的偏微分方程，要求分别对三种类型方程:抛物型、椭圆形和双曲型。
阐述清楚如下方法：1、这里我们先脱离问题所含有的工程背景，分别举三个例子，大致的说明一下如何使用Matlab 的PDE 工具箱来求解三种类型的偏微分方程。
(1)椭圆形方程： 考虑一块圆形金属片,中心挖去一正方形，外边界满足Neumann 条件,内边界满足Dirichlet条件:考虑到入射波以方向，所以上式可以写成这样得到求解这个入射波的定解问题：这里取波长为0。
1。
现在用GUI 来求解上述这个问题，并最终获得其解得图形:图1 初始网格、加密以及网格剖分数据图 2解的三维图形Ω=+∇⋅∇-in f au u c ,)(.),(xik e y x v r π-=-=x -,ikx e r -=),(y x r ⎪⎪⎩⎪⎪⎨⎧-=-=∂∂-=-==+∆--外边界上内边界上在区域内,60,,0602i ik n ree r r k r xi ikx ,60=kλ图3 解的二维动画图附录1：二维动画的Matlab程序(2)抛物型方程：考虑一个圆柱形放射性杆，其左端供热，右端保持常温，侧面与环境有热交换。
偏微分方程的matlab解法讲课文档
现在二十五页，总共二十五页。
Байду номын сангаас
初始条件
现在三页，总共二十五页。
先确定方程大类
现在四页，总共二十五页。
Draw Mode
画图模式，先将处理的区域画出来，二维， 方形，圆形，支持多边形，可以手动更改坐
标，旋转rotate
例如，对于细杆导热，虽然是一维问题， 可以将宽度y虚拟出来，对应于y的边界条 件和初始条件按照题意制定
现在五页，总共二十五页。
°C，板的右边热量从板向环境空气定常流动，其他边及内孔边界保
持绝缘。初始
题;
t t 是板的温度为0 °C ，于是概括为如下定解问 0
d u u0 , t
u 1 0 0 ,在 左 边 界 上
u 1， 在 右 边 界 上 n u = 0， 其 他 边 界 上 n
u t to 0
区域的边界顶点坐标为(-0.5,-0.8), (0.5,-0.8), (-0.5,0.8), (0.5,0.8)。 内边界顶点坐标(-0.05,-0.4), (-0.05,0.4) ,(0.05,-0.4), (0.05,0.4)。
MATLAB形成M文件．
现在十六页，总共二十五页。
现在十七页，总共二十五页。
第四步：设置方程类型
选择PDE菜单中PDE Mode命令，进入PDE模式， 再单击PDE菜单中PDE Secification选项，打开
PDE Secification对话框，设置方程类型．
本例取抛物型方程 du(cu)auf, t
Conditions选项，打开Boundary Conditions对话框， 输入边界件．本例取默认条件，即将全部边界设为齐次
Dirichlet条件，边界显示为红色． 如果想将几何与边界信息存储，可选Boundary菜单 中的Export Decomposed Geometrv．Boundary
MATLAB求解PDE问题
MATLAB求解PDE问题(1）-—概述、例子(转)（2011—07-20 16:48:45）MATLAB PDE T oolbox提供利用有限元方法求解偏微分方程的GUI以及相应的命令行函数。
利用该工具箱可以求解椭圆型方程、抛物型方程、双曲型方程、特征值方程以及非线性方程.PDE Toolbox的功能非常强大，网上有许多利用PDE Toolbox解决各种物理问题的论文，还有专门介绍工具箱的参考书。
网上的例子虽然很多,但是大部分是介绍PDE工具箱自带的一些例子，这些例子中解的区域，边界条件是PDE工具箱已经编写好的，直接调用就可以。
对于该如何自己设定求解区域及边界条件,却很少有人涉及.网上搜索发现只有刘平在博客中详细介绍过求解区域的设定。
下面以一个椭圆型方程的例子来详细说明求解的各个步骤，希望对大家能有所帮助。
设要求如下形式的椭圆方程的解:按照PDE的要求,将方程化为标准形式求解后的图像如下,第一幅图是解的图像,第二幅是计算误差。
从第二幅图可以看到，计算的最大误差是10—3方量级。
通过这个例子我们可以基本掌握PDE求解偏微分方程的步骤和方法，后面我将详细介绍如何设置区域及边界条件。
掌握了区域和边界条件的设定,就可以轻松求解遇到的偏微分方程了。
图后是附带的matlab命令以及注释，并提供m文件附件下载,下载后解压即可.希望能对大家有所帮助。
下面是编写的求解上述方程的matlab语句及说明:g='mygeom'；b=’mybound'；定义区域，边界条件。
mygeom是定义区域的子函数名，函数名可根据自己的需要取定,区域的确定规则由pdegeom函数说明，注意pdegeom函数只是说明如何定义区域,它并不直接确定区域；mybound是定义边界条件的子函数名，与区域类似，边界的确定规则由函数pdebound确定。
后面我会详细介绍区域和边界的取法。
[p,e,t］= initmesh（g）;网格初始化，此处也可以写成[p,e，t］= initmesh（’mygeom')；这样可以省略上面的语句［p,e，t］= refinemesh（g,p，e，t）;[p，e，t］= refinemesh(g,p,e,t）；加密网格两次，需要加密几次重复几次即可，根据具体问题确定加密次数U= assempde(b,p,e，t，1，0，'2*(x+y)-4'）；调用assempde函数计算方程的数值解，assempde函数的详细用法可以参考MATH 网站或者PDE的使用指南。
Matlab偏微分方程求解方法
Matlab 偏微分方程求解方法目录：§1 Function Summary on page 10-87§2 Initial Value Problems on page 10-88§3 PDE Solver on page 10-89§4 Integrator Options on page 10-92§5 Examples” on page 10-93§1 Function Summary1.1 PDE Solver” on page 10-871,2 PDE Helper Functi on” on page 10-871.3 PDE SolverThis is the MATLAB PDE solver.PDE Helper Function§2 Initial Value Problemspdepe solves systems of parabolic and elliptic PDEs in one spatial variable x and time t, of the form)xu ,u ,t ,x (s ))x u ,u ,t ,x (f x (x x t u )x u ,u ,t ,x (c m m ∂∂+∂∂∂∂=∂∂∂∂- (10-2) The PDEs hold for b x a ,t t t f 0≤≤≤≤.The interval [a, b] must be finite. mcan be 0, 1, or 2, corresponding to slab, cylindrical, or spherical symmetry,respectively. If m > 0, thena ≥0 must also hold.In Equation 10-2,)x /u ,u ,t ,x (f ∂∂ is a flux term and )x /u ,u ,t ,x (s ∂∂ is a source term. The flux term must depend on x /u ∂∂. The coupling of the partial derivatives with respect to time is restricted to multiplication by a diagonal matrix )x /u ,u ,t ,x (c ∂∂. The diagonal elements of this matrix are either identically zero or positive. An element that is identically zero corresponds to an elliptic equation and otherwise to a parabolic equation. There must be at least one parabolic equation. An element of c that corresponds to a parabolic equation can vanish at isolated values of x if they are mesh points.Discontinuities in c and/or s due to material interfaces are permitted provided that a mesh point is placed at each interface.At the initial time t = t0, for all x the solution components satisfy initial conditions of the form)x (u )t ,x (u 00= (10-3)At the boundary x = a or x = b, for all t the solution components satisfy a boundary condition of the form0)xu ,u ,t ,x (f )t ,x (q )u ,t ,x (p =∂∂+ (10-4) q(x, t) is a diagonal matrix with elements that are either identically zero or never zero. Note that the boundary conditions are expressed in terms of the f rather than partial derivative of u with respect to x-x /u ∂∂. Also, ofthe two coefficients, only p can depend on u.§3 PDE Solver3.1 The PDE SolverThe MATLAB PDE solver, pdepe, solves initial-boundary value problems for systems of parabolic and elliptic PDEs in the one space variable x and time t.There must be at least one parabolic equation in the system.The pdepe solver converts the PDEs to ODEs using a second-order accurate spatial discretization based on a fixed set of user-specified nodes. The discretization method is described in [9]. The time integration is done with ode15s. The pdepe solver exploits the capabilities of ode15s for solving the differential-algebraic equations that arise when Equation 10-2 contains elliptic equations, and for handling Jacobians with a specified sparsity pattern. ode15s changes both the time step and the formula dynamically.After discretization, elliptic equations give rise to algebraic equations. If the elements of the initial conditions vector that correspond to elliptic equations are not “consistent” with the discretization, pdepe tries to adjust them before eginning the time integration. For this reason, the solution returned for the initial time may have a discretization error comparable to that at any other time. If the mesh is sufficiently fine, pdepe can find consistent initial conditions close to the given ones. If pdepe displays amessage that it has difficulty finding consistent initial conditions, try refining the mesh. No adjustment is necessary for elements of the initial conditions vector that correspond to parabolic equations.PDE Solver SyntaxThe basic syntax of the solver is:sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan)Note Correspondences given are to terms used in “Initial Value Problems” on page 10-88.The input arguments arem: Specifies the symmetry of the problem. m can be 0 =slab, 1 = cylindrical, or 2 = spherical. It corresponds to m in Equation 10-2. pdefun: Function that defines the components of the PDE. Itcomputes the terms f,c and s in Equation 10-2, and has the form[c,f,s] = pdefun(x,t,u,dudx)where x and t are scalars, and u and dudx are vectors that approximate the solution and its partial derivative with respect to . c, f, and s are column vectors. c stores the diagonal elements of the matrix .icfun: Function that evaluates the initial conditions. It has the formu = icfun(x)When called with an argument x, icfun evaluates and returns the initial values of the solution components at x in the column vector u.bcfun:Function that evaluates the terms and of the boundary conditions. Ithas the form[pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t)where ul is the approximate solution at the left boundary xl = a and ur is the approximate solution at the right boundary xr = b. pl and ql are column vectors corresponding to p and the diagonal of q evaluated at xl. Similarly, pr and qr correspond to xr. When m>0 and a = 0, boundedness of the solution near x = 0 requires that the f vanish at a = 0. pdepe imposes this boundary condition automatically and it ignores values returned in pl and ql.xmesh:Vector [x0, x1, ..., xn] specifying the points at which a numerical solution is requested for every value in tspan. x0 and xn correspond to a and b , respectively. Second-order approximation to the solution is made on the mesh specified in xmesh. Generally, it is best to use closely spaced mesh points where the solution changes rapidly. pdepe does not select the mesh in automatically. You must provide an appropriate fixed mesh in xmesh. The cost depends strongly on the length of xmesh. When , it is not necessary to use a fine mesh near to x=0 account for the coordinate singularity.The elements of xmesh must satisfy x0 < x1 < ... < xn.The length of xmesh must be ≥3.tspan:Vector [t0, t1, ..., tf] specifying the points at which a solution is requested for every value in xmesh. t0 and tf correspond tot and f t,respectively.pdepe performs the time integration with an ODE solver that selects both the time step and formula dynamically. The solutions at the points specified in tspan are obtained using the natural continuous extension of the integration formulas. The elements of tspan merely specify where you want answers and the cost depends weakly on the length of tspan.The elements of tspan must satisfy t0 < t1 < ... < tf.The length of tspan must be ≥3.The output argument sol is a three-dimensional array, such that•sol(:,:,k) approximates component k of the solution .•sol(i,:,k) approximates component k of the solution at time tspan(i) and mesh points xmesh(:).•sol(i,j,k) approximates component k of the solution at time tspan(i) and the mesh point xmesh(j).4.2 PDE Solver OptionsFor more advanced applications, you can also specify as input arguments solver options and additional parameters that are passed to the PDE functions.options:Structure of optional parameters that change the default integration properties. This is the seventh input argument.sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options)See “Integrator Options” on page 10-92 for more information.Integrator OptionsThe default integration properties in the MATLAB PDE solver are selected to handle common problems. In some cases, you can improve solver performance by overriding these defaults. You do this by supplying pdepe with one or more property values in an options structure.sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options)Use odeset to create the options structure. Only those options of the underlying ODE solver shown in the following table are available for pdepe.The defaults obtained by leaving off the input argument options are generally satisfactory. “Integrator Options” on page 10-9 tells you how to create the structure and describes the properties.PDE Properties§4 Examples•“Single PDE” on page 10-93•“System of PDEs” on page 10-98•“Additional Examples” on page 10-1031.Single PDE• “Solving the Equation” on page 10-93• “Evaluating the Solution” on page 10-98Solving the Equation. This example illustrates the straightforward formulation, solution, and plotting of the solution of a single PDE222x u t u ∂∂=∂∂π This equation holds on an interval 1x 0≤≤ for times t ≥ 0. At 0t = the solution satisfies the initial condition x sin )0,x (u π=.At 0x =and 1x = , the solution satisfies the boundary conditions0)t ,1(xu e ,0)t ,0(u t =∂∂+π=- Note The demo pdex1 contains the complete code for this example. The demo uses subfunctions to place all functions it requires in a single MATLAB file.To run the demo type pdex1 at the command line. See “PDE Solver Syntax” on page 10-89 for more information. 1 Rewrite the PDE. Write the PDE in the form)xu ,u ,t ,x (s ))x u ,u ,t ,x (f x (x x t u )x u ,u ,t ,x (c m m ∂∂+∂∂∂∂=∂∂∂∂- This is the form shown in Equation 10-2 and expected by pdepe. For this example, the resulting equation is0x u x x x t u 002+⎪⎭⎫ ⎝⎛∂∂∂∂=∂∂π with parameter and the terms 0m = and the term0s ,xu f ,c 2=∂∂=π=2 Code the PDE. Once you rewrite the PDE in the form shown above (Equation 10-2) and identify the terms, you can code the PDE in a function that pdepe can use. The function must be of the form[c,f,s] = pdefun(x,t,u,dudx)where c, f, and s correspond to the f ,c and s terms. The code below computes c, f, and s for the example problem.function [c,f,s] = pdex1pde(x,t,u,DuDx)c = pi^2;f = DuDx;s = 0;3 Code the initial conditions function. You must code the initial conditions in a function of the formu = icfun(x)The code below represents the initial conditions in the function pdex1ic. Partial Differential Equationsfunction u0 = pdex1ic(x)u0 = sin(pi*x);4 Code the boundary conditions function. You must also code the boundary conditions in a function of the form[pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t) The boundary conditions, written in the same form as Equation 10-4, are0x ,0)t ,0(x u .0)t ,0(u ==∂∂+and1x ,0)t ,1(xu .1e t ==∂∂+π- The code below evaluates the components )u ,t ,x (p and )u ,t ,x (q of the boundary conditions in the function pdex1bc.function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t)pl = ul;ql = 0;pr = pi * exp(-t);qr = 1;In the function pdex1bc, pl and ql correspond to the left boundary conditions (x=0 ), and pr and qr correspond to the right boundary condition (x=1).5 Select mesh points for the solution. Before you use the MATLAB PDE solver, you need to specify the mesh points at which you want pdepe to evaluate the solution. Specify the points as vectors t and x.The vectors t and x play different roles in the solver (see “PDE Solver” on page 10-89). In particular, the cost and the accuracy of the solution depend strongly on the length of the vector x. However, the computation is much less sensitive to the values in the vector t.10 CalculusThis example requests the solution on the mesh produced by 20 equally spaced points from the spatial interval [0,1] and five values of t from thetime interval [0,2].x = linspace(0,1,20);t = linspace(0,2,5);6 Apply the PDE solver. The example calls pdepe with m = 0, the functions pdex1pde, pdex1ic, and pdex1bc, and the mesh defined by x and t at which pdepe is to evaluate the solution. The pdepe function returns the numerical solution in a three-dimensional array sol, wheresol(i,j,k) approximates the kth component of the solution,u, evaluated atkt(i) and x(j).m = 0;sol = pdepe(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t);This example uses @ to pass pdex1pde, pdex1ic, and pdex1bc as function handles to pdepe.Note See the function_handle (@), func2str, and str2func reference pages, and the @ section of MATLAB Programming Fundamentals for information about function handles.7 View the results. Complete the example by displaying the results:a Extract and display the first solution component. In this example, the solution has only one component, but for illustrative purposes, the example “extracts” it from the three-dimensional array. The surface plot shows the behavior of the solution.u = sol(:,:,1);surf(x,t,u)title('Numerical solution computed with 20 mesh points') xlabel('Distance x') ylabel('Time t')Distance xNumerical solution computed with 20 mesh points.Time tb Display a solution profile at f t , the final value of . In this example,2t t f ==.figure plot(x,u(end,:)) title('Solution at t = 2') xlabel('Distance x') ylabel('u(x,2)')Solutions at t = 2.Distance xu (x ,2)Evaluating the Solution. After obtaining and plotting the solution above, you might be interested in a solution profile for a particular value of t, or the time changes of the solution at a particular point x. The kth column u(:,k)(of the solution extracted in step 7) contains the time history of the solution at x(k). The jth row u(j,:) contains the solution profile at t(j). Using the vectors x and u(j,:), and the helper function pdeval, you can evaluate the solution u and its derivative at any set of points xout [uout,DuoutDx] = pdeval(m,x,u(j,:),xout)The example pdex3 uses pdeval to evaluate the derivative of the solution at xout = 0. See pdeval for details.2. System of PDEsThis example illustrates the solution of a system of partial differential equations. The problem is taken from electrodynamics. It has boundary layers at both ends of the interval, and the solution changes rapidly for small . The PDEs are)u u (F xu017.0t u )u u (F xu 024.0t u 212222212121-+∂∂=∂∂--∂∂=∂∂ where )y 46.11exp()y 73.5exp()y (F --=. The equations hold on an interval1x 0≤≤ for times 0t ≥.The solution satisfies the initial conditions0)0,x (u ,1)0,x (u 21≡≡and boundary conditions0)t ,1(xu,0)t ,1(u ,0)t ,0(u ,0)t ,0(x u 2121=∂∂===∂∂ Note The demo pdex4 contains the complete code for this example. The demo uses subfunctions to place all required functions in a single MATLAB file. To run this example type pdex4 at the command line.1 Rewrite the PDE. In the form expected by pdepe, the equations are⎥⎦⎤⎢⎣⎡---+⎥⎦⎤⎢⎣⎡∂∂∂∂∂∂=⎥⎦⎤⎢⎣⎡∂∂⎥⎦⎤⎢⎣⎡)u u (F )u u (F )x /u 170.0)x /u (024.0x u u t *.1121212121 The boundary conditions on the partial derivatives of have to be written in terms of the flux. In the form expected by pdepe, the left boundary condition is⎥⎦⎤⎢⎣⎡=⎥⎦⎤⎢⎣⎡∂∂∂∂⎥⎦⎤⎢⎣⎡+⎥⎦⎤⎢⎣⎡00)x /u (170.0)x /u (024.0*.01u 0212and the right boundary condition is⎥⎦⎤⎢⎣⎡=⎥⎦⎤⎢⎣⎡∂∂∂∂⎥⎦⎤⎢⎣⎡+⎥⎦⎤⎢⎣⎡-00)x /u (170.0)x /u (024.0*.1001u 2112 Code the PDE. After you rewrite the PDE in the form shown above, you can code it as a function that pdepe can use. The function must be of the form[c,f,s] = pdefun(x,t,u,dudx)where c, f, and s correspond to the , , and terms in Equation 10-2. 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];3 Code the initial conditions function. The initial conditions function must be of the formu = icfun(x)The code below represents the initial conditions in the function pdex4ic. function u0 = pdex4ic(x);u0 = [1; 0];4 Code the boundary conditions function. The boundary conditions functions must be of the form[pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t)The code below evaluates the components p(x,t,u) and q(x,t) (Equation 10-4) of the boundary conditions in the function pdex4bc.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];5 Select mesh points for the solution. The solution changes rapidly for small t . The program selects the step size in time to resolve this sharp change, but to see this behavior in the plots, output times must be selected accordingly. There are boundary layers in the solution at both ends of [0,1], so mesh points must be placed there to resolve these sharp changes. Often some experimentation is needed to select the mesh that reveals the behavior of the solution.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];6 Apply the PDE solver. The example calls pdepe with m = 0, the functions pdex4pde, pdex4ic, and pdex4bc, and the mesh defined by x and t at which pdepe is to evaluate the solution. The pdepe function returns the numerical solution in a three-dimensional array sol, wheresol(i,j,k) approximates the kth component of the solution, μk, evaluated at t(i) and x(j).m = 0;sol = pdepe(m,@pdex4pde,@pdex4ic,@pdex4bc,x,t);7 View the results. The surface plots show the behavior of the solution components. u1 = sol(:,:,1); u2 = sol(:,:,2); figure surf(x,t,u1) title('u1(x,t)') xlabel('Distance x') 其输出图形为Distance xu1(x,t)Time tfigure surf(x,t,u2) title('u2(x,t)') xlabel('Distance x') ylabel('Time t')Distance xu2(x,t)Time tAdditional ExamplesThe following additional examples are available. Type edit examplename to view the code and examplename to run the example.。
matlab偏微分方程工具箱使用手册
MATLAB偏微分方程工具箱使用手册一、Matlab偏微分方程工具箱介绍Matlab偏微分方程工具箱是Matlab中用于求解偏微分方程（PDE）问题的工具。
它提供了一系列函数和工具，可以用于建立、求解和分析PDE问题。
PDE是许多科学和工程领域中的重要数学模型，包括热传导、扩散、波动等现象的数值模拟、分析和优化。
Matlab偏微分方程工具箱为用户提供了丰富的功能和灵活的接口，使得PDE问题的求解变得更加简单和高效。
二、使用手册1. 安装和启用在开始使用Matlab偏微分方程工具箱前，首先需要确保Matlab已经安装并且包含了PDE工具箱。
确认工具箱已经安装后，可以通过以下命令启用PDE工具箱：```pdetool```这将打开PDE工具箱的图形用户界面，用户可以通过该界面进行PDE 问题的建立、求解和分析。
2. PDE建模在PDE工具箱中，用户可以通过几何建模工具进行PDE问题的建立。
用户可以定义几何形状、边界条件、初值条件等，并选择适当的PDE方程进行描述。
PDE工具箱提供了各种几何建模和PDE方程描述的选项，用户可以根据实际问题进行相应的设置和定义。
3. 求解和分析一旦PDE问题建立完成，用户可以通过PDE工具箱提供的求解器进行求解。
PDE工具箱提供了各种数值求解方法，包括有限元法、有限差分法等。
用户可以选择适当的求解方法，并进行求解。
求解完成后，PDE工具箱还提供了丰富的分析功能，用户可以对结果进行后处理、可视化和分析。
4. 结果导出和应用用户可以将求解结果导出到Matlab环境中，并进行后续的数据处理、可视化和分析。
用户也可以将结果导出到其他软件环境中进行更进一步的处理和应用。
三、个人观点和理解Matlab偏微分方程工具箱是一个非常强大的工具，它为科学和工程领域中的PDE问题提供了简单、高效的解决方案。
通过使用PDE工具箱，用户可以快速建立、求解和分析复杂的PDE问题，从而加快科学研究和工程设计的进程。
利用Matlab中的PDE工具箱进行偏微分方程求解！
利用Matlab中的PDE工具箱进行偏微分方程求解！在科学技术各领域中，有很多问题都可以归结为偏微分方程问题。
在物理专业的力学、热学、电学、光学、近代物理课程中都可遇见偏微分方程。
偏微分方程，再加上边界条件、初始条件构成的数学模型，只有在很特殊情况下才可求得解析解。
随着计算机技术的发展，采用数值计算方法，可以得到其数值解。
下面的几个简单例子，将为大家介绍如何利用Matlab中的PDE 工具箱进行偏微分方程的求解！抛物线型受热金属块的热传导问题：一块受热的有矩形裂纹的金属块。
金属块的左侧被加热到100℃，右侧热量则以恒定速率降低到周围空气的温度，所有其他边界都是独立的，于是边界条件为：初始温度为0℃。
指定起始时间为0，研究5s的热扩散问题。
求解结果如下：具体步骤：1.建模。
先画一个起点为(-0.5,-0.8)，终点为(0.5,0.8)的矩形，再建另一个矩形，起点(-0.05,-0.4)，终点(0.05,0.4)，然后在设置公式中输入R1-R2；2.设置边界条件；3.PDE中设置参数，d=1,c=1,a=0,f=0；4.初始化网格并细化一次；5.在Solve Parameter中，设置u0=0，时间为[0:0.5:5]，然后求解。
注意到金属块的温度升高很快，可试着改变时间列表的表达式为logspace(-2,0.5,20)以便观察。
双曲线型方形薄膜横向波动问题（波动方程）：薄膜左侧和右侧固定（u=0），上端和下端自由：满足边界条件的初值为：求解结果如下：具体步骤：1.建模，画出起点(-1,-1)，终点(1,1)的矩形；2.确定边界条件；3.PDE Specification对话框中选择双曲线型，参数设为c=1,a=0,f=0,d=1；4.初始化网格，并细化一次；5.在Solve Parameters对话框的时间中输入linspace(0,5,31)，u的初值为atan(cos(pi/2*x))，一次偏导的初值为3*sin(pi*x).*exp(sin(pi/2*y))，然后求解。
偏微分方程的MATLAB解法
偏微分⽅程的MATLAB解法引⾔偏微分⽅程定解问题有着⼴泛的应⽤背景。
⼈们⽤偏微分⽅程来描述、解释或者预见各种⾃然现象，并⽤于科学和⼯程技术的各个领域fll。
然⽽，对于⼴⼤应⽤⼯作者来说，从偏微分⽅程模型出发，使⽤有限元法或有限差分法求解都要耗费很⼤的⼯作量，才能得到数值解。
现在，MATLAB PDEToolbox已实现对于空间⼆维问题⾼速、准确的求解过程。
偏微分⽅程如果⼀个微分⽅程中出现的未知函数只含⼀个⾃变量，这个⽅程叫做常微分⽅程，也简称微分⽅程；如果⼀个微分⽅程中出现多元函数的偏导数，或者说如果未知函数和⼏个变量有关，⽽且⽅程中出现未知函数对⼏个变量的导数，那么这种微分⽅程就是偏微分⽅程。
常⽤的⽅法有变分法和有限差分法。
变分法是把定解问题转化成变分问题，再求变分问题的近似解；有限差分法是把定解问题转化成代数⽅程，然后⽤计算机进⾏计算；还有⼀种更有意义的模拟法，它⽤另⼀个物理的问题实验研究来代替所研究某个物理问题的定解。
虽然物理现象本质不同，但是抽象地表⽰在数学上是同⼀个定解问题，如研究某个不规则形状的物体⾥的稳定温度分布问题，由于求解⽐较困难，可作相应的静电场或稳恒电流场实验研究，测定场中各处的电势，从⽽也解决了所研究的稳定温度场中的温度分布问题。
随着物理科学所研究的现象在⼴度和深度两⽅⾯的扩展，偏微分⽅程的应⽤范围更⼴泛。
从数学⾃⾝的⾓度看，偏微分⽅程的求解促使数学在函数论、变分法、级数展开、常微分⽅程、代数、微分⼏何等各⽅⾯进⾏发展。
从这个⾓度说，偏微分⽅程变成了数学的中⼼。
⼀、MATLAB⽅法简介及应⽤1.1 MATLAB简介MATLAB是美国MathWorks公司出品的商业数学软件，⽤于算法开发、数据可视化、数据分析以及数值计算的⾼级技术计算语⾔和交互式环境，主要包括MATLAB和Simulink两⼤部分。
1.2 Matlab主要功能数值分析数值和符号计算⼯程与科学绘图控制系统的设计与仿真数字图像处理数字信号处理通讯系统设计与仿真财务与⾦融⼯程1.3 优势特点1) ⾼效的数值计算及符号计算功能，能使⽤户从繁杂的数学运算分析中解脱出来；2) 具有完备的图形处理功能，实现计算结果和编程的可视化；3) 友好的⽤户界⾯及接近数学表达式的⾃然化语⾔，使学者易于学习和掌握；4) 功能丰富的应⽤⼯具箱(如信号处理⼯具箱、通信⼯具箱等) ，为⽤户提供了⼤量⽅便实⽤的处理⼯具。
MATLAB偏微分方程求解课件
常微分方程：在微分方程中，若自变量的个数只有一个的微分方程。
偏微分方程：自变量的个数有两个或两个以上的微分方程。
求解偏微分方程的方法
• 求解偏微分方程的数值方法： • 1. 有限元法（Finite Element Method, FEM）---
hp-FEM • 2. 有限体积法（Finite Volume Method, FVM） • 3. 有限差分法（Finite Difference Method, FDM）。
• 其它：广义有限元法（Generalized Finite Element Method, FFEM）、扩展有限元法（eXtended Finite Element Method, XFEM）、无网格有限元法（Meshfree Finite Element Method）、离散迦辽金有限元法（Discontinuous Galerkin Finite Element Method, DGFEM）等。
step3:边界条件和初值条件
• 初值条件可以通过【Solve】->【 Parameters…】设置
• 边值条件设置如下 • (1)点击工具栏的第6 个按钮【区域边界】，
显示如下
• (2)【Boundary】->【Remove All Subdomain Borders】移除所有子域的边界， 将得到所有子域合并成一个求解域
是name后边的名字'NumberTitle','off'是关掉默认显示名字。 • subplot(211) • surf(x,t,sol(:,:,1))%sol(:,:,i)表示ui的解 • title('The Solution of u_1') • xlabel('X') • ylabel('T') • zlabel('U') • subplot(212) • surf(x,t,sol(:,:,2))%sol(:,:,i)表示ui的解 • title('The Solution of u_2') • xlabel('X') • ylabel('T') • zlabel('U')
偏微分数值解(2,MATLAB求解方法)
初始条件:
u1 ( x,0) 1,
u2 ( x,0) 0
边界条件:
u1 (0, t ) 0, u1 (1, t ) 1 x u2 u2 (0, t ) 0, (1, t ) 0 x
方程来自电动力学中关于电磁场理论的一个 偏微分方程组。
2.1 用偏微分方程工具箱求解微分方程
直接使用图形用户界面（ Graphical User Interface，简记作GUI）求解．
图22.1 所讨论定解问题的区域
第三步：选取边界 首先选择Boundary菜单中Boundary Mode命 令，进入边界模式．然后单击Boundary菜单中 Remove All Subdomain Borders选项，从而去掉子 域边界，如图22.2．单击Boundary菜单中Specify Boundary Conditions选项，打开Boundary Conditions对话框，输入边界条件．本例取默认条 件，即将全部边界设为齐次Dirichlet条件，边界显 示为红色．如果想将几何与边界信息存储，可选择 Boundary菜单中的Export Decomposed Geometry,Boundary Cond‟s命令，将它们分别存储 在g、b变量中，并通过MATLAB形成M文件．
第八步：若要画等值线图和矢量场图，单击 Plot 菜单中 Parameter 选项，在 Plot selection 对话框中选中 Contour 和 Arrows 两项．然后单击 Plot 按钮，可显示解的等值 线图和矢量场图，如图 2. 6 所示。
图 2.6 解的等值线图和矢量场图
(1) u1=hyperbolic（u0，ut0，tlist，b，p，e，t，c，a，f，d）
《应用matlab的pdetoolbox求解偏微分方程综合实验》
pderect([0 1 0 1],'SQ1')
③单击Boundary菜单下的Boundary Mode选项，进入边界模式.选择Boundary菜单中的Show Subdomain Labels命令，显示区域编号(如图10所示)，再单击所要删去的线段或弧段，执行Boundary菜单中Remove Subdomain Border命令，则可以得到图11.
%取复数解的实部
h=newplot;set(get(h,'Parent'),'Renderer','zbuffer')
pdeplot(p,e,t,'xydata',real(u),'zdata',real(u),...
'mesh','off');
colormap(cool)
clc
%制作反射波的动画程序
②在PDE Toolbox面板中依次拖出圆 双击所绘制的圆，将圆心设置为(0,0),半径分别为1,0.8,0.6,0.5,0.4；依次拖出矩形 双击所绘制的矩形，分别设置Left为-0.2、-0.1，Bottom分别为0.2、0.1,Width分别为0.2、0.2，Height分别为0.9、0.9.这样便可以得到图9所示的平面图.
for j=1:m,...
uu=real(exp(-j*2*pi/m*sqrt(-1))*u);...
fprintf('%d',j);...
pdeplot(p,e,t,'xydata',uu,'colorbar','off',...
matlab用有限元法求解偏微分方程组
matlab用有限元法求解偏微分方程组使用有限元法求解偏微分方程组是一种常见的数值计算方法，它在工程领域和科学研究中广泛应用。
本文将介绍如何利用MATLAB软件进行有限元法求解偏微分方程组的基本步骤和注意事项。
我们需要了解有限元法的基本原理。
有限元法是一种将连续问题离散化为有限个小区域，通过在每个小区域内建立适当的数学模型，然后将这些小区域连接起来形成整个问题的数学模型的方法。
在有限元法中，我们通常将问题的域分割成许多小的有限元，每个有限元都具有简单的几何形状，如线段、三角形或四边形。
然后，在每个有限元上建立适当的近似函数，通过对这些函数的系数进行求解，我们可以得到问题的近似解。
在MATLAB中，有限元法的求解过程可以分为以下几个步骤：1. 离散化域：根据问题的几何形状，将问题的域进行离散化处理。
离散化可以采用三角剖分法或四边形剖分法，将域分割成许多小的有限元。
2. 建立数学模型：在每个有限元上建立适当的数学模型。
这通常涉及选择适当的近似函数，并在每个有限元上求解这些函数的系数。
3. 组装方程：将每个有限元上的数学模型组装成整个问题的数学模型。
这涉及到将有限元之间的边界条件进行匹配，并建立整个问题的刚度矩阵和载荷向量。
4. 求解方程：利用线性代数求解方法，求解得到问题的近似解。
MATLAB提供了各种求解线性方程组的函数，如“\”运算符、LU 分解和共轭梯度法等。
5. 后处理：对求解结果进行后处理，包括绘制解的图形、计算问题的误差等。
在进行有限元法求解偏微分方程组时，需要注意以下几点：1. 网格剖分的合理性：网格剖分的精细程度对结果的精确性有很大影响。
网格过于粗糙可能导致结果的不准确，而网格过于细小则会增加计算的复杂性。
因此，需要根据问题的特点和计算资源的限制选择合适的网格剖分。
2. 近似函数的选择：近似函数的选择直接影响到结果的准确性和计算的效率。
一般情况下，近似函数的阶数越高，结果的准确性越高，但计算的复杂性也越大。
matlab 求解偏微分方程
matlab 求解偏微分方程使用MATLAB求解偏微分方程摘要：偏微分方程(partial differential equation, PDE)是数学中重要的一类方程，广泛应用于物理、工程、经济、生物等领域。
MATLAB 是一种强大的数值计算软件，提供了丰富的工具箱和函数，可以用来求解各种类型的偏微分方程。
本文将介绍如何使用MATLAB来求解偏微分方程，并通过具体案例进行演示。
引言：偏微分方程是描述多变量函数的方程，其中包含了函数的偏导数。
一般来说，偏微分方程可以分为椭圆型方程、双曲型方程和抛物型方程三类。
求解偏微分方程的方法有很多，其中数值方法是最常用的一种。
MATLAB作为一种强大的数值计算软件，提供了丰富的工具箱和函数，可以用来求解各种类型的偏微分方程。
方法：MATLAB提供了多种求解偏微分方程的函数和工具箱，包括pdepe、pdetoolbox和pde模块等。
其中，pdepe函数是用来求解带有初始条件和边界条件的常微分方程组的函数，可以用来求解一维和二维的偏微分方程。
pdepe函数使用有限差分法或有限元法来离散化偏微分方程，然后通过求解离散化后的常微分方程组得到最终的解。
案例演示：考虑一维热传导方程的求解，偏微分方程为：∂u/∂t = α * ∂^2u/∂x^2其中，u(x,t)是温度分布函数，α是热扩散系数。
假设初始条件为u(x,0)=sin(pi*x)，边界条件为u(0,t)=0和u(1,t)=0。
我们需要定义偏微分方程和边界条件。
在MATLAB中，可以使用匿名函数来定义偏微分方程和边界条件。
然后，我们使用pdepe函数求解偏微分方程。
```matlabfunction [c,f,s] = pde(x,t,u,DuDx)c = 1;f = DuDx;s = 0;endfunction u0 = uinitial(x)u0 = sin(pi*x);endfunction [pl,ql,pr,qr] = uboundary(xl,ul,xr,ur,t)pl = ul;ql = 0;pr = ur;qr = 0;endx = linspace(0,1,100);t = linspace(0,0.1,10);m = 0;sol = pdepe(m,@pde,@uinitial,@uboundary,x,t);u = sol(:,:,1);surf(x,t,u);xlabel('Distance x');ylabel('Time t');zlabel('Temperature u');```在上述代码中，我们首先定义了偏微分方程函数pde，其中c、f和s分别表示系数c、f和s。
有限元的MATLAB解法
有限元的MATLAB解法1.打开MATLAB。
2.输入“pdetool”再回车，会跳出PDE Toolbox的窗口（PDE意为偏微分方程，是partial differential equations的缩写），需要的话可点击Options菜单下Grid命令，打开栅格。
3.完成平面几何模型：在PDE Toolbox的窗口中，点击工具栏下的矩形几何模型进行制作模型，可画矩形R，椭圆E，圆C，然后在Set formula栏进行编辑并（如双脊波导R1+R2+R3改为RI-R2-R3,设定a、b、s/a、d/b的值从而方便下步设定坐标）用算术运算符将图形对象名称连接起来，若还需要，可进行储存，形成M文件。
4.用左键双击矩形进行坐标设置：将大的矩形left和bottom都设为0，width是矩形波导的X轴的长度，height是矩形波导的y轴的长度，以大的矩形左下角点为原点坐标为参考设置其他矩形坐标。
5.进行边界设置：点击“Boundary”中的“Boundary Mode”,再点击“Boundary”中的“Specify Boundary Conditions”,选择符合的边界条件，Neumann为诺曼条件,Dirichlet为狄利克雷条件，边界颜色显示为红色。
6.进入PDE模式：点击"PDE"菜单下“PDE Mode”命令，进入PDE 模式，单击“PDE Specification”,设置方程类型，“Elliptic”为椭圆型，“Parabolic”为抛物型，“Hyperbolic”为双曲型，“Eigenmodes”为特征值问题。
7.对模型进行剖分：点击“Mesh”中“Initialize Mesh”进行初次剖分，若要剖的更细，再点击“Refine Mesh”进行网格加密。
8.进行计算：点击“Solve”中“Solve PDE”,解偏微分方程并显示图形解，u值即为Hz或者Ez。
9.单击“Plot”菜单下“Parameters”选项，打开“Plot Selection”对话框。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
M a t l a b PDE工具箱有限元法求解偏微分
方程
在科学技术各领域中，有很多问题都可以归结为偏微分方程问题。
在物理专业的力学、热学、电学、光学、近代物理课程中都可遇见偏微分方程。
偏微分方程，再加上边界条件、初始条件构成的数学模型，只有在很特殊情况下才可求得解析解。
随着计算机技术的发展，采用数值计算方法，可以得到其数值解。
偏微分方程基本形式
而以上的偏微分方程都能利用PDE工具箱求解。
PDE工具箱
PDE工具箱的使用步骤体现了有限元法求解问题的基本思路，包括如下基本步骤：
1) 建立几何模型
2) 定义边界条件
3) 定义PDE类型和PDE系数
4) 三角形网格划分
5) 有限元求解
6) 解的图形表达
以上步骤充分体现在PDE工具箱的菜单栏和工具栏顺序上，如下
具体实现如下。
打开工具箱
输入pdetool可以打开偏微分方程求解工具箱，如下
首先需要选择应用模式，工具箱根据实际问题的不同提供了很多应用模式，用户可以基于适当的模式进行建模和分析。
在Options菜单的Application菜单项下可以做选择，如下
或者直接在工具栏上选择，如下
列表框中各应用模式的意义为：
① Generic Scalar：一般标量模式（为默认选项）。
② Generic System：一般系统模式。
③ Structural Mech.，Plane Stress：结构力学平面应力。
④ Structural Mech.，Plane Strain：结构力学平面应变。
⑤ Electrostatics：静电学。
⑥ Magnetostatics：电磁学。
⑦ Ac Power Electromagnetics：交流电电磁学。
⑧ Conductive Media DC：直流导电介质。
⑨ Heat Tranfer：热传导。
⑩ Diffusion：扩散。
可以根据自己的具体问题做相应的选择，这里要求解偏微分方程，故使用默认值。
此外，
对于其他具体的工程应用模式，此工具箱已经发展到了Comsol Multiphysics软件，它提供了更强大的建模、求解功能。
另外，可以在菜单Options下做一些全局的设置，如下
l Grid：显示网格
l Grid Spacing…：控制网格的显示位置
l Snap：建模时捕捉网格节点，建模时可以打开
l Axes Limits…：设置坐标系范围
l Axes Equal：同Matlab的命令axes equal命令
建立几何模型
使用菜单Draw的命令或使用工具箱命令可以实现简单几何模型的建立，如下
各项代表的意义分别为
l 绘制矩形或方形；
l 绘制同心矩形或方形；
l 绘制椭圆或圆；
l 绘制同心椭圆或圆；
l 绘制多义线。
这里只绘制一个圆如下
定义边界条件
选择Boundary菜单下的Specify Boundary Conditions…，如下
定义PDE类型和PDE系数
选择PDE菜单下的PDE Specifications…，如下
三角形网格划分
选择Mesh菜单下的Initialize Mesh初始化三角形网格，再选择Refine Mesh改进初始网格并细化网格，如下
初始化网格
细化网格
另外还可以进一步选择Jiggle Mesh微调网格。
最后可以选择Display Triangle Quality 显示三角形网格的质量图，其中1表示质量最好，0表示最差，如下
有限元求解
选择Solve菜单下的Solve PDE选项进行PDE问题的求解，如下
解的图形表达
选择Plot菜单下的Parameters…可以设置显示的效果，如下显示结果如下
比较数值解与精确解的误差：可见数值解的精度是很高的。