Matlab求解微分方程组及偏微分方程组Word格式文档下载.docx

1、tspan=t。,t1,t2J|ltf（要求是单调的）.（4）因为没有一种算法可以有效的解决所有的 ODE问题，为此，Matlab提供了多种求解器solver，对于不同的ODE问题，采用不同的solver.表1 Matlab中文本文件读写函数求解器ODE类型特点说明ode45非刚性单步算法：4、5阶Runge-Kutta方程；累计截断误差（Ax）3大部分场合的首选算法ode232、3阶Runge-Kutta使用于精度较低的情形ode113多步法：Adams算法；高低精度可达10“10计算时间比ode45短ode23t适度刚性采用梯形算法适度刚性情形ode15s刚性Gears反向数值微分；精度中

2、等若ode45失效时，可尝试使用ode23s单步法：2阶Rosebrock算法；低精度当精度较低时，计算时间比ode15s短ode23tb梯形算法；低精度ode23、ode45是极其常用的用来求解非刚性的标准形式的一阶微分方程（组）的初值问题的解的Matlab常用程序，其中：ode23采用龙格-库塔2阶算法，用3阶公式作误差估计来调节步长，具有低 等的精度ode45则采用龙格-库塔4阶算法，用5阶公式作误差估计来调节步长，具有 中等的精度3.在matlab命令窗口、程序或函数中创建局部函数时，可用内联函数inline， inline函数形式相当于编写M函数文件，但不需编写 M-文件就可以描述出

3、某种 数学关系 调用inline函数，只能由一个 matlab表达式组成，并且只能返回一个 变量，不允许u,v这种向量形式.因而，任何要求逻辑运算或乘法运算以求得最 终结果的场合，都不能应用inline函数，inline函数的一般形式为：FunctionName=inline（函数内容所有自变量列表例如：（求解F（x）=x2*cos（a*x）-b ,a,b是标量；x是向量）在命令窗口输入：Fofx=inline（ x .A2*cos（a（*X） , x a , b）;g= Fofx（pi/3 pi/3.5,4,1）系统输出为：g=-1.5483 -1.7259注意：由于使用内联对象函数inli

4、ne不需要另外建立m文件，所有使用比较 方便，另外在使用ode45函数的时候，定义函数往往需要编辑一个 m文件来单 独定义，这样不便于管理文件，这里可以使用 inline来定义函数.二实例介绍1.几个可以直接用Matlab求微分方程精确解的实例例1求解微分方程y，2xy=xe程序：syms x y; y=dsolve（ Dy+2*x*y=x*exp（-xA2）例2求微分方程xy y-ex =0在初始条件y=2e下的特解并画出解函数的图形. y=dsolve（ x*Dy+y-exp（1）=0 y（1）=2*exp（1） ezplot（y）dx t5x y = e例3求解微分方程组（dt 在初始条

5、件x|y = 1,y |t曲=0下的特解鱼-x3y=0.dt并画出解函数的图形.syms x y tx,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求解微分方程初值问题 dx = _2y 2 2x的数值解，求解范围为区【y（0）=i间0,0.5.fun=inline（,-2*y+2*xA2+2*x,x,y,）;x,y=ode23（fu n,0

6、,0.5,1）;plot（x,y,o-）例5求解微分方程马 7（1一丫2）巴 y =0, y（0） =1,y（0） =0 的解，并画出dt dt解的图形.分析：这是一个二阶非线性方程，我们可以通过变换，将二阶方程化为一阶 方程组求解令Mg哼7，则捲（0） = 1x2（0） = 0dx!頁dX2 =7（1x2）X2 7, .dt编写M-文件vdp.mfun cti on fy=vdp（t,x） fy=x（2）;7*（1-x（1）A2）*x （2） -x（1）; end在Matlab命令窗口编写程序y0=1;0t,x=ode45（vdp,0,40,y0）;或t,x=ode45（vdp,0,40,y

7、0）; y=x（:,1）;dy=x（:,2）;plot（t,y,t,dy）练习与思考：M-文件vdp.m改写成inline函数程序？3用Euler折线法求解Euler折线法求解的基本思想是将微分方程初值问题f （x, y）.y（xO =y。化成一个代数（差分）方程，主要步骤是用差商y（x h） - y（x）替代微商 鱼，于是 h dx乂二 y（x。记 Xk Xk h, yy（Xk）,从而 yk y（Xk h）,于是y。=y（x）,/ x“=Xk+h, k = 0,1,2j|, n1yk 1 二 yk hf （Xk,yQ.例6用Euler折线法求解微分方程初值问题dy 2xy 2dx y2y（0

8、）=1的数值解（步长h取0.4）,求解范围为区间0,2.本问题的差分方程为k =0,1,2,川，n-1=0, y0 = 1,h = 0.4xk 1 =Xk h,yk 1 二 yk hf （Xk, yj clear f=sym（y+2*x/yA2 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;endszj plot（szj（:,1）,szj（:,2）替换函数subs例如：输入subs（a+b,a,4）意思就是把a

9、用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法的迭代公式为）.Xk 厂Xk h.hyk 1 = yk （Li 2L2 2L3 L4）,6k=0,1,2,|（, n-1Li 二 f （Xk, yk）,h hL2 - f （Xk -,yk -Li）,L f （Xk -,yk - L2）,L4 = f（X

10、k h,yk hL3）.相应的Matlab程序为：y+2*X/yA2l1= subs（f, 替换函数l2=subs（f, ,x+h/2,y+l1*h/2）;l3=subs（f, ,x+h/2,y+l2*h/2）;l4=subs（f, ,x+h,y+l3*h）;y=y+h*（l1+2*l2+2*l3+l4）/6;（1）ode45求解问题并比较差异.利用Matlab求微分方程y-2y y = 0的解.H Q I（3）求解微分方程 y 2（1-y）y y =0,0 乞 x 岂 30,y（0） hy（0） = 0 的特解., , Q H I I（4）利用Matlab求微分方程初值问题（1 x ）y 2

11、xy , y卜卫=1, y 1x=3的解.提醒：尽可能多的考虑解法 三微分方程转换为一阶显式微分方程组Matlab微分方程解算器只能求解标准形式的一阶显式微分方程（组）问题， 因此在使用ODE解算器之前，我们需要做的第一步，也是最重要的一步就是借 助状态变量将微分方程（组）化成 Matlab可接受的标准形式.当然，如果ODEs由一个或多个高阶微分方程给出，则我们应先将它变换成一阶显式常微分方程组 下面我们以两个高阶微分方程组构成的 ODEs为例介绍如何将它变换成一个一阶显式微分方程组.Step 1将微分方程的最高阶变量移到等式左边， 其它移到右边，并按阶次从 低到高排列.形式为：x（m）=f（t,x,x,x,川,x（m“y, y,yJ|,y（na）yn）=g（t,x,x,x,川xZyyyMyg）