我目前正在使用RK4算法对相当深奥的系统进行物理模拟。

它涉及奇点(例如r-> 0),并在一个方程中除以该半径。会有奇点,发生这种情况后,我需要程序停止进一步的计算。一旦程序过于接近奇点,我便建立了一个> break操作,但是有时用户可能没有选择足够好的阈值以在出现无限大之前突破。无法预先确定此阈值的值-系统在动力学方面明显混乱。我决定尝试捕获异常(特别是RuntimeWarning:在double_scalars中遇到的溢出),以便用户知道他们选择的奇点阈值太低。

for i in range(indices - 1):
    try:
        if((theta0 == 0 or theta0 == pi) and i == 0):
            print('Unstable or stable equilibrium chosen as initial value. Please change the initial angle.')
            flag = True
            break
        if(variables[i][0] <= (singularitythreshold * r0)):
            print('The trajectory came within the threshold for identifying a singularity (%.5f%%). The program has finished early to avoid infinities.' % (singularitythreshold * r0 * 100))
            break
        k1 = step * RKaccel(variables[i], times[i])
        k2 = step * RKaccel(variables[i] + k1 / 2, times[i] + step/2)
        k3 = step * RKaccel(variables[i] + k2 / 2, times[i] + step/2)
        k4 = step * RKaccel(variables[i] + k3, times[i] + step)

        variables[i + 1] = variables[i] + k1/6 + k2/3 + k3/3 + k4/6
    except RuntimeWarning:
        print('A Runtime Warning was triggered, indicating infinities as r -> 0. Increase the singularity threshold.')
        flag = True
        print('Plotting procedures have been abandoned to avoid nonsensical data.')
        break


我已经阅读了10多篇有关如何处理此类问题的文章(包括seterr和seterrfunc),但似乎无法正确解决。它永远不会触发异常。

计算是在RKaccel中进行的,我认为此函数与捕获错误无关,因此我没有将其包括在内(很多令人讨厌的方程式)。但是,我包含了要打印的特定警告:

/tmpw0ojdH.py:35: RuntimeWarning: overflow encountered in double_scalars
  radiusdotdot = ((radius / (1 + mu)) * ((thetadot) ** 2)) + (((g * cos(theta)) - (g * mu)) / (1 + mu))
/tmpw0ojdH.py:36: RuntimeWarning: overflow encountered in double_scalars
  thetadotdot  = - ((g * sin(theta)) / radius) - (2 * ((radiusdot) * (thetadot)) / radius)

最佳答案

我已经解决了我的问题。关键是在玩seterr()。通过添加

seterr(all = 'raise')


和改变

except RuntimeWarning:




except FloatingPointError:


我的例外被触发,所有问题都已解决。希望这对以后的人有所帮助。

10-05 21:14