sympy:求解微分方程组

sympy: dsolve for a system of differential equations

如何解一首诗的系统?

这就是我试过的方法:

from sympy import *;
init_printing()

Ex,Ey,Ez,Bx,By,Bz,c,m0,t,u,we,wb,q=symbols('E_x E_y E_z B_x B_y B_z c m_0 t u omega_E omega_B q')
tau=dynamicsymbols('tau')

u0=Function('u^0')(tau) 
u1=Function('u^1')(tau)
u2=Function('u^2')(tau)
u3=Function('u^3')(tau) # edit

u=Matrix([u0,u1,u2,u3])
Emat=Matrix([[0,Ex/c,Ey/c,Ez/c],[Ex/c,0,Bz,-By],[Ey/c,-Bz,0,Bx],[Ez/c,By,-Bx,0]])
RHS=q/m0*Emat*u

odesys=[Eq(u0.diff(),RHS[0]),Eq(u1.diff(),RHS[1]),Eq(u2.diff(),RHS[2]),Eq(u3.diff(),RHS[3])]
odesys
dsolve(odesys,[u0(tau),u1(tau),u2(tau),u3(tau)])

不幸的是,它给了我一个 TypeError:

TypeError                                 Traceback (most recent call last)
<ipython-input-94-f4af1f7e485c> in <module>()
     16 odesys=[Eq(u0.diff(),RHS[0]),Eq(u1.diff(),RHS[1]),Eq(u2.diff(),RHS[2]),Eq(u3.diff(),RHS[3])]
     17 odesys
---> 18 dsolve(odesys,[u0(tau),u1(tau),u2(tau),u3(tau)])

TypeError: 'u^0' object is not callable

我已经尝试了很多东西,但我就是想不出正确的语法来让它起作用。 有什么想法吗?

顺便说一句:如果我使用动态符号而不是函数,它没有帮助:

...
u0,u1,u2,u3=dynamicsymbols('u^0 u^1 u^2 u^3')
...
dsolve(odesys,[u0(t),u1(t),u2(t),u3(t)])

-> 相同类型错误

编辑: 我让 dsolve 变得更容易一些,只要求它求解第一个方程(共 4 个方程),同时保持我的函数定义:

dsolve(Eq(u0.diff(),RHS[0]))  

-> 这确实有效。

下一步:求解第一个和第二个方程:

dsolve([Eq(u0.diff(),RHS[0]),Eq(u1.diff(),RHS[1])])  

->“...未实施错误”

这是物理学家每天都会遇到的问题。最好坚持使用 Mathematica

你写u0=Function('u^0')(tau),意思是u0是函数,符号为"u^0",参数为"tau"。那么u0(t)就没有意义了;您已经说过函数的参数是什么。就像写u^0(tau)(t)一样,没有意义。

因此,要么让 u0=Function('u^0') 不带参数(然后在等式中提供),要么将 u0(t) 更改为 u0

您在 u3 的定义中也有错别字,这使得它与 u0 相同。

最后,关于 tau 是 t 的函数以及微分方程中未知函数的参数的整个事情是模糊的。你是根据一个功能来区分的,这对人们来说已经足够难以理解了,更不用说计算机了。

我建议陈述系统,以便清楚自变量是什么,并且在该变量中采用导数。

在 SymPy 中求解简单 ODE 系统的示例

from sympy import *
t = symbols('t')
x = Function('x')(t)
y = Function('y')(t)
sol = dsolve([x.diff() - y, y.diff() + x])

让你 [Eq(x(t), C1*sin(t) + C2*cos(t)), Eq(y(t), C1*cos(t) - C2*sin(t))]

我意识到这已经非常过时了,但我正致力于将 Maple 中的东西移至 Python 并偶然发现了这个 post。通过一些重新安排,我能够让它为最多前三个方程式提出一个“解决方案”;一旦我放入第四个,Python 似乎挂起:

#%%
from sympy import *
init_printing()

#%%
Ex,Ey,Ez,Bx,By,Bz,c,m0,t,u,we,wb,q=symbols(
    'E_x E_y E_z B_x B_y B_z c m_0 t u omega_E omega_B q', 
    real=True, positive=True)

#%%
u0=Function('u^0')
u1=Function('u^1')
u2=Function('u^2')
u3=Function('u^3')

#%%
u=Matrix([u0(t),u1(t),u2(t),u3(t)])
Emat=Matrix([
    [0,Ex/c,Ey/c,Ez/c],
    [Ex/c,0,Bz,-By],
    [Ey/c,-Bz,0,Bx],
    [Ez/c,By,-Bx,0]])
RHS=q/m0*Emat@u

#%%
odesys=[
    Eq(u0(t).diff(t),RHS[0]),
    Eq(u1(t).diff(t),RHS[1]),
    Eq(u2(t).diff(t),RHS[2]),
    Eq(u3(t).diff(t),RHS[3])]

#%%
for k in range(1, 5):    
    print("{} equations".format(k))
    soln = dsolve(odesys[0:k])
    pretty_print(soln)

我没有对关于真实物质或正物质的断言进行任何严格的检查...此外,前三个的“解决方案”包含遗漏任何变量的积分,因此称它们为“解决方案”是一个伸展但至少 Python 做了 某事 :) 哦,这是多么了不起的事情 - 这是 u^0 三个方程式: