通过 Sympy 求解的常微分方程不满足给定的初始条件



我试图在Sympy中使用dsolve函数解决以下简单的ODE

v = sp.Function('v', real=True)
t, g, m, c = sp.symbols('t g m c', positive=True, real=True)
eq = sp.Eq(v(t).diff(t), g - c/m*v(t)**2)
sol = sp.dsolve(eq, v(t), ics={v(0):0})

结果并不像预期的那么简单,但奇怪的是,如果您检查解决方案以确保给定的初始条件得到满足,它不会返回正确的结果(v(0)=0),正如ics参数在dsolve中传递的那样。

sol.subs(t,0)

替换u=c/m*vc/m*g=a^2得到方程u'=a^2-u^2,这是可分离的,导致

log|u+a|-log|u-a| = 2*a*t+c

应该转换为

(u+a)/(u-a) = C*exp(2*a*t)
u = a*(C*exp(2*a*t)+1)/(C*exp(2*a*t)-1)

,则初始条件导致C=-1u(t)=a*tanh(a*t)

由于sympy解中的常数显然没有从指数中移除,并且绝对值向负号的分解也显然没有遵循,因此sympy最终得到c=log(-1)/2的变化,这在此实数方程的上下文中是没有意义的。此外,-1在复对数的通常分支上,使其更加不确定。在实际的解决方案中,我们发现log((-1)**sqrt(m))不是一个定义良好的表达式。

问题是dsolve找到一般解:

In [41]: sol
Out[41]: 
√g⋅√m               
v(t) = ───────────────────────────────────
⎛             ⎛ C₁⋅√c⋅√g⋅m⎞⎞
⎜          log⎝ℯ          ⎠⎟
⎜√c⋅√g⋅t - ────────────────⎟
⎜                 2        ⎟
√c⋅tanh⎜──────────────────────────⎟
⎝            √m            ⎠

在这里尝试解决C1之前,最好应用一些常量简化:

In [42]: [logterm] = sol.atoms(log)
In [43]: C1 = Symbol('C1')
In [44]: sol2 = sol.subs(logterm, C1)
In [45]: sol2
Out[45]: 
√g⋅√m         
v(t) = ───────────────────────
⎛  C₁          ⎞
⎜- ── + √c⋅√g⋅t⎟
⎜  2           ⎟
√c⋅tanh⎜──────────────⎟
⎝      √m      ⎠
In [46]: [c] = solve(sol2.rhs.subs(t, 0), C1)
In [47]: c
Out[47]: ⅈ⋅π⋅√m
In [48]: sol2.subs(C1, c)
Out[48]: 
√g⋅√m          
v(t) = ─────────────────────────
⎛          ⅈ⋅π⋅√m⎞
⎜√c⋅√g⋅t - ──────⎟
⎜            2   ⎟
√c⋅tanh⎜────────────────⎟
⎝       √m       ⎠
In [49]: sol2.subs(C1, c).simplify()
Out[49]: 
⎛√c⋅√g⋅t⎞
√g⋅√m⋅tanh⎜───────⎟
⎝   √m  ⎠
v(t) = ───────────────────
√c 

也许实际上最好在求解过程的早期阶段求解常数。也许在笔和纸上,当积分发生时,当方程仍然是隐式形式时,解出常数是有意义的。

最新更新