<p>如果你要提前积分,这是个好主意,而且绝对是可行的方法,那么你可以把方程写下来,作为t的函数:</p>
<pre><code>x'' = -kx/m
x'' + (k/m)x = 0
r^2 + k/m = 0
r^2 = -(k/m)
r = i*sqrt(k/m)
x(t) = A*e^(i*sqrt(k/m)t)
= A*cos(sqrt(k/m)t + B) + i*A*sin(sqrt(k/m)t + B)
= A*cos(sqrt(k/m)t + B)
</code></pre>
<p>从初始条件我们知道</p>
<pre><code>x(0) = 12 = A*cos(B)
v(0) = 0 = -sqrt(k/m)*A*sin(B)
</code></pre>
<p>只有当我们选择A=0或B=0或B=Pi时,这些方程中的第二个才是真的</p>
<ul>
<li>如果A=0,则第一个方程没有解</李>
<li>如果B=0,则第一个方程的解A=12</李>
<li>如果B=Pi,第一个方程的解A=-12</李>
</ul>
<p>我们可能更喜欢B=0和A=12。这给</p>
<pre><code>x(t) = 12*cos(sqrt(k/m)t)
v(t) = -12*sqrt(k/m)*sin(sqrt(k/m)t)
a(t) = -12*(k/m)cos(sqrt(k/m)t)
</code></pre>
<p>因此,在任何增量时间t[n+1]=t[n]+dt,我们可以简单地计算t[n]的精确位置、速度和加速度,而不会累积任何漂移或误差</p>
<p>综上所述,如果你对给定一个任意常微分方程的x(t)、v(t)和a(t)的数值计算感兴趣,那么答案就难多了。有很多很好的方法可以被称为数值积分。欧拉方法是最简单的:</p>
<pre><code>// initial conditions
t[0] = 0
x[0] = …
x'[0] = …
…
x^(n-1)[0] = …
x^(n)[0] = 0
// iterative step
x^(n)[k+1] = f(x^(n-1)[k], …, x'[k], x[k], t[k])
x^(n-1)[k+1] = x^(n-1)[k] + dt * x^(n)[k]
…
x'[k+1] = x'[k] + dt * x''[k]
x[k+1] = x[k] + dt * x'[k]
t[k+1] = t[k] + dt
</code></pre>
<p>选择的dt值越小,运行固定时间所需的时间就越长,但得到的结果越精确。这基本上是对函数及其所有导数进行黎曼和,直到ODE中涉及的最高导数</p>
<p>更准确的版本是辛普森规则,它做的是相同的事情,但取最后一个时间量的平均值(而不是任一端点的值;上面的示例使用间隔的开始)。保证区间内的平均值比任一端点更接近区间内的真值(除非函数在该区间内为常数,在这种情况下,辛普森至少也一样好)</p>
<p>对于常微分方程来说,最好的标准数值积分方法(假设不需要像蛙跳法这样的方法来提高稳定性)可能是龙格-库塔法。一个足够阶数的自适应时间步长龙格库塔方法通常可以做到这一点,并给出准确的答案。不幸的是,解释Runge-Kutta方法的数学可能过于先进和耗时,无法在这里介绍,但您可以在网上或数字配方(Numeric Recipes)中找到有关这些和其他先进技术的信息,这是一系列关于数字方法的书籍,其中包含许多非常有用的代码示例</p>
<p>即使是龙格-库塔方法也基本上是通过在时间量子内对函数值的精确猜测来工作的。他们只是用更复杂的方法来做,这可以证明减少了每一步的误差</p>