(整理)实验2 利用matlab解(非)线性、微分方程(组).
matlab求解微分方程组用matlab对微分方程求解实验报告
matlab求解微分方程组用matlab对微分方程求解实验报告导读：就爱阅读网友为您分享以下“用matlab对微分方程求解实验报告”的资讯，希望对您有所帮助，感谢您对的支持!o 《高等数学》上机作业（三）一、上机目的1、学会用Matlab 求简单微分方程的解析解。
2、学会用Matlab 求微分方程的数值解。
二、上机内容1、求简单微分方程的解析解.2、求微分方程的数值解.3、数学建模实例.4、上机作业. 三、上机作业1. 求微分方程：xy ' ? y ? ex?0在初值条件y (1 ) ? 2 e 下的特解，并画出解函数的图形. 命令>> y=dsolve('x*Dy+y-exp(x)=0','y(1)=2*exp(1) ','x') 运行结果：y = 1/x*exp(x)+1/x*exp(1)函数图象：2. 求微分方程的特解.?d2ydy?4?5y?0?2dx ?dx?y(0)?0,y'(1)?10?命令>> y=dsolve('D2y+4*Dy-5*y=0','y(0)=0,Dy(1)= 10','x') 运行结果：y=10/(exp(1)+5*exp(-5))*exp(x)-10/(exp(1)+5*exp(-5))*exp(-5*x)3. 鱼雷追击问题一敌舰在某海域内沿着正北方向航行时，我方战舰恰好位于敌舰的正西方向1 公里处．我舰向敌舰发射制导鱼雷，敌舰速度为0.42 公里/分，鱼雷速度为敌舰速度的2倍。
试问敌舰航行多远时将被击中？M文件x0=0; xf=0.9999999999999;[x,y]=ode15s('eq1',[x0 xf],[0 0]); plot(x,y(:,1),'b.') hold on; y=0:0.1:1;plot(1,y, '*') 运行结果图像：结论：大概在y=0.67处击中敌方舰艇！(选做)一个慢跑者在平面上沿椭圆以恒定的速率v=1跑步,设椭圆方程为: x=10+20cost, y=20+5sint. 突然有一只狗攻击他. 这只狗从原点出发,以恒定速率w跑向慢跑者,狗的运动方向始终指向慢跑者.分别求出w=20,w=5时狗的运动轨迹.W=20M文件代码function dy=eq3(t,y)dy=zeros(2,1);dy(1)=20*(10+20*cos(t)-y(1))/sqrt((10+20*cos(t)-y(1))+(20 +15*sin(t)-y(2)));dy(2)=20*(20+15*sin(t)-y(2))/sqrt((10+20*cos(t)-y(1))+(20+ 15*sin(t)-y(2)));运行命令t0=0;tf=10;[t,y]=ode45('eq3',[t0 tf],[0 0]);T=0:0.1:2*pi;X=10+20*cos(T);Y=20+15*sin(T);plot(X,Y,'-')hold onplot(y(:,1),y(:,2),'r*')运行结果：利用二分法更改tf tf=5时tf=2.5时tf=3.15时：所以在t=3.15时刻恰好追上！W=5M文件代码function dy=eq4(t,y)dy=zeros(2,1);dy(1)=5*(10+20*cos(t)-y(1))/sqrt((10+20*cos(t)-y(1))+(20+15*sin(t)-y(2)));dy(2)=5*(20+15*sin(t)-y(2))/sqrt((10+20*cos(t)-y(1))+(20+1 5*sin(t)-y(2)));命令：t0=0;tf=10;[t,y]=ode45('eq4',[t0 tf],[0 0]);T=0:0.1:2*pi;X=10+20*cos(T);Y=20+15*sin(T);plot(X,Y,'-')hold onplot(y(:,1),y(:,2),'*')运行结果更改tf=20运行结果Tf=40所以永远追不上！四、上机心得体会高等数学是工科学生的主干科目，它应用于生产生活的方方面面，通过建模，计算可以求出实际问题的最优化问题！因此我们需要掌握建模和利用专业软件处理实际问题的能力！百度搜索“就爱阅读”,专业资料,生活学习,尽在就爱阅读网,您的在线图书馆。
实验2 利用matlab解(非)线性、微分方程(组)
实验2 利用matlab 解（非）线性、微分方程（组）一、实验目的1、线性方程组的解法：直接求解法和迭代法；2、非线性方程以及非线性方程组的求解；3、微分方程的数值解。
二、实验内容1、对于下列线性方程组：⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡66136221143092x （1） 请用直接法求解；>> A=[2 9 0;3 4 11;2 2 6];>> b=[13 6 6]';>> x=A\bx =7.4000-0.2000-1.4000（2） 请用LU 分解方法求解；>> [L,U]=lu(A);>> x1=U\(L\b)x1 =7.4000-0.2000-1.4000>> [L,U,P]=lu(A);>> x2=U\(L\P*b)x2 =7.4000-0.2000-1.4000（3） 请用QR 分解方法求解；>> [Q,R]=qr(A);>> x1=R\(Q\b)x1 =7.4000-0.2000-1.4000>> [Q,R,E]=qr(A);>> x2=E*(R\(Q\b))x2 =7.4000-0.2000-1.4000（4） 请用Cholesky 分解方法求解。
>> R=chol(A)Error using cholMatrix must be positive definite.因此系数矩阵A 不是正定的，故不能用Cholesky 分解法2、设迭代精度为10-6，分别用Jacobi 迭代法、Gauss-Serdel 迭代法求解下列线性方程组，并比较此两种迭代法的收敛速度。
⎪⎩⎪⎨⎧=+-=-+-=-510272109103232121x x x x x x x 解，Ax=b,新建函数如下>> A=[10 -1 0;-1 10 -2;0 -2 10];>> b=[9 7 5]';>> eps=10e-6;>> [x,n]=jacobi(A,b,[0,0,0]',eps)x =0.99370.93680.6874n =9>> [x,n]=gauseidel(A,b,[0,0,0]',eps)x =0.99370.93680.6874n =6故本例中Gauss-Serdel 迭代法优于Jacobi 迭代法3、求解非线性方程010=-+x xe x 在2附近的根。
用matlab对微分方程求解实验报告.
o 《高等数学》上机作业(三一、上机目的1、学会用 M a t l a b 求简单微分方程的解析解。
2、学会用 M a t l a b 求微分方程的数值解。
二、上机内容1、求简单微分方程的解析解.2、求微分方程的数值解.3、数学建模实例.4、上机作业. 三、上机作业1. 求微分方程:在初值条件下的特解,并画出解函数的图形. 命令>> y =d s o l v e ('x *D y +y -e x p (x =0','y (1=2*e x p (1','x ' 运行结果:y = 1/x *e x p (x +1/x *e x p (1'xxy y e +-=12(y e =函数图象:2. 求微分方程的特解.22450(00,'(110d y dyy dx dx y y ?+-=???==?命令>> y=dsolve('D2y+4*Dy-5*y=0','y(0=0,Dy(1=10','x' 运行结果:y=10/(exp(1+5*exp(-5*exp(x-10/(exp(1+5*exp(-5*exp(-5*x3. 鱼雷追击问题一敌舰在某海域内沿着正北方向航行时,我方战舰恰好位于敌舰的正西方向 1 公里处.我舰向敌舰发射制导鱼雷,敌舰速度为0.42 公里/分,鱼雷速度为敌舰速度的2倍。
试问敌舰航行多远时将被击中?M文件x0=0; xf=0.9999999999999; [x,y]=ode15s('eq1',[x0 xf],[0 0]; plot(x,y(:,1,'b.'hold on;y=0:0.1:1;plot(1,y, '*'运行结果图像:结论:大概在y=0.67处击中敌方舰艇!(选做一个慢跑者在平面上沿椭圆以恒定的速率v=1跑步,设椭圆方程为:x=10+20cost, y=20+5sint. 突然有一只狗攻击他. 这只狗从原点出发,以恒定速率w跑向慢跑者,狗的运动方向始终指向慢跑者.分别求出w=20,w=5时狗的运动轨迹.W=20M文件代码function dy=eq3(t,ydy=zeros(2,1;dy(1=20*(10+20*cos(t-y(1/sqrt((10+20*cos(t-y(1^2+(2 0+15*sin(t-y(2^2;dy(2=20*(20+15*sin(t-y(2/sqrt((10+20*cos(t-y(1^2+(2 0+15*sin(t-y(2^2;运行命令t0=0;tf=10;[t,y]=ode45('eq3',[t0 tf],[0 0];T=0:0.1:2*pi;X=10+20*cos(T;Y=20+15*sin(T;plot(X,Y,'-'hold onplot(y(:,1,y(:,2,'r*'运行结果:利用二分法更改tf tf=5时tf=2.5时tf=3.15时:所以在t=3.15时刻恰好追上!W=5M文件代码function dy=eq4(t,ydy=zeros(2,1;dy(1=5*(10+20*cos(t-y(1/sqrt((10+20*cos(t-y(1^2+(20 +15*sin(t-y(2^2; dy(2=5*(20+15*sin(t-y(2/sqrt((10+20*cos(t-y(1^2+(20 +15*sin(t-y(2^2; 命令:t0=0;tf=10;[t,y]=ode45('eq4',[t0 tf],[0 0]; T=0:0.1:2*pi;X=10+20*cos(T; Y=20+15*sin(T; plot(X,Y,'-'hold onplot(y(:,1,y(:,2,'*' 运行结果更改tf=20运行结果Tf=40 11所以永远追不上！四、上机心得体会高等数学是工科学生的主干科目，它应用于生产生活的方方面面，通过建模，计算可以求出实际问题的最优化问题！因此我们需要掌握建模和利用专业软件处理实际问题的能力！ 12。
数学建模实验二：微分方程模型Matlab求解与分析
实验二： 微分方程模型Matlab 求解与分析一、实验目的[1] 掌握解析、数值解法，并学会用图形观察解的形态和进行解的定性分析； [2] 熟悉MATLAB 软件关于微分方程求解的各种命令；[3] 通过范例学习建立微分方程方面的数学模型以及求解全过程； [4] 熟悉离散 Logistic 模型的求解与混沌的产生过程。
二、实验原理1. 微分方程模型与MATLAB 求解解析解用MATLAB 命令dsolve(‘eqn1’,’eqn2’, ...) 求常微分方程（组）的解析解。
其中‘eqni'表示第i 个微分方程，Dny 表示y 的n 阶导数,默认的自变量为t 。
（1） 微分方程 例1 求解一阶微分方程 21y dxdy+= (1) 求通解 输入：dsolve('Dy=1+y^2')输出：ans =tan(t+C1)（2）求特解 输入：dsolve('Dy=1+y^2','y(0)=1','x')指定初值为1，自变量为x 输出：ans =tan(x+1/4*pi)例2 求解二阶微分方程 221()04(/2)2(/2)2/x y xy x y y y πππ'''++-=='=-原方程两边都除以2x ，得211(1)04y y y x x'''++-= 输入:dsolve('D2y+(1/x)*Dy+(1-1/4/x^2)*y=0','y(pi/2)=2,Dy(pi/2)=-2/pi','x')ans =- (exp(x*i)*(pi/2)^(1/2)*i)/x^(1/2) +(exp(x*i)*exp(-x*2*i)*(pi/2)^(3/2)*2*i)/(pi*x^(1/2))试试能不用用simplify 函数化简 输入: simplify(ans)ans =2^(1/2)*pi^(1/2)/x^(1/2)*sin(x) （2）微分方程组例3 求解 d f /d x =3f +4g ; d g /d x =-4f +3g 。
使用Matlab进行微分方程求解的方法
使用Matlab进行微分方程求解的方法引言微分方程是数学中非常重要的一部分，广泛应用于物理、经济、工程等领域。
对于大部分微分方程的解析解往往难以求得，而数值解法则成为了一种常用的解决手段。
Matlab作为一种强大的科学计算软件，也提供了丰富的工具和函数用于求解微分方程，本文将介绍一些常见的使用Matlab进行微分方程求解的方法。
一、数值求解方法1. 欧拉方法欧拉方法是最简单的一种数值求解微分方程的方法，它将微分方程的微分项用差分的方式进行近似。
具体的公式为：y(n+1) = y(n) + hf(x(n), y(n))其中，y(n)表示近似解在第n个点的值，h为步长，f(x, y)为微分方程的右端项。
在Matlab中使用欧拉方法进行求解可以使用ode113函数，通过设定不同的步长，可以得到不同精度的数值解。
2. 中点法中点法是较为精确的一种数值求解微分方程的方法，它的计算公式为：k1 = hf(x(n), y(n))k2 = hf(x(n) + h/2, y(n) + k1/2)y(n+1) = y(n) + k2中点法通过计算两个斜率的平均值来得到下一个点的值，相较于欧拉方法，中点法能提供更精确的数值解。
3. 4阶龙格库塔法龙格库塔法是一类高阶数值求解微分方程的方法，其中4阶龙格库塔法是最常用的一种。
它的计算公式为：k1 = hf(x(n), y(n))k2 = hf(x(n) + h/2, y(n) + k1/2)k3 = hf(x(n) + h/2, y(n) + k2/2)k4 = hf(x(n) + h, y(n) + k3)y(n+1) = y(n) + (k1 + 2k2 + 2k3 + k4)/64阶龙格库塔法通过计算多个斜率的加权平均值来得到下一个点的值，相较于欧拉方法和中点法，它的精度更高。
二、Matlab函数和工具除了可以使用以上的数值方法进行微分方程求解之外，Matlab还提供了一些相关的函数和工具，方便用户进行微分方程的建模和求解。
实验二_基于Matlab的微分方程数值解法
实验二微分方程数值解法一．实验原理及实验内容：对微分方程描述的控制系统，利用欧拉法、二阶龙格-库塔法、四阶龙格-库塔法分别编写M文件，进行数值计算和作图。
1．分别用欧拉法、二阶龙格-库塔法、四阶龙格-库塔法求下面系统的输出响应y(t)在0≤t≤1上，h=0.1时的数值解。
'2,(0)1=-=y y y要求保留4位小数，并将三种方法的结果与真解2=进行比较。
()ty t e-2.若为如2y y y==何编程计算？',(0)1二．实验仪器：计算机Matlab软件三．实验数据记录：程序一：disp('欧拉算法');y=1;h=0.1;for i=0:0.1:1disp(y);y=y+h*(-2*y);enddisp('欧拉算法');ydisp('精确解');yy=exp(-2*t)h=0.1;disp('函数的2阶数值解为');disp('y=');y=1;for t=0:h:1;disp(y);k1=-2*y;k2=-2*(y+k1*h);y(i+1)=y(i)+(k1+k2)*h*1/2;endh=0.1;disp('函数的4阶数值解为');disp('y=');y=1;for t=0:h:1;disp(y);k1=-2*y;k2=-2*(y+k1*h*1/2);k3=-2*(y+k2*h*1/2);k4=-2*(y+k3*h);y=y+h*1/6*(k1+2*k2+2*k3+k4); end>>程序2：t=0:0.1:1;n=length(t);y(1)=1;h=0.1;for i=1:n-1y(i+1)=y(i)+h*(y(i)*y(i)); enddisp('欧拉算法');ydisp('精确解');yy=exp(-2*t)h=0.1;disp('函数的2阶数值解为');disp('y=');y=1;for t=0:h:1;disp(y);k1=y*y;k2=(y+k1*h)^2;y=y+(k1+k2)*h*1/2;endh=0.1;disp('函数的4阶数值解为');disp('y=');y=1;disp(y);k1=y*y;k2=(y+k1*h*1/2)^2;k3=(y+k2*h*1/2)^2;k4=(y+k3*h)^2;y=y+h*1/6*(k1+2*k2+2*k3+k4); end。
用MATLAB求解微分方程及微分方程组
注意:
1、在解n个未知函数的方程组时，x0和x均为n维向量， m-文件中的待解方程组应以x的分量形式写成. 2、使用Matlab软件求数值解时，高阶微分方程必须 等价地变换成一阶微分方程组.
例4
d 2x dx 2 1000 (1 x ) x 0 2 dt dt x (0 ) 2; x ' (0 ) 0
y '1 y 2 1 2 y2 ' 1 y 1 /( 1 x ) 5
解法三(建立参数方程求数值解)
设时刻t乙舰的坐标为(X(t)，Y(t))，导弹的坐标为(x(t)，y(t)). dx 2 dy 2 2 ( ) ( ) w 1．设导弹速度恒为 w ，则 （1） dt dt
（2）
（3）
3．因乙舰以速度v0沿直线x=1运动，设v0=1，则w=5，X=1，Y=t
因此导弹运动轨迹的参数方程为：
5 dx (1 x ) dt 2 2 (1 x ) ( t y ) 5 dy (t y ) 2 2 (1 x ) ( t y ) dt x (0) 0, y (0) 0
解: 令 y1=x，y2=y1’
则微分方程变为一阶微分方程组：
y1 ' y 2 2 y 2 ' 1000 (1 y 1 ) y 2 y 1 y (0) 2, y (0 ) 0 1 2
3、结果如图
2 1.5 1 0.5 0
1、建立m-文件vdp1000.m如下： function dy=vdp1000(t,y) dy=zeros(2,1); dy(1)=y(2); dy(2)=1000*(1-y(1)^2)*y(2)-y(1);
实验二MATLAB数值计算常微分方程（组）的求解
实验⼆MATLAB数值计算常微分⽅程（组）的求解实验⼆ MATLAB 数值计算：常微分⽅程(组)的求解⼀、实验⽬的在物理学和⼯程技术上，很多问题都可以⽤⼀个或⼀组常微分⽅程来描述，因此要解决相应的实际问题往往需要⾸先求解对应的微分⽅程。
在⼤多数情况下这些微分⽅程通常是⾮线性的或者是超越⽅程(⽐如范德堡⽅程，波导本征值⽅程等)，因此往往需要使⽤计算机数值求解。
MATLAB 作为⼀种强⼤的科学计算语⾔，其在数值计算和数据的可视化⽅⾯具有⽆以伦⽐的优势。
在解决常微分⽅程问题上，MATLAB 就提供了多种可适⽤于不同场合(如刚性和⾮刚性问题)下的求解器(Solver)，例如ode45，ode15s ，ode23，ode23s 等等。
本次实验将以范德堡⽅程的计算和地球卫星的运⾏轨道的仿真为例，练习使⽤MATLAB 的常微分⽅程求解器，以期达到如下⼏个⽬的：1. 熟悉常微分⽅程的求解⽅法，了解状态⽅程的概念；2. 能熟练使⽤dsolve 函数解析求解常微分⽅程；3. 能熟练运⽤ode45、ode15s 求解器分别数值求解⾮刚性和刚性常微分⽅程；4. 学习⽤求解器来绘制相图的⽅法。
⼆、实验的预备知识1．微分⽅程的概念未知的函数以及它的某些阶的导数连同⾃变量都由⼀已知⽅程联系在⼀起的⽅程称为微分⽅程。
如果未知函数是⼀元函数，称为常微分⽅程（Ordinary differential equations ，简称odes ）。
n 阶常微分⽅程的⼀般形式（隐式）为：0),,",',,()(=n y y y y t F （1）其中t 为⾃变量。
如果未知函数是多元函数，成为偏微分⽅程。
联系⼀些未知函数的⼀组微分⽅程组称为微分⽅程组。
微分⽅程中出现的未知函数的导数的最⾼阶解数称为微分⽅程的阶。
若⽅程中未知函数及其各阶导数都是⼀次的，称为线性常微分⽅程，⼀般表⽰为)()(')()(1)1(1)(t b y t a y t a y t a y n n n n =++++--若上式中的系数a i (t), i =1,2,…,n 均与t ⽆关，称之为常系数。
实验二MATLAB数值计算常微分方程的求解
实验二MATLAB数值计算常微分方程的求解引言:微分方程是描述物理问题和工程问题的重要工具，也是数学和工程学科中的重要课题。
解析方法可以求得一些简单微分方程的解析解，但对于复杂问题，常常无法得到解析解。
此时，数值求解方法成为一种有效的工具。
MATLAB是一个功能强大的数值计算软件，对于求解常微分方程组也提供了多种数值方法。
本实验将介绍MATLAB中求解常微分方程(组)的方法，并通过实例来演示这些方法的应用。
一、MATLAB的常微分方程(组)求解函数MATLAB提供了多种函数用于求解常微分方程(组)，其中最常用的函数是ode45、ode23和ode15s。
这些函数采用不同的数值方法，精度和效率也不同。
下面分别对这些函数进行简单介绍：1. ode45函数:ode45函数采用了一种自适应的步长控制算法，可以在一个时间段内自动调整步长。
这个函数的语法是：[t,y] = ode45(odefun,tspan,y0,options)其中，odefun是一个函数句柄，表示待解的常微分方程组，tspan是求解区间，y0是初始条件，options是一个结构体，包含其他的参数和选项。
2. ode23函数:ode23函数也采用了自适应的步长控制算法，但相对于ode45，它是一个简化版本，只允许使用比较简单的问题。
这个函数的语法与ode45相似：[t,y] = ode23(odefun,tspan,y0,options)3. ode15s函数:ode15s函数采用了一种隐式数值方法，适用于比较刚性（stiff）的常微分方程。
这个函数的语法与前两个函数也相似：[t,y] = ode15s(odefun,tspan,y0,options)以上是MATLAB提供的三种解常微分方程(组)的函数，根据问题的特点和求解的需求选择相应的函数。
二、实例演示在本实验中，我们将通过一个实例来演示如何使用MATLAB解常微分方程组。
假设有一个简单的线性常微分方程组：dy1/dt = -2y1 + y2dy2/dt = y1 - 2y2给定初始条件为y1(0)=1，y2(0)=0，求在t=[0,5]时，y1和y2的解。
实验二MATLAB 求微分方程的解
实验二 微分方程求解一、问题背景与实验目的实际应用问题通过数学建模所归纳而得到的方程，绝大多数都是微分方程，真正能得到代数方程的机会很少．另一方面，能够求解的微分方程也是十分有限的，特别是高阶方程和偏微分方程（组）．这就要求我们必须研究微分方程（组）的解法，既要研究微分方程（组）的解析解法（精确解），更要研究微分方程（组）的数值解法（近似解）．对微分方程（组）的解析解法(精确解)，Matlab 有专门的函数可以用，本实验将作一定的介绍．本实验将主要研究微分方程(组)的数值解法（近似解），重点介绍 Euler 折线法．二、相关函数（命令）及简介1．dsolve('equ1','equ2',…)：Matlab 求微分方程的解析解．equ1、equ2、…为方程（或条件）．写方程（或条件）时用 Dy 表示y 关于自变量的一阶导数，用用 D2y 表示 y 关于自变量的二阶导数，依此类推．2．simplify(s)：对表达式 s 使用 maple 的化简规则进行化简． 例如： syms xsimplify(sin(x)^2 + cos(x)^2) ans=13．[r,how]=simple(s)：由于 Matlab 提供了多种化简规则，simple 命令就是对表达式 s 用各种规则进行化简，然后用 r 返回最简形式，how 返回形成这种形式所用的规则．例如： syms x[r,how]=simple(cos(x)^2-sin(x)^2) r = cos(2*x) how = combine4．[T,Y] = solver(odefun,tspan,y 0) 求微分方程的数值解． 说明：(1) 其中的 solver 为命令 ode45、ode23、ode113、ode15s 、ode23s 、ode23t 、ode23tb 之一．(2) odefun是显式常微分方程：⎪⎩⎪⎨⎧==00)(),(yt y y t f dt dy(3) 在积分区间 tspan =],[0f t t 上，从0t 到f t ，用初始条件0y 求解．(4) 要获得问题在其他指定时间点 ,210,,t t t 上的解，则令 tspan =],,,[,210f t t t t （要求是单调的）．(5) 因为没有一种算法可以有效地解决所有的 ODE 问题，为此，Matlab 提供了多种求解器 Solver ，对于不同的ODE 问题，采用不同的Solver ．(6) 要特别的是：ode23、ode45 是极其常用的用来求解非刚性的标准形式的一阶常微分方程(组)的初值问题的解的 Matlab 的常用程序，其中：ode23 采用龙格-库塔2 阶算法，用3 阶公式作误差估计来调节步长，具有低等的精度．ode45 则采用龙格-库塔4 阶算法，用5 阶公式作误差估计来调节步长，具有中等的精度．5．ezplot(x,y ,[tmin,tmax])：符号函数的作图命令．x,y 为关于参数t 的符号函数，[tmin,tmax] 为 t 的取值范围．6．inline()：建立一个内联函数．格式：inline('expr', 'var1', 'var2',…) ，注意括号里的表达式要加引号．例：Q = dblquad(inline('y*sin(x)'), pi, 2*pi, 0, pi)三、实验内容1. 几个可以直接用 Matlab 求微分方程精确解的例子： 例1：求解微分方程22xxexy dxdy -=+，并加以验证．求解本问题的Matlab 程序为：syms x y %line1 y=dsolve('Dy+2*x*y=x*exp(-x^2)','x') %line2diff(y ,x)+2*x*y-x*exp(-x^2) %line3 simplify(diff(y ,x)+2*x*y-x*exp(-x^2)) %line4 说明：(1) 行line1是用命令定义x,y 为符号变量．这里可以不写，但为确保正确性，建议写上；(2) 行line2是用命令求出的微分方程的解：1/2*exp(-x^2)*x^2+exp(-x^2)*C1(3) 行line3使用所求得的解．这里是将解代入原微分方程，结果应该为0，但这里给出：-x^3*exp(-x^2)-2*x*exp(-x^2)*C1+2*x*(1/2*exp(-x^2)*x^2+exp(-x^2)*C1)(4) 行line4 用 simplify() 函数对上式进行化简，结果为 0， 表明)(x y y =的确是微分方程的解．例2：求微分方程0'=-+x e y xy 在初始条件e y 2)1(=下的特解，并画出解函数的图形．求解本问题的 Matlab 程序为： syms x yy=dsolve('x*Dy+y-exp(x)=0','y(1)=2*exp(1)','x') ezplot(y)微分方程的特解为：y=1/x*exp(x)+1/x* exp (1) (Matlab 格式)，即xe e y x+=，解函数的图形如图 1：图1例3：求微分方程组⎪⎪⎩⎪⎪⎨⎧=--=++035y x dt dy e y x dtdx t在初始条件0|,1|00====t t y x 下的特解，并画出解函数的图形．求解本问题的 Matlab 程序为： 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 auto微分方程的特解(式子特别长)以及解函数的图形均略． 2. 用ode23、ode45等求解非刚性的标准形式的一阶常微分方程(组)的初值问题的数值解(近似解)．例4：求解微分方程初值问题⎪⎩⎪⎨⎧=++-=1)0(2222y x x y dxdy 的数值解，求解范围为区间[0, 0.5]．fun=inline('-2*y+2*x^2+2*x','x','y'); [x,y]=ode23(fun,[0,0.5],1); x'; y';plot(x,y ,'o-') >> x' ans =0.0000 0.0400 0.0900 0.1400 0.1900 0.2400 0.2900 0.3400 0.3900 0.4400 0.4900 0.5000 >> y' ans =1.0000 0.9247 0.8434 0.7754 0.7199 0.6764 0.6440 0.6222 0.6105 0.6084 0.6154 0.6179 图形结果为图 2．图2例 5：求解描述振荡器的经典的 V er der Pol 微分方程.7,0)0(',1)0(,0)1(222====+--μμy y y dtdy y dty d分析：令,,121dtdx x y x ==则.)1(,1221221x x x dtdx x dtdx --==μ先编写函数文件verderpol.m ： function xprime = verderpol(t,x) global mu;xprime = [x(2);mu*(1-x(1)^2)*x(2)-x(1)]; 再编写命令文件vdp1.m ： global mu; mu = 7; y0=[1;0][t,x] = ode45('verderpol',[0,40],y0); x1=x(:,1);x2=x(:,2); plot(t,x1)图形结果为图3．图33. 用 Euler 折线法求解前面讲到过，能够求解的微分方程也是十分有限的．下面介绍用 Euler 折线法求微分方程的数值解(近似解)的方法．Euler 折线法求解的基本思想是将微分方程初值问题⎪⎩⎪⎨⎧==00)(),,(yx y y x f dx dy化成一个代数方程，即差分方程，主要步骤是用差商hx y h x y )()(-+替代微商dxdy ,于是：⎪⎩⎪⎨⎧==-+)()),(,()()(00x y y x y x f h x y h x y k k k k 记)(,1k k k k x y y h x x =+=+，从而)(1h x y y k k +=+，则有1,,2,1,0).,(,),(1100-=⎪⎩⎪⎨⎧+=+==++n k y x hf y yh x x x y y k k k k k k 例 6：用 Euler 折线法求解微分方程初值问题⎪⎩⎪⎨⎧=+=1)0(,22y y x y dxdy 的数值解(步长h 取0.4)，求解范围为区间[0,2]．解：本问题的差分方程为1,,2,1,0).2),( ),(,,4.0,1,021100-=⎪⎪⎪⎩⎪⎪⎪⎨⎧+=+=+====++n k y x y y x f y x hf y y h x x h y x k k k k k k （其中： 相应的Matlab 程序见附录 1． 数据结果为：0 1.0000 0.4000 1.4000 0.8000 2.1233 1.2000 3.1145 1.6000 4.4593 2.0000 6.3074图形结果见图4：图4特别说明：本问题可进一步利用四阶 Runge-Kutta 法求解，读者可将两个结果在一个图中显示，并和精确值比较，看看哪个更“精确”？（相应的 Matlab 程序参见附录 2）．四、自己动手1. 求微分方程0sin 2')1(2=-+-x xy y x 的通解．2. 求微分方程x e y y y x sin 5'2''=+-的通解．3. 求微分方程组⎪⎪⎩⎪⎪⎨⎧=-+=++00y x dtdy y x dtdx在初始条件0|,1|00====t t y x 下的特解，并画出解函数()y f x =的图形． 4. 分别用 ode23、ode45 求上述第 3 题中的微分方程初值问题的数值解(近似解)，求解区间为[0,2]t ∈．利用画图来比较两种求解器之间的差异．5. 用 Euler 折线法求解微分方程初值问题⎪⎩⎪⎨⎧=-=1)0(,12'32y y xy y 的数值解(步长h 取0.1)，求解范围为区间[0,2]．6. 用四阶 Runge-Kutta 法求解微分方程初值问题⎩⎨⎧=-=1)0(,cos 'y x e y y x 的数值解(步长h 取0.1)，求解范围为区间[0,3]．四阶 Runge-Kutta 法的迭代公式为(Euler 折线法实为一阶 Runge-Kutta 法)：1,,2,1,0),()2,2()2,2(),()22(6,),(342312143211100-=⎪⎪⎪⎪⎪⎪⎪⎩⎪⎪⎪⎪⎪⎪⎪⎨⎧++=++=++==++++=+==++n k hL y h x f L L h y h x f L L h y h x f L y x f L L L L L hy y h x x x y y k k k k k k k k k k k k 相应的 Matlab 程序参见附录 2．试用该方法求解第5题中的初值问题． 7. 用 ode45 方法求上述第 6 题的常微分方程初值问题的数值解(近似解)，从而利用画图来比较两者间的差异．五、附录附录 1：(fulu1.m)clearf=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}); x=x+h;szj=[szj;x,y]; end szjplot(szj(:,1),szj(:,2))附录 2：(fulu2.m)clearf=sym('y-exp(x)*cos(x)'); a=0; b=3; h=0.1;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];endszjplot(szj(:,1),szj(:,2))。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
实验2 利用matlab 解（非）线性、微分方程（组）
一、实验目的
1、线性方程组的解法：直接求解法和迭代法；
2、非线性方程以及非线性方程组的求解；
3、微分方程的数值解。
二、实验内容
1、对于下列线性方程组：
⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡66136221143092x
（1）
（2） 请用直接法求解；
（3）
（4） 请用LU 分解方法求解；
（5） 请用QR 分解方法求解；
（6）
（7） 请用Cholesky 分解方法求解。
2、设迭代精度为10-6，分别用Jacobi 迭代法、Gauss-Serdel 迭代法求解下列线性方程组，并比较此两种迭代法的收敛速度。
⎪⎩
⎪⎨⎧=+-=-+-=-510272109103232121x x x x x x x 3、求解非线性方程010=-+-x xe x 在2附近的根。
4、求下列非线性方程组在(0.5,0.5) 附近的数值解。
⎪⎩⎪⎨⎧=-+=-+02)sin(02)cos(y x xe y ye x
5、通过画图方法描述某非刚性体的运动方程的微分方程⎪⎩⎪⎨⎧=-==2
1131232151.0y y y y y y
y y y ，其初始条件为⎪⎩⎪⎨⎧===1)0(1)0(0)0(3
21y y y 。
6、求二阶微分方程)2sin(3t y e y t y t =-'+''， 1)0(=y ，1)0(-='y 在20≤≤t 时的数值图解。