第12章 常微分方程(组)数值求解方程与方程组的数值解

合集下载
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
一、 微分代数方程(DAE)举例
【例12.4-1】的结果图:
方 法 1计 算 结 果 图 1 1 方 法 2计 算 结 果 图 1 方 法 3计 算 结 果 图
0.8
0.8
0.8
0.6
0.6
0.6
0.4
0.4
0.4
0.2
0.2
0.2
2. 边值问题
两个求解函数函数bvp4c和bvp5c,后者求解精度要比前者好。以 bvpsolver表示bvp4c或者bvp5c,那么这两个函数有着统一的调用格 式:
solinit = bvpinit(x, yinit, params) sol = bvpsolver(odefun,bcfun,solinit) sol = bvpsolver(odefun,bcfun,solinit,options)
0
5
10
15
20
0
5
10
15
20
t
2017/7/20
t
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
第三节 隐式微分方程(组)求解
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
一、 概述
一些 微分方程组在初始给出的时候是不容易显示的 表示成上面提到的标准形式的。这时候就需要想办法表 示成上述的形式。一般来说有三种思路,一种是利用 solve函数符号求解出高阶微分的显式表达式,一种是利 用fzero/fsolve函数求解状态变量的微分值,还有一种是 利用MATLAB自带的ode15i函数 。
常微分方程数值解
三、 利用fzero/fsolve函数
下图是【例12.3-3】结果图像。 【例12.3-4】是利用ode15i求解【例12.33】算例,速度明显增快,结果一致,图像也是下图。
45 y1 (t) 40 35 30 25 20 15 10 5 0 y2 (t) y3 (t) y4 (t)
ode113 非刚性
低到高
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
二、初值问题求解函数
ode15s 刚性 低到中 基于数值差分公式(后向差分公式,BDFs 也叫Gear方法),因此效率不是很高。同 ode113一样,ode15s也是一个多步计算器。 当ode45求解失败,或者非常慢,并且怀 疑问题是刚性的,或者求解微分代数问题 时可以考虑用ode15s 低 基于修正的二阶Rosenbrock公式。由于是 单步解算器,当精度要求不高时,它效率 可能会高于ode15s。它可以解决一些 ode15s求解起来效率不太高的刚性问题。 ode23t 适度刚性 低 低 ode23t可以用来求解微分代数方程。 当方程是刚性的,并且求解要求精度不高 时可以使用。
如果质量矩阵奇异的话,(1)称为微分代数方程组(differential algebraic equations, DAEs.),可以利用求解刚性微分方程的函数如ode15s,ode23s 等来求解,从输入形式上看,求解DAEs和求解普通的ODE很类似,主 要区别是需要给微分方程求解器指定质量矩阵。
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
第一节数值求解常微分方程(组) 函数概述
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
一、 概述
第9章介绍了符号求解各类型的微分方程组,但是 能够求得解析解的微分方程往往只是出现在大学课 堂上,在实际应用中,绝大多数微分方程(组)无 法求得解析解。这就需要利用数值方法求解。 MATLAB以数值计算见长,提供了一系列数值求解 微分方程的函数。 这些函数可以求解非刚性问题,刚性问题,隐式 微分方程,微分代数方程等初值问题,也可以求解 延迟微分方程,以及边值问题等。
常微分方程数值解
常微分方程(组) 数值求解
吴鹏(rocwoods)
rocwoods@126.com
MATLAB从零到进阶
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
主要内容
数值求解常微分方程(组)函数概述 非刚性/刚性常微分方程问题求解 隐式微分方程(组)求解 微分代数方程(DAE)与延迟微分方程(DDE) 求解 边值问题求解
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
二、 非刚性问题举例
问题见书中【例12.2-1】,左图微分方程的解,右图平面相轨迹
x(t) 2.5 2 1.5 1
6
4
=1 =2 =3
2
0.5 0
速度
0 5 10 15 t 20 25 30
0
-0.5 -1
-2
-1.5
-4
-2
-2.5
-6 -3
-2
-1
0 位移
1
2
3
2017/7/20
©
吴鹏, 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
0
0
0
-0.2
y1 (t) y2 (t) y3 (t)
-0.2
y1 (t) y2 (t) y3 (t)
-0.2
y1 (t) y2 (t) y3 (t)
-0.4
来自百度文库
0
5
10
15
20
-0.4
0
5
10
15
20
-0.4
0
2017/7/20
t
t
©
5
10
15
20
吴鹏 t , MATLAB从零到进阶.
常微分方程数值解
t
2017/7/20
0
5
10
0
5
10
0
5
10
t
2017/7/20
t
t
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
三、 刚性问题举例
下图是【例12.2-3】得到解的图像,以及两个解的和的图像
10 y1 (t) 8 6 4 4 2 2 0 0 -2 -4 -6 -2 -4 -6 y2 (t) 10 8 6 12
F(t) = y1 (t) + y2 (t) y = 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
2017/7/20
©
吴鹏, MATLAB从零到进阶.
2017/7/20
©
吴鹏, 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自动选择积分步长。前者 返回实际求解过程中所有求解的时间点上的解,而后者只返回设定 的时间点上的解。后者对计算效率没有太大影响,但是求解大型问 题时,可以减少内存存储。
0
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
t
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
第四节 微分代数方程(DAE)与延 迟微分方程(DDE)求解
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
一、 微分代数方程(DAE)举例
DAE的求解一般有三种办法,一种是变量替换法,一种是 用ode15s函数还有一种是用12.3节中提到的ode15i函数 【例12.4-1】是利用上述三种方法求解的普通微分代数方 程 。 【例12.4-2】是变量替换后用fsolve函数求解出每一计算 节点的值,然后再调用ode45、ode23tb等函数求解,另一种 方法就是直接利用ode15i函数求解 。
2017/7/20
©
吴鹏, 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。
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
二、初值问题求解函数
2. 函数介绍
函数 ode45 ode23 问题类型 精确度 非刚性 非刚性 中等 低 说明 采用算法为4-5阶Runge-Kutta法,大多数 情况下首选的函数 基于 Bogacki-Shampine 2-3阶Runge-Kutta 公式,在精度要求不高的场合,以及对于 轻度刚性方程,ode23的效率可能好于 ode45。 基于变阶次Adams-Bashforth-Moutlon PECE算法。在对误差要求严格的场合或 者输入参数odefun代表的函数本身计算量 很大情况下比ode45效率高。ode113可以看 成一个多步解算器,因为它会利用前几次 时间节点上的解计算当前时间节点的解。 因此它不适应于非连续系统。
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
二、初值问题求解函数
y0 :微分方程(组)的初值,即所有状态变量在t0时刻的值。 options 结构体,通过odeset设置得到的微分优化参数。
返回参数说明: T:时间点组成的列向量 Y:微分方程(组)的解矩阵,每一行对应相应T的该行上时间点的微 分方程(组)的解。 sol:以结构体的形式返回解。
隐式微分方程无法写成(1)或者(2)的形式,其求解方法本章也有讨论。
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
第二节 非刚性/刚性常微分方程 初值问题求解
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
一、 概述
所谓刚性、非刚性问题最直观的判别方法就是从 解 在某段时间区间内的变化来看。非刚性问题变化相对缓 慢,而刚性问题在某段时间内会发生剧烈变化,即很短 的时间内,解的变化巨大。对于刚性问题不适合用 ode45来求解,如果硬要用ode45来求解的话,达到指定 精度所耗费的时间往往会非常长 。
ode23s
刚性
ode23tb 刚性
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
三、 延迟问题以及边值问题求解函数
1. 延迟问题
MATLAB提供了dde23和ddesd函数用来求解。前者用来求解状态变量 延迟为常数的微分方程(组),后者用来求解状态变量延迟不为常数 的微分方程(组)。调用格式以及参数意义大部分类似ode系列求解函 数,不同的是要输入延迟参数以及系统在时间小于初值时的状态函数。
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
二、利用solve函数
问题见书中【例12.3-1】,下图是求解出的结果曲线
3.5
3
2.5
y 2 y
1
2
(t) (t)
1.5
1
0.5 0 5 10 15 20 25 30
t
2017/7/20
©
吴鹏, MATLAB从零到进阶.
常微分方程数值解
一、 微分代数方程(DAE)举例
【例12.4-2】的结果图:
方 法 1计 算 结 果 图 6 6 方 法 2计 算 结 果 图
4
4
2
2
0
0
-2
-2
-4 y1 (t) -6 y2 (t) y3 (t) y4 (t) -8 0 1 2 3 4 5
-4 y1 (t) -6 y2 (t) y3 (t) y4 (t) -8 0 1 2 3 4 5
相关文档
最新文档