1、求解器ODE类型特点说明ode45非刚性单步算法:4、5阶Runge-Kutta方程;累计截断误差大部分场合的首选算法ode232、3阶Runge-Kutta方程;使用于精度较低的情形ode113多步法:Adams算法;高低精度可达计算时间比ode45短ode23t适度刚性采用梯形算法适度刚性情形ode15s刚性Gears反向数值微分;精度中等若ode45失效时,可尝试使用ode23s单步法:2阶Rosebrock算法;低精度当精度较低时,计算时间比ode15s短ode23tb梯形算法;ode23、ode45是极其常用的用来求解非刚性的标准形式的一阶微分方程(组)的初值问题的解的Matlab常
2、用程序,其中:ode23采用龙格-库塔2阶算法,用3阶公式作误差估计来调节步长,具有低等的精度.ode45则采用龙格-库塔4阶算法,用5阶公式作误差估计来调节步长,具有中等的精度.3在matlab命令窗口、程序或函数中创建局部函数时,可用内联函数inline,inline函数形式相当于编写M函数文件,但不需编写M-文件就可以描述出某种数学关系.调用inline函数,只能由一个matlab表达式组成,并且只能返回一个变量,不允许u,v这种向量形式.因而,任何要求逻辑运算或乘法运算以求得最终结果的场合,都不能应用inline函数,inline函数的一般形式为:FunctionName=inline
3、(函数内容, 所有自变量列表)例如:(求解F(x)=x2*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 求解微分方程程序:syms x
4、y; y=dsolve(Dy+2*x*y=x*exp(-x2),x)例2 求微分方程在初始条件下的特解并画出解函数的图形. y=dsolve(x*Dy+y-exp(1)=0,y(1)=2*exp(1),x);ezplot(y)例 3 求解微分方程组在初始条件下的特解并画出解函数的图形.syms x y t x,y=dsolve(Dx+5*x+y=exp(t),Dy-x-3*y=0x(0)=1y(0)=0t)simple(x);simple(y)ezplot(x,y,0,1.3);axis auto2.用ode23、ode45等求解非刚性标准形式的一阶微分方程(组)的初值问题的数值解(近似解)例
5、 4 求解微分方程初值问题的数值解,求解范围为区间0,0.5.fun=inline(-2*y+2*x2+2*xxy);x,y=ode23(fun,0,0.5,1);plot(x,y,o-例 5 求解微分方程的解,并画出解的图形.分析:这是一个二阶非线性方程,我们可以通过变换,将二阶方程化为一阶方程组求解.令,则编写M-文件vdp.mfunction fy=vdp(t,x)fy=x(2);7*(1-x(1)2)*x(2)-x(1);end在Matlab命令窗口编写程序y0=1;0t,x=ode45(vdp,0,40,y0);或t,x=ode45(vdp,0,40,y0);y=x(:,1);dy=
6、x(:,2);plot(t,y,t,dy)练习与思考:M-文件vdp.m改写成inline函数程序?3.用Euler折线法求解Euler折线法求解的基本思想是将微分方程初值问题化成一个代数(差分)方程,主要步骤是用差商替代微商,于是记从而于是例 6 用Euler折线法求解微分方程初值问题的数值解(步长取0.4),求解范围为区间0,2.本问题的差分方程为 clear f=sym(y+2*x/y2 a=0; b=2; h=0.4; n=(b-a)/h+1; x=0; y=1; szj=x,y;%数值解 for i=1:n-1 y=y+h*subs(f,x,y);%subs,替换函数 x=x+h;
7、szj=szj;x,y; endszj 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法的迭代公式为相应的Matlab程序为: l1=subs(f, 替换函数 l2=subs(f, ,x+h/2
8、,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求解问题并比较差异.(2)利用Matlab求微分方程的解.(3)求解微分方程的特解.(4)利用Matlab求微分方程初值问题的解.提醒:尽可能多的考虑解法三微分方程转换为一阶显式微分方程组 Matlab微分方程解算器只能求解标准形式的一阶显式微分方程(组)问题,因此在使用ODE解算器之前,我们需要做的第一步,也是最重要的一步就是借助状态变量将微分方程(组)化成Matlab可接受的标准形式.当然,如
9、果ODEs由一个或多个高阶微分方程给出,则我们应先将它变换成一阶显式常微分方程组.下面我们以两个高阶微分方程组构成的ODEs为例介绍如何将它变换成一个一阶显式微分方程组.Step 1 将微分方程的最高阶变量移到等式左边,其它移到右边,并按阶次从低到高排列.形式为:Step 2 为每一阶微分式选择状态变量,最高阶除外ODEs中所有是因变量的最高阶次之和就是需要的状态变量的个数,最高阶的微分式不需要给它状态变量.Step 3 根据选用的状态变量,写出所有状态变量的一阶微分表达式(1)求解微分方程组其中(2)求解隐式微分方程组提示:使用符号计算函数solve求,然后利用求解微分方程的方法四偏微分方程
10、解法Matlab提供了两种方法解决PDE问题,一是使用pdepe函数,它可以求解一般的PDEs,具有较大的通用性,但只支持命令形式调用;二是使用PDE工具箱,可以求解特殊PDE问题,PDEtoll有较大的局限性,比如只能求解二阶PDE问题,并且不能解决片微分方程组,但是它提供了GUI界面,从复杂的编程中解脱出来,同时还可以通过FileSave As直接生成M代码.1.一般偏微分方程(组)的求解(1)Matlab提供的pdepe函数,可以直接求解一般偏微分方程(组),它的调用格式为:sol=pdepe(m,pdefun,pdeic,pdebc,x,t)pdefun是PDE的问题描述函数,它必须换
11、成标准形式:这样,PDE就可以编写入口函数:c,f,s=pdefun(x,t,u,du),m,x,t对应于式中相关参数,du是u的一阶导数,由给定的输入变量可表示出c,f,s这三个函数.pdebc是PDE的边界条件描述函数,它必须化为形式:于是边值条件可以编写函数描述为:pa,qa,pb,qb=pdebc(x,t,u,du),其中a表示下边界,b表示上边界.pdeic是PDE的初值条件,必须化为形式:,故可以使用函数描述为:u0=pdeic(x)sol是一个三维数组,sol(:,:,i)表示的解,换句话说,对应x(i)和t(j)时的解为sol(i,j,k),通过sol,我们可以使用pdeval函数直接计算某个点的函数值.(2)实例说明求解偏微分其中,且满足初始条件及边界条件解:(1)对照给出的偏微分方程和pdepe函数求解的标准形式,原方程改写为可见%目标PDE函数funct
copyright@ 2008-2022 冰豆网网站版权所有
经营许可证编号:鄂ICP备2022015515号-1