所以基本上我使用 odeint 来解决我拥有的 ode 并且我希望将它用作另一个方程中的变量(最终实际上是我的解的方程)。
所以在这里我解决了第一个 ode(我的代码实际上是在解决一个 ODE 矩阵,但同样的问题仍然适用)
Omeg = odeint(f, zeros, sMesh, args=(dOmeg,om,))
我的 f() 函数定义如下:
def f(y, t, dh, v):
"""Returns dh/ds with initial values plugged in"""
# Create vector for solved DEs
solvedEqs = range(len(dh))
# Establish a dictionary for values to substitue
subValues = {}
for i in range(0,len(v)):
subValues.update({v[i]:y[i]}) #add {symbol:value} to dictionary
# Set each element of dh with initial values subed in to return array
for i in range(0,len(solvedEqs)):
solvedEqs[i] = dh[i].subs(subValues)
return solvedEqs
所以在所有这些好东西之后,我得到了一个数组或数字,点,我可以绘制以显示解决方案等。但我的问题是,我需要使用这个解决方案(为了简单起见,假设它是一个而不是解决方案矩阵)插入我实际上试图找到的功能:
Hs = H0
for i in range(1,order): # length of order
Hs += ((bernoulli(i-1)*1.)/math.factorial(i-1))*(adjointOp(H0,Omeg,i))
其中 bernoulli() 只是一个返回该订单的伯努利数的函数,而 adjointOp() 只是给定订单的伴随运算符(交换子关系)。
我首先在 Mathematica 中做到了这一点,其中 Omeg 的解决方案是 InterpolationFunction,但这只是作为变量工作并相乘等。但是,我不知道如何在 python 中处理这种情况。