尝试使用solve_ivp的矢量化选项，奇怪的是它会抛出一个错误，即y0必须是一维的。微量元素：
from scipy.integrate import solve_ivp
import numpy as np
import math
def f(t, y):
theta = math.pi/4
ham = np.array([[1,0],[1,np.exp(-1j*theta*t)]])
return-1j * np.dot(ham,y)
def main():
y0 = np.eye(2,dtype= np.complex128)
t0 = 0
tmax = 10**(-6)
sol=solve_ivp( lambda t,y :f(t,y),(t0,tmax),y0,method='RK45',vectorized=True)
print(sol.y)
if __name__ == '__main__':
main()
调用签名是 fun(t, y)。这里 t 是一个标量，ndarray y 有两个选项：它可以具有形状 (n,)；那么 fun 必须返回形状为 (n,) 的 array_like。或者它可以具有形状（n，k）；那么 fun 必须返回一个形状为 (n, k) 的 array_like，即每一列对应于 y 中的单个列。两个选项之间的选择由向量化参数决定（见下文）。矢量化实现允许通过有限差分更快地逼近雅可比行列式（刚性求解器所需）。
错误 ：
ValueError：
y0必须是一维的。
Python 3.6.8
scipy。版本 “1.2.1”