2

我目前正在探索 Lorenz 系统，Python并R注意到ode包中的细微差别。odeintfromPython和ode都说他们lsoda用来计算他们的导数。但是，lsoda对两者使用该命令似乎会给出截然不同的结果。我已经尝试ode45让ode函数 inR得到更相似的东西，Python但我想知道为什么我不能得到完全相同的结果：

from scipy.integrate import odeint
def lorenz(x, t):
    return [
        10 * (x[1] - x[0]),
        x[0] * (28 - x[2]) - x[1],
        x[0] * x[1] - 8 / 3 * x[2],
    ]


dt = 0.001
t_train = np.arange(0, 0.1, dt)
x0_train = [-8, 7, 27]
x_train = odeint(lorenz, x0_train, t_train)


x_train[0:5, :]
array([[-8.        ,  7.        , 27.        ],
       [-7.85082366,  6.98457874, 26.87275343],
       [-7.70328919,  6.96834721, 26.74700467],
       [-7.55738803,  6.95135316, 26.62273959],
       [-7.41311133,  6.93364263, 26.49994363]])

library(deSolve)
n <- round(100, 0)
# Lorenz Parameters: sigma, rho, beta
parameters <- c(s = 10, r = 28, b = 8 / 3)
state <- c(X = -8, Y = 7, Z = 27) # Initial State
# Lorenz Function used to generate Lorenz Derivatives
lorenz <- function(t, state, parameters) {
  with(as.list(c(state, parameters)), {
    dx <- parameters[1] * (state[2] - state[1])
    dy <- state[1] * (parameters[2] - state[3]) - state[2]
    dz <- state[1] * state[2] - parameters[3] * state[3]
    list(c(dx, dy, dz))
  })
}
times <- seq(0, ((n) - 1) * 0.001, by = 0.001)
# ODE45 used to determine Lorenz Matrix
out <- ode(y = state, times = times,
           func = lorenz, parms = parameters, method = "ode45")[, -1]
out[1:nrow(out), , drop = FALSE]
             X        Y        Z
 [1,] -8.00000000 7.000000 27.00000
 [2,] -7.85082366 6.984579 26.87275
 [3,] -7.70328918 6.968347 26.74700
 [4,] -7.55738803 6.951353 26.62274
 [5,] -7.41311133 6.933643 26.49994

我不得不打电话out[1:nrow(out), , drop = FALSE]来获得完全提供的小数位，似乎head四舍五入到最接近的五分之一。我知道这是非常微妙的，但希望得到完全相同的结果。有谁知道这是否不仅仅是和之间的舍入R问题Python？

提前致谢。

4

1 回答 1

2

所有求解 ODE 的数值方法都是达到给定精度的近似值。atol=1e-6, rtol=1e-6默认情况下，deSolve 求解器的精度设置atol为绝对容差和rtol相对容差。此外，ode45还有一些额外的参数来微调自动步长算法，它可以利用插值。

要增加容差，例如设置：

out <- ode(y = state, times = times, func = lorenz, 
  parms = parameters, method = "ode45", atol = 1e-10, rtol = 1e-10)

最后，我建议使用 odepack 求解器，例如lsodaorvode而不是经典的ode45. 更多详细信息可以在ode和lsoda帮助页面以及ode45?rkMethod 帮助页面中找到。

也可能存在类似的参数odeint。

最后一点：由于 Lorenz 是一个混沌系统，局部误差会由于误差放大而导致发散行为。这是混沌系统的一个基本特征，从理论上讲，从长远来看，混沌系统是不可预测的。因此，无论您做什么，以及您设置了多少精度，模拟轨迹都不是“真实的”，它们只是显示出类似的模式。

于 2021-05-20T16:40:58.173 回答