从系列:在MATLAB中求解ode
克里夫硅藻土,MathWorks
ODE23通过比较二阶和三阶方法,自动选择步长并保持一定的精度。它是最简单的MATLAB求解器,具有自动误差估计和连续插值等现代特征。ODE23适用于计算机图形等粗精度要求。
实现现代数值方法的软件有两个在ODE4和古典龙格-库塔(Runge-Kutta)等代码中没有的特性。软件中的方法可以估计误差,并提供自动步长控制。您不需要指定步长h。您需要指定您想要的精度。该方法对误差进行估计,并相应地调整步长。它们提供了一个完全精确的连续插值。它们不只是提供离散点集的解。它们提供了一个函数,定义了区间内所有地方的解。你可以绘制它,找到函数的0,提供一个叫做事件处理的工具,等等。
Larry Shampine是常微分方程数值解方面的权威。他是这本用MATLAB求解ode的教科书的主要作者。他现在是达拉斯南卫理公会大学的名誉教授。他长期担任MathWorks的顾问,负责ODE套件的开发。Shampine和他的学生Przemyslaw Bogacki在1989年发表了这种方法。它是ODE23的基础,我们将在MATLAB ODE套件中使用的第一个方法。
基本的方法是三阶。误差估计是基于三阶法和二阶法之间的差异。有四个斜坡。
第一个是函数在区间开始处的值。但这是基于FSAL的,第一步和最后一步一样,斜率很可能是前一步剩下的。如果上一步成功,则该函数值与上一步的最后一个函数值相同。
这个斜率用于步进到区间的中间,函数在这里求值。这个斜率用来跨越区间的3/4然后得到第三个斜率。然后用这三个值来执行步骤。Yn + 1是这三个函数值的线性组合。然后对函数求值,在区间的末尾得到第四个斜率。然后,用这四个斜率来估计误差。
这里的误差估计值是yn + 1和另一个解的估计值之间的差这个解是通过二阶方法得到的我们实际上不计算。我们只需要这个方法和yn + 1的差来估计误差。
这个估计误差与用户提供的公差进行比较。如果估计的误差小于一个公差,那么这一步就是成功的。第四个斜率,s4,变成了下一步的s1。
如果结果大于公差,则误差可作为调整步长的依据。在任何一种情况下,误差估计都是调整下一步步长的基础。这是boacki - shampine阶3(2)方法,它是ODE 23的基础。
首先,让我们看看ODE23的一些非常简单的用法。我要取微分方程y ' = y,我要计算e ^ t,调用ODE23在0到1的区间上,初始值为1。没有输出参数。如果我叫它ODE23,它只是画出了解。在这儿。它只是生成了一个图。它选择一个步长,从0到1,然后得到e的最终值,2.7左右。
如果我提供输出参数。我说(t, y) = ODE23,它返回的是t和y的值,ODE23选择它想要的t的值。这是一个微不足道的问题。它最终选择的步长为0.1。在它开始之后,它选择一个初始步长为。08的容错值。y的最终值是2.718,也就是e的值。
这就是ODE23的两种简单用法。如果不提供任何输出参数,它会绘制一个图。如果你提供输出参数,t和y,它会返回t和y的值选择t的值来满足错误。默认的误差是10的- 3次方。所以这个值精确到三位数。果然,这就是我们得到的。
现在让我们尝试一些更有挑战性的东西,看看自动错误控制步长选择是如何运行的。设a等于四分之一。然后令y0 = 15.9。如果设它为16,也就是1 / a²,就会遇到奇点。微分方程是y ' = 2 (a - t)乘以y²。我要把它和ODE23积分在0到1的区间上从y0开始,把结果保存在t和y中,然后画出来。这是我的绘图命令,这是解。
在a点有一个近似奇点,它几乎要爆炸了。然后它又平静下来。所以当你上升到奇点,然后下降到奇点时,这些点聚在一起,但随着解的稳定,它们之间的距离会越来越远。ODE解算器可以采取更大的步骤。
为了看看实际采取了什么步骤,我们来计算t的差,然后画出来。这是所取的步长。我们可以看到,在几乎是奇点的0.25处,我们取了一个小的步长。当我们接近这个区间的末尾时,步长会变大。最后,到达区间末端的步长作为最后一步。这就是ODE23的自动步长选择。
BS23有一个很好的自然插值它实际上已经被知道100多年了。它被称为埃尔米特三次插值。我们知道两个点决定一条直线。两个点和两个斜率决定了一个三次方程。
在每个区间都有y和yn + 1的值。也有两个斜率,就是这个。我们有端点的导数,yn '和yn + 1 ',这是微分方程在这些点上的值。这四个值决定了一个经过这两点的三次曲线并且有这两个斜率。
这个立方允许评估解决方案的软件在任何时候在定义的时间间隔没有额外的成本函数f的评估。这可以用来画图形解决方案的,平滑的图形解决方案,找到解决方案的0,事件处理,等等。ODE23提供的另一个特性。