matlab常微分方程和常微分方程组的求解
重要：MATLAB常微分方程(组)数值解法
Matlab常微分方程求解问题分类
边值问题：
初值问题：
• 定解附加条件在自变量 的一端
• 一般形式为： y' f (x, y)
y(a)
y0
• 初值问题的数值解法一 般采用步进法，如 Runge-Kutta法
➢ 在自变量两端均给定附加 条件
y' f (x, y)
➢ 一般形式：y(a)y1, y(b)y2
1.根据常微分方程要求的求解精度与速度要求
求解初值问题：
y
'
y
2x y
y ( 0 ) 1
(0x1)
比较ode45和ode23的求解精度和速度
ode45和ode23的比较-1
function xODE clear all clc
format long
y0 = 1; [x1,y1] = ode45(@f,[0,1],y0); [x2,y2] = ode23(@f,[0,1],y0); plot(x1,y1,'k-',x2,y2,'b--') xlabel('x') ylabel('y')
rD = k(3)*C(2)-k(5)*C(4);
rE = k(4)*C(3)+k(5)*C(4);
% Mass balances dCdt = [rA; rB; rC; rD; rE];
三个串联的CSTR等温反应器（例4-3）
function IsothermCSTRs clear all clc CA0 = 1.8; % kmol/m^3 CA10 = 0.4; % kmol/m^3 CA20 = 0.2; % kmol/m^3 CA30 = 0.1; % kmol/m^3 k = 0.5; % 1/min tau = 2; stoptime = 2.9; % min [t,y] = ode45(@Equations,[0 stoptime],[CA10 CA20 CA30],[],k,CA0,tau); disp(' Results:') disp(' t CA1 CA2 CA3') disp([t,y]) plot(t,y(:,1),'k--',t,y(:,2),'b:',t,y(:,3),'r-') legend('CA_1','CA_2','CA_3') xlabel('Time (min)') ylabel('Concentration') % -----------------------------------------------------------------function dydt = Equations(t,y,k,CA0,tau) CA1 = y(1); CA2 = y(2); CA3 = y(3); dCA1dt = (CA0-CA1)/tau - k*CA1; dCA2dt = (CA1-CA2)/tau - k*CA2; dCA3dt = (CA2-CA3)/tau - k*CA3; dydt = [dCA1dt; dCA2dt; dCA3dt];
matlab-常微分方程
Events
含义 为‘on’时，控制解向量 有效值： 范数的相对误差，使每 on、off 步计算中，满足： 缺省值： norm(e)<=max(RelTol*n off orm(y),AbsTol) 有效值： 为‘on’时，返回相应的 on、off 事件记录
取值
参数设置
属性名 含义 若无输出参量，则solver 将执行下面操作之一： 有效值： 画出解向量中各元素随 odeplot、 时间的变化； odephas2、画出解向量中前两个分 odephas3、量构成的相平面图； odeprint 画出解向量中前三个分 缺省值： 量构成的三维相空间图 odeplot ； 随计算过程，显示解向 量 取值
使用于精度较低 的情形
OD 求解器 E类 Solver 型 非 ode113 刚 性 适 度 ode23t 刚 性 刚 ode15s 性
特点
说明
多步法；Adams算 计算时间比ode45 法；高低精度均可 短 -3～10-6 到10 采用梯形算法 适度刚性情形
多步法；Gear’s反 若ode45失效时， 向数值微分；精度 可尝试使用 中等
• (6)若没有给定输出参量，则在命令窗口显 示解列表。若该命令找不到解析解，则返 回一警告信息，同时返回一空的sym对象。 这时，用户可以用命令ode23或ode45求解 方程组的数值解。
y ′′ = −a y ′ y ( 0) = 1 y ′(π / a ) = 0
2
例1
例2
• >> [u,v] = dsolve('Du=v,Dv=u') u= C1*exp(-t)+C2*exp(t) V= -C1*exp(-t)+C2*exp(t)
第六章用MATLAB求解微分方程及微分方程组
运行结果：u = tan(t-c)
例 2 求微分方程的特解.
d 2 y dy 2 4 29 y 0 dx dx y (0) 0, y ' (0) 15
解 输入命令: y=dsolve('D2y+4*Dy+29*y=0','y(0)=0,Dy(0)=15','x') 运行结果为 : y =3e-2xsin（5x）
不同求解器Solver的特点
求解器 Solver ode23s
OD E类 型
刚 性
特点
说明
刚 ode23tb 性
一步法；2阶 当精度较低时， Rosebrock算法； 计算时间比 低精度 ode15s短 当精度较低时， 梯形算法；低精 计算时间比 度 ode15s短
d 2x dx 2 1000(1 x 2 ) x 0 例 4 dt dt x(0) 2; x' (0) 0 解: 令 y1=x，y2=y1’
y Me k1t 即
模型2：快速饮酒后，体液中酒精含量的变化率
dx k2 x k1Mek1t dt x(0) 0
用Matlab求解模型2：
syms x y k1 k2 M t x=dsolve('Dx+k2*x=k1*M*exp(-k1*t)','x(0)=0','t') 运行结果： M*k1/(-k1+k2)*exp(-k2*t+t*(-k1+k2))-exp(-k2*t)*M*k1/(-k1+k2) 即：
（1）
2. 由于弹头始终对准乙舰，故导弹的速度平行于乙舰与导弹头位置的差向量，
如何使用MATLAB求解微分方程(组)ppt课件
差,输出参数,事件等，可缺省。 9
使用ODE?时如何编 写微分方程 ？
方式一：带额外参数，使用时需对参数进行赋值
function odefun(t,x,flag,R,L,C) %用flag说明R、L、C为变 量
xdot=zeros(2,1);
%表明其为列向量
xdot(1)=-R/L*x(1)-1/L*x(2)+1/L*f(t);
2
Where ?
工程控制
ODE
医学生理
航空航天
金融分析
3
Where ?
算法开发 数据分析
数值计算 MAT LAB
数据可视化
4
When ?
当对问题进行建模后，有常微分方程 需要求解时。
在生物建模中，经常需要求解常微分 方程。如药物动力学的房室模型的建模 仿真。
5
方法 ode23
ode45
数 ode113
当无法藉由微积分技巧求 得解析解时，这时便只能利 用数值分析的方式来求得其 数值解了。实际情况下，常 微分方程往往只能求解出其
数值解。
在数学中，刚性方程是指一 个微分方程，其数值分析的解 只有在时间间隔很小时才会稳 定，只要时间间隔略大，其解 就会不稳定。
目前很难去精确地去定义哪 些微分方程是刚性方程，但是 大体的想法是：这个方程的解
y(1)=x(2);
y1
y2
y(2)= -t*x(1)+exp(t)*x(2)+3*sin(2*t);
end
1000
求解程序ห้องสมุดไป่ตู้键步骤
[t,y]=ode45('odefun',[3.9 4.0],[2 8])
y
matlab_常微分方程数值解法
dt 2
简朴问题可以求得解析解，多数实际问题靠数值求解 。
第4页
一阶常微分方程(ODE )初值问题 ： ODE ：Ordinary Differential Equation
dy
f
(x,
y)
dx
x0 x xn
y(x0 ) y0
数值解法就是求y(x)在某些分立旳节点 xn 上旳近似值 yn，用以近似y(xn)
x0
y0
x1 f y(x), x dx
x0
x2 f y(x), x dx
x1
y(x1) f y(x1), x1 h
第17页
同样，在[x0,xn+1] ，积分采用矩形近似，得：
y(xn1) y0
f xn1
x0
y(x), x dx
y(xn ) f y(xn ), xn h
yn y(xn )
第5页
2、欧拉近似办法
2.1 简朴欧拉(L.Euler, 1707-1783)办法。
dy
dx
f
(y, x)
y(x0 ) y0
欧拉数值算法就是由初值通过递推求解，递推求解
就是从初值开始，后一种函数值由前一种函数值得到。核 心是构造递推公式。
y0 y1 y2 yn
第6页
i 1,2,...
第36页
没有一种算法可以有效地解决所有旳 ODE 问题，因此 MATLAB 提供了多种ODE函数。
函数 ODE类
特点
阐明
型
ode45
非刚性 单步法；4，5 阶 R-K 措施；合计 大部分场合旳首选措施
截断误差为 (△x)3
ode23
非刚性 单步法；2，3 阶 R-K 措施；合计 使用于精度较低旳情形
matlab梯形法求常微分方程
近年来，随着科技的迅猛发展，人们对数学问题的求解需求也越来越迫切。
在数值分析中，常微分方程的求解一直是一个备受关注的领域。
而在这个领域中，matlab梯形法求解常微分方程成为了一种被广泛应用的方法。
那么，什么是matlab梯形法？它又是如何应用于求解常微分方程的呢？让我们来深入了解matlab梯形法。
在matlab中，梯形法是一种常用的数值求解方法，它可以用于求解常微分方程。
该方法的基本思想是将微分方程中的导数用差分代替，从而将微分方程转化为代数方程组，再利用matlab进行求解。
通过该方法，我们可以得到微分方程的数值解，从而更好地理解和分析问题。
现在，让我们来探讨matlab梯形法在求解常微分方程中的应用。
假设我们需要求解如下的一阶常微分方程：\[ \frac{dy}{dt} = f(t,y) \]其中，\( f(t,y) \) 是关于\( t \)和\( y \)的函数。
我们需要将微分方程离散化，即用差分代替导数。
通过将时间区间\( [a, b] \)进行均匀划分，我们可以得到：\[ t_0 = a, t_1, t_2, ..., t_n = b \]\[ y_0 = \alpha, y_1, y_2, ..., y_n \]\[ h = \frac{b-a}{n} \]其中，\( t_i \) 是时间节点，\( y_i \) 是对应的近似解，\( h \) 是时间步长。
接下来，我们可以利用梯形法进行求解。
梯形法的迭代公式为：\[ y_{i+1} = y_i + \frac{h}{2}[f(t_i, y_i) + f(t_{i+1}, y_{i+1})] \]通过不断迭代，我们可以得到微分方程的数值解。
在实际应用中，matlab梯形法可以很好地处理各种类型的常微分方程。
无论是线性方程还是非线性方程，matlab梯形法都能提供较为准确的数值解。
该方法还可以用于求解初值问题和边值问题，具有较好的通用性和适用性。
matlab用欧拉法求解常微分方程
matlab用欧拉法求解常微分方程在数学和科学领域中，常微分方程是一种非常有用的工具，用于描述许多自然和物理现象。
MATLAB是一种强大的数学软件，可以用来解决许多数学问题。
本文将介绍如何使用欧拉法在MATLAB中求解常微分方程。
欧拉法是一种基本的数值方法，用于近似解决微积分方程问题。
该方法使用离散时间步长，将微积分方程转换成差分方程，并不断迭代求解。
欧拉法的实现非常简单，因此它很适合初学者。
下面是使用欧拉法在MATLAB中求解常微分方程的步骤：1. 定义常微分方程以 y' = -0.5y + 3sin(t) 为例，我们先定义常微分方程。
在MATLAB中，可以使用 anonymous functions 实现：dydt = @(t,y) -0.5*y + 3*sin(t);2. 定义时间范围和时间步长我们需要定义时间范围和时间步长，以便在一定时间范围内求解差分方程。
在这个例子中，我们定义时间范围为 0 到 10，并定义时间步长为 0.1：tspan = [0 10];h = 0.1;3. 定义初始条件我们需要定义初始条件，即 y(0) 的值。
在这个例子中，我们假设 y(0) = 1：y0 = 1;4. 求解差分方程现在我们可以使用欧拉法求解差分方程了。
在MATLAB中，可以使用 odeEuler 函数（需要自己编写）：[t,y] = odeEuler(dydt,tspan,y0,h);5. 可视化结果最后，我们可以将结果可视化，以便更好地理解求解过程。
在这个例子中，我们可以用 plot 函数将求解结果绘制出来：plot(t,y)xlabel('Time')ylabel('y(t)')title('Solution of y'' = -0.5y + 3sin(t) using Euler''s method')以上就是使用欧拉法在MATLAB中求解常微分方程的基本步骤。
matlab求解常微分方程
matlab求解常微分⽅程本⽂主要介绍matlab中求解常微分⽅程（组）的dsolve和ode系列函数，并通过例⼦加深读者的理解。
⼀、符号介绍D: 微分符号；D2表⽰⼆阶微分，D3表⽰三阶微分，以此类推。
⼆、函数功能介绍及例程1、dsolve 函数dsolve函数⽤于求常微分⽅程组的精确解，也称为常微分⽅程的符号解。
如果没有初始条件或边界条件，则求出通解；如果有，则求出特解。
1)函数格式Y = dsolve(‘eq1,eq2,…’ , ’cond1,cond2,…’ , ’Name’)其中，‘eq1,eq2,…’:表⽰微分⽅程或微分⽅程组;’cond1,cond2,…’:表⽰初始条件或边界条件;‘Name’:表⽰变量。
没有指定变量时，matlab默认的变量为t；2)例程例1.1(dsolve 求解微分⽅程)求解微分⽅程：dsolve('Dy=3*x^2','x')例1.2（加上初始条件）求解微分⽅程：例2(dsolve 求解微分⽅程组)求解微分⽅程组：由于x,y均为t的导数，所以不需要在末尾添加’t’。
2、ode函数在上⽂中我们介绍了dsolve函数。
但有⼤量的常微分⽅程，虽然从理论上讲，其解是存在的，但我们却⽆法求出其解析解，此时，我们需要寻求⽅程的数值解。
ode是Matlab专门⽤于解微分⽅程的功能函数。
该求解器有变步长（variable-step）和定步长（fixed-step）两种类型。
不同类型有着不同的求解器，具体说明如下图。
其中，ode45求解器属于变步长的⼀种，采⽤Runge-Kutta算法；其他采⽤相同算法的变步长求解器还有ode23。
ode45表⽰采⽤四阶-五阶Runge-Kutta算法，它⽤4阶⽅法提供候选解，5阶⽅法控制误差，是⼀种⾃适应步长（变步长）的常微分⽅程数值解法，其整体截断误差为(Δx)^5。
解决的是Nonstiff(⾮刚性)常微分⽅程。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
下载的，感觉不错，共享一下
常微分方程和常微分方程组的求解
一、实验目的：
熟悉Matlab 软件中关于求解常微分方程和常微分方程组的各种命令，掌握利用Matlab 软件进行常微分方程和常微分方程组的求解。
二、相关知识
在MATLAB 中，由函数dsolve()解决常微分方程（组）的求解问题，其具体格式如下：
X=dsolve(‘eqn1’,’eqn2’,…)
函数dsolve 用来解符号常微分方程、方程组，如果没有初始条件，则求出通解，如果有初始条件，则求出特解。
例1：求解常微分方程1dy dx x y =
+的MATLAB 程序为：dsolve('Dy=1/(x+y)','x')，
注意，系统缺省的自变量为t ，因此这里要把自变量写明。
结果为：-lambertw(-C1*exp(-x-1))-x-1
其中：Y=lambertw(X)表示函数关系Y*exp(Y)=X 。
例2：求解常微分方程2
'''0yy y -=的MATLAB 程序为：
Y2=dsolve('y*D2y-Dy^2=0’,’x’) 结果为：
Y2 =[ exp((x+C2)/C1)]
[ C2]
我们看到有两个解，其中一个是常数。
例3：求常微分方程组253t
t dx x y e dt dy x y e dt ⎧++=⎪⎪⎨
⎪--=⎪⎩通解的MATLAB 程序为： [X,Y]=dsolve('Dx+5*x+y=exp(t),Dy-x-3*y=exp(2*t)','t')
例4：求常微分方程组020
210cos ,224,0
t t t dx dy x t x dt dt dx dy y e y dt dt =-=⎧+-==⎪⎪⎨
⎪++==⎪⎩通解的MATLAB 程序为：
[X,Y]=dsolve('Dx+2*x-Dy=10*cos(t),Dx+Dy+2*y=4*exp(-2*t)','x(0)=2','y(0)=0')
以上这些都是常微分方程的精确解法，也称为常微分方程的符号解。
但是，我们知道，有大量的常微分方程虽然从理论上讲，其解是存在的，但我们却无法求出其解析解，此时，我们需要寻求方程的数值解，在求常微分方程数值解方面，MATLAB 具有丰富的函数，我们将其统称为solver ，其一般格式为：
[T,Y]=solver(odefun,tspan,y0)
该函数表示在区间tspan=[t0,tf]上，用初始条件y0求解显式常微分方程
'(,)y f t y =。
solver 为命令ode45，ode23，ode113，ode15s ，ode23s ，ode23t ，ode23tb 之一，这些命令各有特点。
我们列表说明如下： 求解器 ODE 类型
特点
说明
ode45
非刚性
一步算法，4,5阶Runge-Kutta 方法累积截断误差3
()x ∆
大部分场合的首选算
法
ode23 非刚性
一步算法，2,3阶Runge-Kutta 方法累积截断误差3()x ∆
使用于精度较低的情
形
ode113
非刚性
多步法，Adams 算法，高低精
度均可达到3
6
10~10--
计算时间比ode45短 ode23t 适度刚性 采用梯形算法 适度刚性情形 ode15s 刚性 多步法，Gear’s 反向 数值积分，精度中等 若ode45失效时， 可尝试使用 ode23s
刚性
一步法，2阶Rosebrock 算法，
低精度。
当精度较低时， 计算时间比ode15s 短
odefun 为显式常微分方程'(,)y f t y =中的(,)f t y tspan 为求解区间，要获得问题在其他指定点012,,,
t t t 上的解，则令
012[,,,
,]f tspan t t t t =（要求i t 单调），
y0初始条件。
例5：求解常微分方程
2
'222y y x x =-++，00.5x ≤≤，(0)1y =的MATLAB 程序如下：fun=inline('-2*y+2*x*x+2*x')；[x,y]=ode23(fun,[0,0.5],1) 结果为：
x =
0,0.0400,0.0900,0.1400,0.1900,0.2400,0.2900,0.3400,0.3900,0.4400,0.4900,0.5000 y =
1.0000,0.9247,0.8434,0.7754,0.7199,0.6764,0.6440,0.6222,0.6105,0.6084,0.6154,0.6179
例6：求解常微分方程222
(1)0,(0)1,'(0)0d y dy y y y y dt dt μ--+===的解，并画
出解的图形。
分析：这是一个二阶非线性方程，用现成的方法均不能求解，但我们可以通过下面的变换，将二阶方程化为一阶方程组，即可求解。
令：1x y =，
2dy
x dt =
，7μ=，则得到：
接着，编写vdp.m 如下： function fy=vdp(t,x)
fy=[x(2);7*(1-x(1)^2)*x(2)-x(1)]; 再编写m 文件sy12_6.m 如下： y0=[1;0]
[t,x]=ode45(@vdp,[0,40],y0); y=x(:,1);dy=x(:,2); plot(t,y,t,dy) 三、实验内容
1．利用MATLAB 求常微分方程的初值问题38dy
y dx +=，02x y ==的解。
2．利用MATLAB 求常微分方程的初值问题2
(1)''2'x y xy +=，01x y ==，0'3
x y ==的解。
3．利用MATLAB 求常微分方程(4)
2'''''0y
y y -+=的解。
4．利用MATLAB 求常微分方程组003
24,2
30,0
t
t t dx dy x y e x dt dt dx x y y dt ==⎧++-==
⎪⎪⎨
⎪++==⎪
⎩的特解。
5．求解常微分方程2
''2(1)'0y y y y --+=，030x ≤≤，(0)1y =，'(0)0y =的特解，
并做出解函数的曲线图。
6．完成实验报告。