Python 复杂的耦合 ODE 错误

Python complex coupled ODEs error

目前,我正在尝试求解具有复数项的耦合 ODE 系统。我正在使用 scipy.integrate.ODE,我已经成功地解决了之前涉及仅实数项的耦合 ODE 系统的问题。在那种情况下我使用了odeint,它不适合我现在面临的问题。该系统有 9 个耦合 ODE,这是我的代码:

from pylab import *
import numpy as np
from scipy.integrate import ode
#Here I am calculating a parameter H of the system
h=1
w1=0.1
w2=0.12
g=1
H=array([[h/2*(w1+w2),h/2*w1+h*g,h*w2/2,0], [h/2*w1+h*g,h/2*(w1-w2),0,-h*w2/2], [h*w2/2,0,-h/2*(w1-w2),-h/2*w1-h*g],[0,-h*w2/2,-h/2*w1-h*g,-h/2*(w1+w2)]])

#this is the system of ODEs
def resol(p1111, p1212, p2121, p1112, p1121, p1122, p1221, p1222, p2122,t,par):
     gamma, h, H = par
     i=1j
     dp1111 = (H[0,1]*conjugate(p1112)+H[0,2]*conjugate(p1121)+H[0,3]*conjugate(p1122)-H[1,0]*p1112-H[2,0]*p1121-H[3,0]*p1122)*i+delta*(p1212-p1111)
     dp1112 = (H[0,0]*p1112+H[0,1]*p1212+H[0,2]*conjugate(p1221)-H[0,1]*p1111-H[1,1]*p1112-H[3,1]*p1122)*i+delta*(conjugate(p1112)-p1112)
     dp1121 = (H[0,0]*p1121+H[0,1]*p1221+H[0,2]*p2121-H[0,2]*p1111-H[2,2]*p1121-H[3,2]*p1122)*i+delta*(p1222-p1121)
     dp1122 = (H[0,0]*p1122+H[0,1]*p1222+H[0,2]*p2122-H[1,3]*p1112-H[2,3]*p1121-H[3,3]*p1122)*i+delta*(p1221-p1122)
     dp1212 = (H[1,0]*p1112+H[1,3]*conjugate(p1222)-H[0,1]*conjugate(p1112)-H[3,1]*p1222)*i+delta*(p1111-p1212)
     dp1221 = (H[1,0]*p1121+H[1,1]*p1221+H[1,3]*conjugate(p2122)-H[0,2]*conjugate(p1112)-H[2,2]*p1221-H[3,2]*p1222)*i+delta*(p1122-p1221)
     dp1222 = (H[1,0]*p1122+H[1,1]*p1222+H[1,3]*p2122-H[1,3]*p1212-H[2,3]*p1221-H[3,3]*p1222)*i+delta*(p1121-p1222)
     dp2121 = (H[2,0]*p1121+H[2,3]*conjugate(p2122)-H[0,2]*conjugate(p1121)-H[3,2]*p2122)*i+delta*(1-p1111-p1212-2*p2121)
     dp2122 = (H[2,0]*p1122+H[2,2]*p2122+H[2,3]*(1-p1111-p1212-p2121)-H[1,3]*conjugate(p1221)-H[2,3]*p2121-H[3,3]*p2122)*i+delta*(conjugate(p2122)-p2122)
     sol = [dp1111, dp1212, dp2121, dp1112, dp1121, dp1122, dp1221, dp1222, dp2122]

     return sol
#I set the Initial values of echa term, note that there are 9 because I have 9 terms
condin=np.array([1/3,0,2/3,0,sqrt(2)/3,0,0,0,0])
#I set the parameters
par=[g, h, H]
#I set the integrator, I don't use jacobian because I it is tedious to calculate it for my system
r=ode(resol).set_integrator('zvode', method='bdf',with_jacobian=False)
#Just defining the time and the steps
t_start = 0.0
t_final = 10.0
delta_t = 0.1
num_steps = np.floor((t_final - t_start)/delta_t) + 1

r.set_initial_value(condin, t_start) #I give initial values to the integrator
r.set_f_params(par) #the parameter aswell

t = np.zeros((num_steps, 1)) #I define a (x,1) dim matrix  to save each term of the time
resultados = np.zeros((num_steps, 9)) # the same to save the results but with a (x,9) dim matrix
t[0] = t_start
resultados[0] = condin
k = 1
while r.successful() and k < num_steps:
    r.integrate(r.t + delta_t)
    t[k] = r.t #save the time and the results
    resultados[k] = r.y
    k += 1

如您所见,我将我的函数定义为 resol(p1111, p1212, p2121, p1112, p1121, p1122, p1221, p1222, p2122,t,par):。使用 odeint 时我用不同的方式写了它:

def resol(var,t,par):
    p1111, p1212, p2121, p1112, p1121, p1122, p1221, p1222, p2122 = var
    gamma, h, H = par

但是当使用 integrate.ode 时,它给了我这个错误:

p1111, p1212, p2121, p1112, p1121, p1122, p1221, p1222, p2122 = var
TypeError: 'float' object is not iterable

反正我按照我发布的方式定义函数解决了它。问题是,当我 运行 代码时,我总是得到这个错误,我不知道为什么:

create_cb_arglist: Failed to build argument list (siz) with enough arguments (tot-opt) required by user-supplied function (siz,tot,opt=3,11,0).
Traceback (most recent call last):
File "Untitled.py", line 48, in <module>
r.integrate(r.t + delta_t)
File "/Users/kilian/src/scipy/scipy/integrate/_ode.py", line 333, in integrate
self.f_params, self.jac_params)
File "/Users/kilian/src/scipy/scipy/integrate/_ode.py", line 760, in run
args[5:]))
vode.error: failed in processing argument list for call-back f.

我不知道这个错误是不是因为我在我的函数中使用函数 conjugate() 来计算热力学的复共轭。但是我觉得不是这个原因。

如果有人能帮我解决这个问题,我将不胜感激。提前谢谢你

这应该可以解决这个问题:

def resol(t, y, par):
    p1111, p1212, p2121, p1112, p1121, p1122, p1221, p1222, p2122 = y
    gamma, h, H = par

你会遇到更多问题,例如 delta 需要被调用 delta_t(需要重命名)和 ComplexWarning: Casting complex values to real discards the imaginary part 但它运行。