第12章 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_常微分方程数值解法
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提供了多种函数用于求解常微分方程(组),其中最常用的函数是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是一种用于科学计算和工程应用的高级编程语言和交互式环境。
它在数学建模、模拟和分析等方面有着广泛的应用。
在MATLAB 中,常微分方程的数值求解是一个常见的应用场景。
在实际工程问题中,通常需要对常微分方程进行数值求解来模拟系统的动态行为。
本文将介绍MATLAB中对常微分方程进行数值求解的快速方法。
1. 基本概念在MATLAB中,可以使用ode45函数来对常微分方程进行数值求解。
ode45是一种常用的Runge-Kutta法,它可以自适应地选取步长,并且具有较高的数值精度。
使用ode45函数可以方便地对各种类型的常微分方程进行求解,包括一阶、高阶、常系数和变系数的微分方程。
2. 函数调用要使用ode45函数进行常微分方程的数值求解,需要按照以下格式进行函数调用:[t, y] = ode45(odefun, tspan, y0)其中,odefun表示用于描述微分方程的函数,tspan表示求解的时间跨度,y0表示初值条件,t和y分别表示求解得到的时间序列和对应的解向量。
3. 示例演示为了更好地理解如何使用ode45函数进行常微分方程的数值求解,下面我们以一个具体的例子来进行演示。
考虑如下的一阶常微分方程:dy/dt = -2*y其中,y(0) = 1。
我们可以编写一个描述微分方程的函数odefun:function dydt = odefun(t, y)dydt = -2*y;按照上述的函数调用格式,使用ode45函数进行求解:tspan = [0 10];y0 = 1;[t, y] = ode45(odefun, tspan, y0);绘制出解曲线:plot(t, y);4. 高级用法除了基本的函数调用方式外,MATLAB中还提供了更多高级的方法来对常微分方程进行数值求解。
可以通过设定选项参数来控制数值求解的精度和稳定性,并且还可以对刚性微分方程进行求解。
5. 性能优化在实际工程应用中,常常需要对大规模的常微分方程进行数值求解。
matlab求解微分方程数值解与解析解
matlab求解微分方程数值解与解析解微分方程是数学中的重要内容,它描述了物理、工程、经济等领域中的许多现象和问题。
在实际应用中,我们经常需要求解微分方程的解析解或数值解。
本文将以Matlab为工具,探讨如何求解微分方程并对比解析解与数值解的差异。
一、引言微分方程是描述自然界中许多现象和问题的数学语言,它包含了未知函数及其导数与自变量之间的关系。
微分方程的求解可以帮助我们了解问题的性质和变化规律,并为实际应用提供参考。
在许多情况下,微分方程的解析解很难求得,这时我们可以利用计算机进行数值求解。
二、微分方程的数值解法1.欧拉法欧拉法是最简单的数值求解微分方程的方法之一。
它通过将微分方程转化为差分方程,然后利用离散的点逼近连续的解。
具体步骤如下:(1)将微分方程转化为差分方程,即用近似的导数代替真实的导数;(2)选择初始条件,即确定初始点的值;(3)选择步长和求解区间,即确定求解的范围和步长;(4)使用迭代公式计算下一个点的值;(5)重复步骤(4),直到达到指定的求解区间。
2.改进的欧拉法欧拉法存在精度较低的问题,为了提高精度,可以使用改进的欧拉法。
改进的欧拉法是通过使用两次导数的平均值来计算下一个点的值,从而提高了数值解的精度。
3.龙格-库塔法龙格-库塔法是一种常用的数值求解微分方程的方法,它通过使用多个点的导数来逼近连续解。
龙格-库塔法的步骤如下:(1)选择初始条件和步长;(2)使用迭代公式计算下一个点的值;(3)计算下一个点的导数;(4)根据导数的值和步长计算下一个点的值;(5)重复步骤(3)和(4),直到达到指定的求解区间。
龙格-库塔法的精度较高,适用于求解一阶和高阶微分方程。
三、微分方程的解析解解析解是指能够用公式或函数表示的方程的解。
有些微分方程具有解析解,可以通过数学方法求得。
例如,一阶线性常微分方程和某些特殊类型的二阶微分方程等。
解析解的优势在于精确性和直观性,能够帮助我们深入理解问题的本质。
常微分方程(ODEs)的MATLAB数值解法
Matlab 中常微分方程数值解法讲解©
作者:dynamic 时间:2008.12.10 版权:All Rights Reserved By
★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★★ Matlab Sky 联盟----打造最优秀、专业和权威的 Matlab 技术交流平台! 网址: /com/org/net 邮箱:matlabsky@ QQ 群:23830382 40510634 16233891(满了)
然而现实总是残酷的,有时方程偏偏是隐式的,没法写成上面的样子,不用担心 Matlab 早就为我们想到了,这 个留在后面的隐式微分方程数值求解中再详细讲解!
step2.为每一阶微分式选择状态变量,最高阶除外
x1 = x, x2 = x ', x3 = x '', L, xm = x ( m −1) xm+1 = y, xm+1 = y ', xm + 2 = y '',L , xm + n = y ( n−1)
4
打造最优秀、专业和权威的 Matlab 技术交流平台!
2.微分方程转换
好, 上面我们把 Matlab 中的常微分方程(ODE)的解算器讲解的差不多了, 下面我们就具体开始介绍如何使用上面 的知识吧! 现实总是残酷的,要得到就必须先付出,不可能所有的 ODE 一拿来就可以直接使用,因此,在使用 ODE 解算 器之前,我们需要做的第一步,也是最重要的一步就是将微分方程组化成 Matlab 可接受的标准形式! 如果 ODEs 由一个或多个高阶微分方程给出,则我们应先将其变换成一阶显式常微分方程组! 下面我们以两个高阶微分方程构成的 ODEs 为例介绍如何将之变换成一个一阶显式常微分方程组。
matlab求微分方程组的解析解
matlab求微分方程组的解析解(实用版)目录1.引言2.MATLAB 求微分方程组的解析解的方法3.示例:求解一阶微分方程组4.示例:求解二阶微分方程组5.总结正文一、引言微分方程组在数学建模和实际问题中有着广泛的应用,求解微分方程组对于理解问题的内在机制和预测未来发展趋势具有重要意义。
在众多数学软件中,MATLAB 凭借其强大的数值计算和图形绘制功能,成为求解微分方程组的常用工具。
本文将介绍如何使用 MATLAB 求解微分方程组的解析解。
二、MATLAB 求微分方程组的解析解的方法MATLAB 求解微分方程组的解析解主要依赖于符号计算函数和数值计算函数。
其中,符号计算函数主要用于求解微分方程组的解析解,数值计算函数则用于求解微分方程组的数值解。
在使用这些函数时,需要确保符号计算和数值计算的顺序,以避免计算错误。
三、示例:求解一阶微分方程组考虑如下一阶微分方程组:```dy/dx = x + ydz/dx = x - z```我们可以使用 MATLAB 的符号计算函数`symfun`和`symvar`来求解该方程组。
首先,定义符号变量 x、y、z 和 p(表示参数),然后使用`symfun`函数创建微分方程组的符号表达式。
接着,利用`symvar`函数求解微分方程组,并将结果转换为数值形式。
最后,使用`plot`函数绘制解的图形。
四、示例:求解二阶微分方程组考虑如下二阶微分方程组:```x"" + 3x" + 2x = 0y"" + 3y" + 2y = 0```我们可以使用 MATLAB 的符号计算函数`symfun`和`symvar`来求解该方程组。
首先,定义符号变量 x、y 和 p(表示参数),然后使用`symfun`函数创建微分方程组的符号表达式。
接着,利用`symvar`函数求解微分方程组,并将结果转换为数值形式。
最后,使用`plot`函数绘制解的图形。
matlab求微分方程组的解析解
MATLAB求微分方程组的解析解引言在科学与工程领域,微分方程组是一种常见的数学模型,用于描述各种物理现象和工程问题。
解析解是指能够用公式表达出来的精确解。
本文将介绍如何使用M ATL A B求解微分方程组的解析解。
1.方程组的建立首先,我们需要确定待求解的微分方程组。
假设我们有一个由n个微分方程组成的方程组,可以写为如下形式:d y1/dt=f1(t,y1,y2,...,yn)d y2/dt=f2(t,y1,y2,...,yn)......d y n/dt=f n(t,y1,y2,...,yn)其中`t`是自变量,`y1,y2,...,y n`是因变量,`f1,f2,...,fn`是给定的函数关系。
我们的目标是求解`y1(t),y2(t),...,yn(t)`的解析解。
2.使用MAT LAB求解M A TL AB提供了强大的求解微分方程组的工具,我们可以使用其中的函数来求解上述方程组的解析解。
首先,我们需要在MA T LA B中定义方程组的函数形式。
可以通过定义一个函数或者使用匿名函数来实现。
例如,我们可以定义一个名为`m yE qu at io ns`的函数,其输入参数为`t`和一个向量`y`,输出为一个向量`d y`,代表方程组的左侧和右侧的变量分别。
函数示例如下:f u nc ti on dy=m yE qua t io ns(t,y)%定义方程组d y=z er os(n,1);d y(1)=f1(t,y(1),y(2),...,y(n));d y(2)=f2(t,y(1),y(2),...,y(n));......d y(n)=fn(t,y(1),y(2),...,y(n));e n d然后,我们可以使用M AT LA B的`d so lv e`函数来求解微分方程组的解析解。
示例如下:s y ms ty1(t)y2(t)...yn(t)a s su me(t,'re al')a s su me(y1(t),'rea l')a s su me(y2(t),'rea l')......a s su me(y n(t),'rea l')e q n1=d if f(y1(t),t)==f1(t,y1(t),y2(t),...,y n(t));e q n2=d if f(y2(t),t)==f2(t,y1(t),y2(t),...,y n(t));......e q nn=d if f(yn(t),t)==fn(t,y1(t),y2(t),...,y n(t));e q ns=[eq n1,e qn2,...,eq nn];S=ds ol ve(e qn s);`S`即为方程组的解析解集合。
matlab 求微分方程组数值解
matlab 求微分方程组数值解使用Matlab求解微分方程组是一种常见的数值方法。
微分方程组是描述自然界中许多现象的数学模型,它们可以用一组关于未知函数及其导数的方程来表示。
通过求解微分方程组,我们可以得到未知函数在给定条件下的数值解。
在Matlab中,求解微分方程组可以使用ode45函数。
该函数是一个常用的求解常微分方程初值问题的函数,它使用四阶龙格-库塔法(RK4)进行数值求解。
使用ode45函数求解微分方程组的步骤如下:定义微分方程组。
在Matlab中,可以使用匿名函数或函数句柄的方式定义微分方程组。
例如,对于一个二阶微分方程组:dy1/dt = f1(t, y1, y2)dy2/dt = f2(t, y1, y2)可以定义一个匿名函数:f = @(t, y) [f1(t, y(1), y(2)); f2(t, y(1), y(2))]其中,t是自变量,y是未知函数的向量。
接下来,指定求解的时间区间和初值条件。
时间区间可以通过指定起始时间和结束时间来确定。
初值条件是指在起始时间处未知函数的值。
初值条件可以通过一个向量来表示。
例如,对于一个二阶微分方程组,初值条件可以表示为一个长度为2的向量。
然后,调用ode45函数进行求解。
ode45函数的输入参数包括定义的微分方程组、时间区间和初值条件。
该函数会返回数值解和对应的时间点。
可以通过绘制图形或打印数值解来展示结果。
Matlab提供了丰富的绘图函数,可以方便地将数值解可视化。
需要注意的是,求解微分方程组时,应选择合适的数值方法和步长,以保证数值解的精度和稳定性。
对于复杂的微分方程组,可能需要进行参数调整和迭代求解,以得到满意的结果。
使用Matlab求解微分方程组是一种便捷而有效的数值方法。
通过定义微分方程组、指定时间区间和初值条件,调用ode45函数进行求解,可以得到微分方程组的数值解。
这种方法在科学研究和工程实践中具有广泛的应用,可以帮助我们更好地理解和分析自然界中的现象。
Matlab学习——求解微分方程(组)
Matlab学习——求解微分⽅程(组)介绍:1.在 Matlab 中,⽤⼤写字母 D 表⽰导数,Dy 表⽰ y 关于⾃变量的⼀阶导数,D2y 表⽰ y 关于⾃变量的⼆阶导数,依此类推.函数 dsolve ⽤来解决常微分⽅程(组)的求解问题,调⽤格式为X=dsolve(‘eqn1’,’eqn2’,…)如果没有初始条件,则求出通解,如果有初始条件,则求出特解系统缺省的⾃变量为 t。
2.函数 dsolve 求解的是常微分⽅程的精确解法,也称为常微分⽅程的符号解.但是,有⼤量的常微分⽅程虽然从理论上讲,其解是存在的,但我们却⽆法求出其解析解,此时,我们需要寻求⽅程的数值解,在求常微分⽅程数值解⽅⾯,MATLAB 具有丰富的函数,将其统称为 solver,其⼀般格式为:[T,Y]=solver(odefun,tspan,y0)说明:(1)solver 为命令 ode45、ode23、ode113、ode15s、ode23s、ode23t、ode23tb、ode15i 之⼀.(2)odefun 是显⽰微分⽅程y ' = f (t , y) 在积分区间 tspan = [t 0 , t f ] 上从t0 到t f⽤初始条件 y0求解.(3)如果要获得微分⽅程问题在其他指定时间点t 0 , t1 , t 2 , , t f上的解,则令tspan = [t 0 , t1 , t 2 , t f ](要求是单调的).(4)因为没有⼀种算法可以有效的解决所有的 ODE 问题,为此,Matlab 提供了多种求解器 solver,对于不同的 ODE 问题,采⽤不同的 solver3.在 matlab 命令窗⼝、程序或函数中创建局部函数时,可⽤内联函数 inline,inline 函数形式相当于编写 M 函数⽂件,但不需编写 M-⽂件就可以描述出某种数学关系.调⽤ inline 函数,只能由⼀个 matlab 表达式组成,并且只能返回⼀个变量,不允许[u,v]这种向量形式.因⽽,任何要求逻辑运算或乘法运算以求得最终结果的场合,都不能应⽤ inline 函数,inline 函数的⼀般形式为:FunctionName=inline(‘函数内容’, ‘所有⾃变量列表’)例如:(求解 F(x)=x^2*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 来定义函数。
第12章 MATLAB常微分方程(组)数值求解方程与方程组的数值解.
2. 边值问题
两个求解函数函数bvp4c和bvp5c,后者求解精度要比前者好。以 bvpsolver表示bvp4c或者bvp5c,那么这两个函数有着统一的调用格 式:
solinit = bvpinit(x, yinit, params) sol = bvpsolver(odefun,bcfun,solinit) sol = bvpsolver(odefun,bcfun,solinit,options)
-2.5
-6 -3
-2
-1
0 位移
1
2
3
2018/10/9
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
三、 刚性问题举例
问题见书中【例12.2-2】, 【例12.2-3】。下图是【例12.2-2】不同 求解器得到解的图像对比。
子 函 数 形 式 /ode23 100 90 80 70 60 50 40 30 20 10 0 子 函 数 形 式 /ode15s 匿 名 函 数 形 式 /ode15s 100 100 90 80 70 60 50 40 30 20 10 0 90 80 70 60 50 40 30 20 10 0
三、 利用fzero/fsolve函数
问题见书中【例12.3-2】, 【例12.3-3】, 【例12.3-4】 。下图是 【例12.3-2】结果图像。
1.8
1.6
1.4
1.2
1
0.8
0.6
0.4
0.2
0 0 2 4 6 8 10 12 14 16 18 20
t
2018/10/9
©
吴鹏, MATLAB从零到进阶.
2018/10/9
©
MATLAB求解常微分方程数值解
MATLAB求解常微分方程数值解利用MATLAB求解常微分方程数值解1.内容简介把《高等工程数学》看了一遍,增加对数学内容的了解,对其中数值解法比较感兴趣,这大概是因为在其它各方面的学习和研究中经常会遇到数值解法的问题。
理解模型然后列出微分方程,却对着方程无从下手,无法得出精确结果实在是让人难受的一件事情。
实际问题中更多遇到的是利用数值法求解偏微分方程问题,但考虑到先从常微分方程下手更为简单有效率,所以本文只研究常微分方程的数值解法。
把一个工程实际问题弄出精确结果远比弄清楚各种细枝末节更有意思,因此文章中不追求非常严格地证明,而是偏向如何利用工具实际求解出常微分方程的数值解,力求将课程上所学的知识真正地运用到实际方程的求解中去,在以后遇到微分方程的时候能够熟练运用MATLAB得到能够在工程上运用的结果。
文中求解过程中用到MATLAB进行数值求解,主要目的是弄清楚各个函数本质上是如何对常1微分方程进行求解的,对各种方法进行MATLAB 编程求解,并将求得的数值解与精确解对比,其中源程序在附录中。
最后考察MATLAB中各个函数的适用范围,当遇到实际工程问题时能够正确地得到问题的数值解。
2.Euler Method(欧拉法)求解Euler法求解常微分方程主要包括3种形式,即显式Euler法、隐式Euler法、梯形公式法,本节内容分别介绍这3种方法的具体内容,并在最后对3种方法精度进行对比,讨论Euler法的实用性。
本节考虑实际初值问题使用解析法,对方程两边同乘以得到下式两边同时求积分并采用分部积分得到解析2解:本节后面将对此方程进行求解,并与精确解进行对比,分析Euler的可行性。
2.1.显式Euler法和隐式Euler法显式和隐式Euler法都属于一阶方法,显式Euler法的迭代公式简单,如下所示:对过上述公式对式进行迭代,其中步长3图2.1 显式Euler法精确解和数值解图像从图2.1中可以看出,显式Euler法在斜率很大的时候存在非常大的误差。
实验二 MATLAB数值计算_常微分方程(组)的求解
y ( n ) = f (t , y' , y" , ⋯, y ( n −1 ) )
(2)
设 y1 = y, y2 = y ', y3 = y '',...., yn = y( n −1) ,可将上式化为 n 元一阶常微分方程组:
⎧ y1 ' = y 2 ⎪ y '= y 2 3 ⎪ ⎪ ⎨⋯ ⎪ y' = y n ⎪ n −1 ⎪ y ' = f ( t , y1 , y 2 ,⋯ , y n ) ⎩ n
方程(4)式要有确定的解必须给定初值条件(t0 为初始时刻):
⎛ y10 ⎞ ⎜y ⎟ y0 = y ( t0 ) = ⎜ 20 ⎟ (7) ⎜ ⋮ ⎟ ⎟ ⎜ ⎝ y n0 ⎠
所谓数值解法, 就是在求解区间[t0 ,tf ]上寻求 y (t)在一系列离散节点 t 0 < t 1 < t 2 < ⋯ < t m ≤ t f 上的近似值 y ( t k ), k = 0,1,..., m 。称 hk = t k +1 − t k 为步长,通常取为常量 h。最简单的数 值解法是 Euler 方法。 Euler 法的思路极其简单:在节点处用差商近似代替(4)式中的导数
值,即 f1 , f2 ,…, fn 。故 odefun 函数的返回值一般是列向量,其最简单的编写格式为: F = odefun (t, y) t y F 为标量(自变量) ,代表迭代点(迭代的当前时刻,如 t0 , t1 等) 列向量(见 5 式), 因计算 t 时刻的函数表达式 F 需要知道该时刻的 y 值, 见式(8) 返回值,代表 t 时刻微分方程组右边的函数值,为列向量(见 6 式) 。有了向量 F 后,若再给出初始条件,根据(8)或(9)式就可迭代出方程的解。 求解器在计算时将会不断地在各个时刻调用 odefun 的值。此外,MATLAB 还专门提供了一 个 ode 函数的编写模板,给出了该函数的详细编写格式,可使用命令 : open odefile.m 打开, 附录 2 给出了 ode 函数模板的中文说明。 tspan 指定积分迭代的区间 [t0 , tf ],t0 是迭代的初始时刻。 (i). 若 tspan = [ t0, tf] 为 2
matlab常微分方程的数值解
plot(z(:,1),z(:,2),'m-')
-1
0
1
2
3
4
5
6
第5页/共8页
x 2ty 0
例4
y
2tx
0
H=[-30,30]
t 0,[x0, y0 ] [0,0],[x0 , y0 ] [1,0]
(1)转化方程。
(2)建立文件,设文件名为dzdt04.m。
function dz=f(t,z) dz(1)=z(2);
odeij ode45 ode23
问题类型 非刚性 非刚性
精度 中等 较低
ode113 非刚性
低到高
ode15s 刚性
低到高
ode23s 刚性
低
ode23t 适度刚性 低
ode23tb 刚性
低
适用对象
多数情况下可优先选用,但不能用来 解刚性问题 可用来解中等刚性问题,或误差允许 范围较宽的问题
不能解刚性问题,当误差容限要求严 格时效果较ode45好
(1)转化方程。
(2)建立文件,设文件名为dzdt03.m。
function dz=f(t,z)
2
a=-0.2;
1.5
dz(1)=z(2);
1
dz(2)=a*z(2)-sin(t);
dz=[dz(1);dz(2)];
0.5
(3)调用求解。
0
H=[0,100];z0=[0;2]; -0.5
[t,z]=ode45('dzdt03',H,z0);
可用于解刚性问题,当采用ode45失败 或效果很差时,可考虑使用
可用于解刚性问题,当误差容限较宽 时效果比ode15s好
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
如果质量矩阵奇异的话,(1)称为微分代数方程组(differential algebraic equations, DAEs.),可以利用求解刚性微分方程的函数如ode15s,ode23s 等来求解,从输入形式上看,求解DAEs和求解普通的ODE很类似,主 要区别是需要给微分方程求解器指定质量矩阵。
2018/10/9
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
二、初值问题求解函数
y0 :微分方程(组)的初值,即所有状态变量在t0时刻的值。 options 结构体,通过odeset设置得到的微分优化参数。
返回参数说明: T:时间点组成的列向量 Y:微分方程(组)的解矩阵,每一行对应相应T的该行上时间点的微 分方程(组)的解。 sol:以结构体的形式返回解。
常微分方程数值解
常微分方程(组) 数值求解
吴鹏(rocwoods)
rocwoods@
MATLAB从零到进阶
2018/10/9
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
主要内容
数值求解常微分方程(组)函数概述 非刚性/刚性常微分方程问题求解 隐式微分方程(组)求解 微分代数方程(DAE)与延迟微分方程(DDE) 求解 边值问题求解
ode113 非刚性
低到高
2018/10/9
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
二、初值问题求解函数
ode15s 刚性 低到中 基于数值差分公式(后向差分公式,BDFs 也叫Gear方法),因此效率不是很高。同 ode113一样,ode15s也是一个多步计算器。 当ode45求解失败,或者非常慢,并且怀 疑问题是刚性的,或者求解微分代数问题 时可以考虑用ode15s 低 基于修正的二阶Rosenbrock公式。由于是 单步解算器,当精度要求不高时,它效率 可能会高于ode15s。它可以解决一些 ode15s求解起来效率不太高的刚性问题。 ode23t 适度刚性 低 低 ode23t可以用来求解微分代数方程。 当方程是刚性的,并且求解要求精度不高 时可以使用。
2018/10/9
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
二、初值问题求解函数
2. 函数介绍
函数 ode45 ode23 问题类型 精确度 非刚性 非刚性 中等 低 说明 采用算法为4-5阶Runge-Kutta法,大多数 情况下首选的函数 基于 Bogacki-Shampine 2-3阶Runge-Kutta 公式,在精度要求不高的场合,以及对于 轻度刚性方程,ode23的效率可能好于 ode45。 基于变阶次Adams-Bashforth-Moutlon PECE算法。在对误差要求严格的场合或 者输入参数odefun代表的函数本身计算量 很大情况下比ode45效率高。ode113可以看 成一个多步解算器,因为它会利用前几次 时间节点上的解计算当前时间节点的解。 因此它不适应于非连续系统。
ode23s
刚性
ode23tb 刚性
2018/10
三、 延迟问题以及边值问题求解函数
1. 延迟问题
MATLAB提供了dde23和ddesd函数用来求解。前者用来求解状态变量 延迟为常数的微分方程(组),后者用来求解状态变量延迟不为常数 的微分方程(组)。调用格式以及参数意义大部分类似ode系列求解函 数,不同的是要输入延迟参数以及系统在时间小于初值时的状态函数。
2. 边值问题
两个求解函数函数bvp4c和bvp5c,后者求解精度要比前者好。以 bvpsolver表示bvp4c或者bvp5c,那么这两个函数有着统一的调用格 式:
solinit = bvpinit(x, yinit, params) sol = bvpsolver(odefun,bcfun,solinit) sol = bvpsolver(odefun,bcfun,solinit,options)
2018/10/9
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
第一节数值求解常微分方程(组) 函数概述
2018/10/9
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
一、 概述
第9章介绍了符号求解各类型的微分方程组,但是 能够求得解析解的微分方程往往只是出现在大学课 堂上,在实际应用中,绝大多数微分方程(组)无 法求得解析解。这就需要利用数值方法求解。 MATLAB以数值计算见长,提供了一系列数值求解 微分方程的函数。 这些函数可以求解非刚性问题,刚性问题,隐式 微分方程,微分代数方程等初值问题,也可以求解 延迟微分方程,以及边值问题等。
2018/10/9
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
二、初值问题求解函数
1. 提供的函数
ode23, ode45, ode113, ode15s, ode23s, ode23t, ode23tb,这些函数统一 的调用格式如下: [T,Y] = solver(odefun,tspan,y0) [T,Y] = solver(odefun,tspan,y0,options) sol = solver(odefun,[t0 tf],y0...) 输入参数说明: odefun 表示微分方程(组)的句柄。 tspan 微分方程(组)的求解时间区间,有两种格式[t0,tf]或者 [t0,t1,…,tf],两者都以t0为初值点,根据tf自动选择积分步长。前者 返回实际求解过程中所有求解的时间点上的解,而后者只返回设定 的时间点上的解。后者对计算效率没有太大影响,但是求解大型问 题时,可以减少内存存储。
2018/10/9
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
四、 求解前准备工作
微分方程的形式是多种多样的,一般来说,很多高阶微分方程可以通过 变量替换转化成一阶微分方程组,即可以写成下面的形式: M t, y y ' F t , y ( 1) M t , y 称为质量矩阵,如果其非奇异的话,上式可以写成: y ' M 1 t, y F t, y ( 2) 将(2)式右半部分用odefun表示出来(具体表现形式可以采用匿名函数、 子函数、嵌套函数、单独m文件等形式),就是ode45,ode23等常微分 方程初值问题求解的输入参数odefun。