实验四-求微分方程的解--数学软件与数学实验-教学课件.ppt
实验四求微分方程的解数学实验q 自牛顿发明微积分以来,微分方程在描述事物运动规律上已发挥了重要的作用。实际应用问题通过数学建模所得到的方程,绝大多数是微分方程。q 由于实际应用的需要,人们必须求解微分方程。然而能够求得解析解的微分方程十分有限,绝大多数微分方程需要利用数值方法来近似求解。q 本实验主要研究如何用 Matlab 来计算微分方程(组)的数值解,并重点介绍一个求解微分方程的基本数值解法Euler折线法。问题背景和实验目的q 考虑一维经典初值问题u 基本思想:用差商代替微商根据 Talyor 公式,y(x)在点 xk 处有Euler 折线法初值问题的Euler折线法q 具体步骤:等距剖分:步长:u 分割求解区间u 差商代替微商得方程组:分割求解区间,差商代替微商,解代数方程 为分割点k=0,1,2,.,n-1yk 是 y(xk)的近似Euler 折线法源程序clearf=sym(y+2*x/y2);a=0;b=2;h=0.4;n=(b-a)/h+1;%n=(b-a)/h;x=0;y=1;szj=x,y;for i=1:n-1%i=1:n y=y+h*subs(f,x,y,x,y);x=x+h;szj=szj;x,y;endszjplot(szj(:,1),szj(:,2),or-)Euler折线法举例(续)解析解:解析解近似解y=1/3*(-18-54*x+45*exp(3*x)(1/3)Runge-Kutta 方法q 为了减小误差,可采用以下方法:u 让步长 h 取得更小一些;u 改用具有较高精度的数值方法:q 龙格-库塔方法Runge-Kutta(龙格-库塔)方法u 是一类求解常微分方程的数值方法u 有多种不同的迭代格式四阶 R-K 方法源程序clear;f=sym(y+2*x/y2);a=0;b=2;h=0.4;n=(b-a)/h+1;%n=(b-a)/h;x=0;y=1;szj=x,y;for i=1:n-1%i=1:n l1=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;endplot(szj(:,1),szj(:,2),dg-)Runge-Kutta 方法Matlab 解初值问题q 用 Maltab自带函数 解初值问题u 求解析解:dsolveu 求数值解:ode45、ode23、ode113、ode23t、ode15s、ode23s、ode23tbdsolve 的使用q 几点说明l 如果省略初值条件,则表示求通解;l 如果省略自变量,则默认自变量为 t dsolve(Dy=2*x,x);dy/dx=2xdsolve(Dy=2*x);dy/dt=2xl 若找不到解析解,则返回其积分形式。l 微分方程中用 D 表示对 自变量 的导数,如:Dy y;D2y y;D3y ydsolve 举例例 2:求微分方程 在初值条件 下的特解,并画出解函数的图形。y=dsolve(x*Dy+y-exp(x)=0,y(1)=2*exp(1),x)ezplot(y);v dsolve 举例例:x,y=dsolve(Dx+5*x=0,Dy-3*y=0,.x(0)=1,y(0)=1,t)r=dsolve(Dx+5*x=0,Dy-3*y=0,.x(0)=1,y(0)=1,t)这里返回的 r 是一个 结构类型 的数据r.x%查看解函数 x(t)r.y%查看解函数 y(t)只有很少一部分微分方程(组)能求出解析解。大部分微分方程(组)只能利用数值方法求数值解。dsolve的输出个数只能为一个 或 与方程个数相等。Matlab函数数值求解T,Y=solver(odefun,tspan,y0)其中 y0 为初值条件,tspan为求解区间;Matlab在数值求解时自动对求解区间进行分割,T(列向量)中返回的是分割点的值(自变量),Y(数组)中返回的是这些分割点上的近似解,其列数等于因变量的个数。solver 为Matlab的ODE求解器(可以是 ode45、ode23、ode113、ode15s、ode23s、ode23t、ode23tb)没有一种算法可以有效地解决所有的 ODE 问题,因此MATLAB 提供了多种ODE求解器,对于不同的ODE,可以调用不同的求解器。