从系列中:在MATLAB中求解ode
克利夫·莫勒
包含高阶导数的常微分方程被改写为只包含一阶导数的向量系统。给出了典型的范德堡尔非线性振子的例子。随着参数的增加,VdP方程变得僵硬。
许多数学模型涉及高阶导数。但是MATLAB ODE解算器仅适用于一阶常微分方程组。所以我们必须重写模型,只考虑一阶导数。让我们看看如何用一个非常简单的模型,谐振子来实现这一点。
x双素数加x等于0。这涉及到一个二阶导数。为了把它写成一阶系统,我们引入了向量y。这是一个包含两个分量x和x的导数的向量。
我们只是改变了符号让y有两个分量,x和x '那么y的导数就是向量x和x '微分方程变成y2 ' + y1 = 0。看到我们是如何重写这个微分方程的了吗。所以y2 '是x ' '的函数?
一旦你这样做了,其他的事情就很容易了。向量系统现在是y1, y2 '是y2 - y1。第一个分量是y1 ' = y2。也就是说第一分量的导数等于第二分量。这是微分方程本身。Y2 '是- y1是谐振子微分方程。
当我们为MATLAB编写一个自治函数时,这里是自治函数。f是t和y的自治函数,它不依赖于t。首先是向量,现在是列向量。f的第一个分量是y2。第二个分量是-y1。
这里的第一个组件只是一个符号问题。所有的内容都在第二部分,它表示微分方程。
对于某些初始条件,假设初始条件是x (0) = 0, x '(0) = 1。用向量y表示,就是y1 (0) y的第一个分量是0。第二个分量是1。
或者在向量项中,初始向量是0,1。这意味着它们的解是正弦t和余弦t。当我们在MATLAB中编写初始条件时,它是列向量0,1。
让我们打开MATLAB命令窗口。这是微分方程。y1素数是y2。y2素数是-y1。这是谐振子。
我们将从0积分到2pi,因为它们是trig函数。我将要求以2π乘以36的步长输出,这相当于每10度,就像机场的跑道一样。
我需要一个初始条件。y0不是0。对于这两个分量,我需要一个列向量0,1。我将使用ODE45,如果我在没有输出参数的情况下调用它,微分方程f的ODE45,t跨越时间间隔,y0是初始条件。
如果我不带输出参数调用ODE45,它就会自动绘制解决方案。这里我们得到了cos t从1开始,sin t从0开始的图像。
现在如果我回到命令窗口,要求捕获t和y的输出,然后我得到输出的向量。37步,向量t,和两个分量y,包含正弦和余弦的两列。现在我可以把它们画在相平面上。
阴谋对付y2。如果我调整轴,我会得到一个很好的图,每10度有一个小圆圈,就像我说的,就像机场的跑道一样。
范德堡尔振荡器由荷兰电气工程师于1927年引入,用来模拟真空管电路中的振荡。它有一个非线性阻尼项。从那时起,它就被用来模拟各种领域的现象,包括地质学和神经学。
它表现出混乱的行为。我们对它感兴趣是为了进行数值分析,因为随着参数mu的增加,问题变得越来越棘手。为了将它写成一阶系统用于MATLAB ODE求解器,我们引入向量y,包含x和x '。
y '等于x '和x ' '然后微分方程写成y '的第一个分量是y2。然后微分方程写在y的第二分量中。
这是非线性阻尼项减去y1。当= 0时,这是谐振子。这里是匿名函数。
让我们用范德波尔振荡器做一些实验。首先,我必须定义μ的值。我将选择一个适度的mu值。μ等于100。现在使用mu set,我可以定义匿名函数。
它包含一个mu的值。这是范德波尔方程,在f的第二分量。我将收集ODE解算器工作有多努力的统计数据。为此,我将使用ODE集,并告诉它我想打开统计。
我需要一个初始条件。现在我要运行ODE45来解决这个问题。我只指定了一个起始值t和一个最终值t。ODE45将选择自己的时间步长。我知道t final等于460,它将在两个周期的解上积分。
现在我们可以看一会儿了。它采取了很多步骤。而且它开始减速,因为它需要越来越多的步骤。现在它开始变得非常缓慢,因为它变得僵硬。
我把摄像机关了一会儿,这样你就不用看这些步骤了。我们正努力提高到460度。快结束的时候我会把它打开。
好吧,摄像机已经关了大约三分钟了。你可以看到我们已经走了多远。我们离终点还很远。所以我要按下这里的停止按钮。我们将返回到命令窗口。
哦,我们还没讲到最后。让我回到时间间隔,试一下这个值。看看这是怎么回事。所以这仍然需要很多步骤。
但我们可以-,这将持续一个时期。我们将在这里结束。
我会让照相机开着直到我们结束。好吧,这花了不到一分钟的时间。它走了近15000步。让我们用一个僵硬的解算器来运行它。
那里所以只花了半秒钟,只走了500步。所以这里有一个适度的刚度例子。让我们用刚性解算器来研究范德波尔方程。让我们捕获输出并绘制它。
因为那个情节不是很有趣。我想用两种不同的方式来描绘它。再一次,我想回到-,捕捉几个周期。
让我们绘制一个当前组件作为时间的函数。就在这里。这是范德波尔方程。
你可以看到它有一个初始的瞬态,然后它进入这个周期振荡这里有一些非常陡峭的峰值。甚至这个僵硬的解算者也在这些快速变化中努力工作。这里有相当多的点,因为它选择了步长。
现在,让我们回到命令行,做一个相平面图。这是带有阻尼的振荡器的相平面图。它不在圆附近,如果是0就会在圆附近。
这就是范德堡尔振荡器的特性。
你也可以从以下列表中选择一个网站:
选择中国网站(中文或英文)以获得最佳网站性能。其他MathWorks国家站点没有针对您所在位置的访问进行优化。