SymPy/SciPy:求解具有不同变量的常微分方程组

Question

我是SymPy和Python的新手,我目前正在使用Python 2.7和SymPy 0.7.5,其目标是:a)从文本文件中读取微分方程组b)解决系统问题

我已经读过这个问题和另外一个问题了,它们几乎就是我要找的,但我还有一个额外的问题:我事先并不知道方程组的形式,所以我无法在def里面创建相应的函数脚本,如本例所示.整个过程必须在运行时进行管理.

所以,这里是我的代码的一些片段.假设我有一个包含以下内容的文本文件system.txt:

dx/dt = 0.0387*x - 0.0005*x*y
dy/dt = 0.0036*x*y - 0.1898*y

我所做的是:

# imports
import sympy
import scipy
import re as regex

# define all symbols I am going to use
x = sympy.Symbol('x')
y = sympy.Symbol('y')
t = sympy.Symbol('t')

# read the file
systemOfEquations = []
with open("system.txt", "r") as fp :
   for line in fp :
            pattern = regex.compile(r'.+?\s+=\s+(.+?)$')
            expressionString = regex.search(pattern, line) # first match ends in group(1)   
            systemOfEquations.append( sympy.sympify( expressionString.group(1) ) )

此时,我仍然坚持使用systemOfEquation列表中的两个符号表达式.如果我可以从另一个文件中读取ODE系统的初始条件,为了使用scipy.integrate.odeint,我必须将系统转换为Python可读的函数,如:

def dX_dt(X, t=0):
return array([ 0.0387*X[0] - 0.0005*X[0]*X[1] ,
              -0.1898*X[1] + 0.0036*X[0]*X[1] ])

有没有一种很好的方法在运行时创建它？例如,将函数写入另一个文件,然后将新创建的文件作为函数导入？(也许我在这里很傻,但请记住我对Python比较新:-D)

我已经看到了sympy.utilities.lambdify.lambdify它有可能将符号表达式转换为lambda函数,但我想知道这是否可以帮助我... lambdify当时似乎与一个表达式一起工作,而不是系统.

提前感谢您的任何建议:-)

编辑:

经过极少的修改,沃伦的答案完美无瑕.我有listOfSymbols中所有符号的列表; 此外,它们的显示顺序与odeint将使用的数据X列的顺序相同.所以,我使用的功能是

def dX_dt(X, t):
    vals = dict()
    for index, s in enumerate(listOfSymbols) :
            if s != time :
                    vals[s] = X[index]
    vals[time] = t
    return [eq.evalf(subs=vals) for eq in systemOfEquations]

我只是在我的具体问题中对变量'time'做了一个例外.再次感谢!:-)

Answer 1

如果您要在读取文件的同一脚本中求解系统（因此systemOfEquations可以作为全局变量使用），并且如果中使用的唯一systemOfEquations变量是x，y并且可能是t，您可以dX_dt在同一个文件中定义，如下所示：

def dX_dt(X, t):
    vals = dict(x=X[0], y=X[1], t=t)
    return [eq.evalf(subs=vals) for eq in systemOfEquations]

dX_dt可以用在odeint. systemOfEquations在下面的 ipython 会话中，我已经运行了创建和定义的脚本dX_dt：

In [31]: odeint(dX_dt, [1,2], np.linspace(0, 1, 5))
Out[31]: 
array([[ 1.        ,  2.        ],
       [ 1.00947534,  1.90904183],
       [ 1.01905178,  1.82223595],
       [ 1.02872997,  1.73939226],
       [ 1.03851059,  1.66032942]]