SciPy odeint 通过矩阵乘法给出不同的结果

计算科学 Python scipy
2021-12-26 01:56:42

我在stackoverflow上问过这个问题,但也许这个社区会对答案有更好的了解。

odeint我目前正在尝试开发一个执行矩阵乘法的函数,同时用in扩展一个微分方程,Python我看到了奇怪的结果。

我有下面的值矩阵和函数,它采用该矩阵的特定值:

from scipy.integrate import odeint
x0_train = [2,0]
dt = 0.01
t = np.arange(0, 1000, dt)
matrix_a = np.array([-0.09999975, 1.999999, -1.999999, -0.09999974])
# Function to run odeint with
def f(x, t, a):
    return [
        a[0] * x[0] + a[1] * x[1],
        a[2] * x[0] - a[3] * x[1]
    ]
odeint(f, x0_train, t, args=(matrix_a,))

>>> array([[ 2.        ,  0.        ],
       [ 1.99760115, -0.03999731],
       [ 1.99440529, -0.07997867],
       ...,
       [ 1.69090227,  1.15608741],
       [ 1.71199436,  1.12319701],
       [ 1.73240339,  1.08985846]])

这似乎是对的,但是当我创建自己的函数来执行乘法/回归时,我看到数组底部的结果完全不同。我有两个稀疏数组,它们提供相同的条件,matrix_a但它们周围有零。

from sklearn.preprocessing import PolynomialFeatures
new_matrix_a = array([[ 0.        , -0.09999975,  1.999999  ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ],
       [ 0.        , -1.999999  , -0.09999974,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ,  0.        ,  0.        ,  0.        ,  0.        ,
         0.        ]])
# New function
def f_new(x, t, parameters):
    polynomials = PolynomialFeatures(degree=5)
    x = np.array(x).reshape(-1,2)
    #x0_train_array_reshape = x0_train_array.reshape(1,2)
    polynomial_transform = polynomials.fit(x)
    polynomial_features = polynomial_transform.fit_transform(x).T
    x_ode = np.matmul(parameters[0],polynomial_features)
    y_ode = np.matmul(parameters[1],polynomial_features)
    return np.concatenate((x_ode, y_ode), axis=None).tolist()

odeint(f_new, x0_train, t, args=(new_matrix_a,))

>>> array([[ 2.00000000e+00,  0.00000000e+00],
       [ 1.99760142e+00, -3.99573216e-02],
       [ 1.99440742e+00, -7.98188169e-02],
       ...,
       [-3.50784051e-21, -9.99729456e-22],
       [-3.50782881e-21, -9.99726119e-22],
       [-3.50781711e-21, -9.99722781e-22]])

如您所见,我在数组末尾得到了完全不同的值。我一直在运行我的代码,但似乎找不到它们不同的原因。有没有人有明确的原因或者我做错了什么f_new理想情况下,我想开发一个可以在其中获取任何值的函数matrix_a,这就是我尝试创建这个新函数的原因。

提前致谢。

1个回答

我最初专注于第二个功能,假设第一个是正确的,但我意识到问题实际上出在第一个功能上。

您的矩阵乘法在第一个函数中不正确。它应该如下所示:

def f(x, t, a):
    return [
        a[0] * x[0] + a[1] * x[1],
        a[2] * x[0] + a[3] * x[1]
    ]

您是在减去a[3] * x[1]而不是添加它。这正是-0.2*x0_train[1]我在评论中指出的;自从a30.1,并且您减去它而不是添加它,这将导致输出的第二部分2*a[3]*x[1]在每一步都关闭。进行此更改时,两者在传递到 时都会产生相同的结果odeint