尝试求解 ODE 系统时出错

Getting an Error When Trying to Solve Systems of ODEs

这是我的代码:

import numpy as np
from scipy.integrate import odeint
import matplotlib.pyplot as plt
import math

def VectorField(w, t, p):
    """
    Defines the differential equations to find a in terms of eta.

    Arguments:
        w :  vector of the state variables:
                 w = [eta, a]
        t :  time
        p :  vector of the parameters:
                 p = [H_0, m, r, d]
    """
    eta, a = w
    H_0, m, r, d = p

    # Create f = (eta', a'):
    f = [a, H_0 * (m*a^(-1) + r*a^(-2) + d* a^2)^(1/2)]
    return f

# Parameter values
H_0 = 0.000000000000000002184057032
m = 0.315
r = 0.0000926
d = 0.685

# Initial Conditions
a0 = eta_0 = 0

# ODE solver parameters
abserr = 1.0e-8
relerr = 1.0e-6
stoptime = 10.0
numpoints = 250

#x axis values
t = [stoptime * float(i) / (numpoints - 1) for i in range(numpoints)]

# Pack up the parameters and initial conditions:
initials = [eta_0, a0]
p = [H_0, m, r, d]

# Call the ODE solver.
wsol = odeint(VectorField, initials, t, args=(p,), atol=abserr, rtol=relerr)

但是,当我运行代码时,我可以得到错误信息

TypeError: ufunc 'bitwise_xor' 不支持输入类型,并且根据转换规则''safe''

无法将输入安全地强制转换为任何受支持的类型

代码指向我的向量f。我不确定在声明这个向量 f 时我做错了什么。谁能告诉我我做错了什么?我怎样才能摆脱这个错误?谢谢

当我在不同语言之间来回切换时,这个有时会把我绊倒。

这一行,例如:

f = [a, H_0 * (m*a^(-1) + r*a^(-2) + d* a^2)^(1/2)]

应该是

f = [a, H_0 * (m*a**(-1) + r*a**(-2) + d* a**2)**(1/2)]

因为**是幂运算符,不是^