从系列:在MATLAB常微分方程求解
克里夫·莫勒尔,MathWorks公司
ODE4采用了经典的龙格 - 库塔法,这是在过去的100年微分方程使用最广泛的数值方法。它的主要缺点是缺乏一个错误的估计。火焰的增长,一个简单的模型是在这里,并在以后的影片中使用的例子。
这是经典的龙格 - 库塔方法。这是由远及远,超过100年的手工计算,在20世纪上半叶世界上最流行的数值计算方法,然后在20世纪后半叶数字计算机计算。我怀疑今天它仍然在使用。
您评估功能四次每一步,先在区间的开始。然后使用该步入区间,中间以获得S2。然后使用S2步入区间中部再次。
并再次评估功能那里得到S3。然后使用S3跨越间隔步骤清晰,并获得S4。然后取这四个斜坡的组合,加权中间的两个沉重,把你的最后一步。这是经典的龙格 - 库塔方法。
这是我们的MATLAB实现。我们将其称之为ODE4,因为它每一步的计算结果为函数的四倍。相同的论点,向量y出来。现在,我们有四个slopes-- S1之初,S2中途在中间,S3再次在中间,然后在S4右手。S1的1/6,S2的1/3,S3的1/3,S4的1/6给你的最后一步。这是经典的龙格 - 库塔方法。
卡尔·龙格是一个相当著名的德国数学家和物理学家,谁发表了这个方法,几个人一起,在1895年他制作了一些其他的数学论文,并相当著名。
马丁·库塔独立发现了这个方法,并将其发布在1901年,他是不为别的等近众所周知的。
我想追求一个简单的燃烧模型。因为模型有一些重要的数值性质。如果你点燃一根火柴,火球会迅速增长,直到达到临界大小。然后,由于球体内部燃烧所消耗的氧气量与通过表面可获得的氧气量保持了平衡,所以球体就保持了那个尺寸。
这里的无量纲模型。本场比赛是一个球体,y是它的半径。在Y立方项是球体的体积。和Y立方占了内部的燃烧。
球体的表面是平方成正比收率而y平方项占通过表面的可用氧气。关键参数,重要参数,是初始半径,Y0,Y前功尽弃。
在Y0半径启动和增长,直到在y立方和y平方项相互平衡,此时的生长速度为0,半径不会再增长。我们整合了很长一段时间。我们整合在一个时间是反比于初始半径。这是模型。
这里有一个动画。我们开始用小火焰在这里,小球火焰。你看到的将是一个小半径那里。时间和半径在该图的顶部示出。它开始成长。
当时间到达50时,我们一半。火焰排序爆炸,然后起身半径1,此时两个方面相互平衡。和火焰停止增长。它仍然略有增长这里,虽然你不能看到它在这个这个规模。
让我们为龙格-库塔建立这个。微分方程是y ' = y²- y²。从0开始,临界初始半径是0。01。这意味着我们要积分到2 / y0到200。
我要选择的步长走500步。我只是要挑选有些武断。好了,现在我已经准备好使用ODE4。我会存储在y中的结果。
它上升到1。我要绘制的结果。因此,这里是I的T需要的值。而这里的情节。
现在,你可以看到火焰开始增长。它增长得相当缓慢。在时间间隔的中间,它会爆炸并迅速上升,直到半径为1,然后停留在这里。
现在,这个过渡期是相当狭窄。而且我们将继续研究这个问题。并且是要为数值方法提供了一个挑战,这个过渡区域。
现在这里,我们只是通过它去。我们有一个步长小时,我们挑选漂亮的任意。我们只是产生这些值。我们还真有点想法,这些数字如何准确的。
他们看起来OK。但是,如何准确的他们是谁?这是关于对经典龙格 - 库塔方法的关键问题。如何可靠的是价值观,我们必须在我们这里的图形?
我有四个练习供你考虑。如果微分方程不包含y,那么这个解就是一个积分。龙格-库塔法成为一种经典的数值积分方法。如果你学习过这种方法,那么你应该能够识别这种方法。
数。two--找到Ÿ黄金的精确解等于1加y平方,以0等于零年。然后看到ODE4会发生什么,当你试图解决它在t从0到2的时间间隔。
第三,如果间隔的长度不能被步长整除会发生什么?例如,如果t的最终值是,步长是0。1。不要试图解决这个问题。这只是步长固定的危害之一。
最后,运动four--用的1/1000的初始半径调查火焰的问题。对于T的什么价值呢半径达到90%其最终价值?