1

我正在使用 sympy.dsolve 来解决衰减链的简单 ODE。对于不同的衰减率（例如 lambda_1 > lambda_2），我得到的答案是错误的。代入 C1=0 后，我得到一个简单的指数

-N_0*lambda_1*exp(-lambda_1*t)/(lambda_1 - lambda_2)

而不是正确的答案：

(exp(-lambda_1*t)-exp(-lambda_2*t)).

我究竟做错了什么？这是我的代码

sp.var('lambda_1,lambda_2 t')
sp.var('N_0')
daughter = sp.Function('N2',Positive=True)(t)

stage1 = N_0*sp.exp(-lambda_1*t)
eq = sp.Eq(daughter.diff(t),stage1*lambda_1 - daughter*lambda_2)

sp.dsolve(eq,daughter)
4

2 回答 2

0

您的微分方程是（使用较短的变量标识符）

y' = A*N0*exp(-A*t) - B*y

应用积分因子exp(B*t)得到等效值

(exp(B*t)*y(t))' = A*N0*exp((B-A)*t)

整合得到

exp(B*t)*y(t) = A*N0*exp((B-A)*t)/(B-A) + C

y(t) = A*N0*exp(-A*t)/(B-A) + C*exp(-B*t)

这正是求解器计算的结果。

于 2018-01-16T22:35:22.673 回答
0

您是否插入了 C1 = 0，期望得到 y(0) = 0 的解？这不是它的工作原理。C1 是公式中的任意常数，将其设置为 0 并不能保证当 t = 0 时表达式的计算结果为 0。

这是一个正确的方法，一步一步

sol1 = sp.dsolve(eq, daughter)

这将返回 Piecewise，因为不知道两个 lambda 是否相等：

Eq(N2(t), (C1 + N_0*lambda_1*Piecewise((t, Eq(lambda_2, lambda_1)), (-exp(lambda_2*t)/(lambda_1*exp(lambda_1*t) - lambda_2*exp(lambda_1*t)), True)))*exp(-lambda_2*t))

我们可以澄清它们不是：

sol2 = sol1.subs(Eq(lambda_2, lambda_1), False)

得到

Eq(N2(t), (C1 - N_0*lambda_1*exp(lambda_2*t)/(lambda_1*exp(lambda_1*t) - lambda_2*exp(lambda_1*t)))*exp(-lambda_2*t))

接下来，我们需要 C1 使得当 t = 0 时右侧变为 0。因此，取右侧，将 0 代入 t，求解 C1：

C1 = Symbol('C1')
C1val = solve(sol2.rhs.subs(t, 0), C1, dict=True)[0][C1]

（没有必要包含dict=True，但我喜欢它，因为它强制求解的统一输出：它是一个字典数组。）

顺便说一句，C1val 现在是N_0*lambda_1/(lambda_1 - lambda_2). 放进去：

sol3 = sol2.subs(C1, C1val).simplify()

你有它：

Eq(N2(t), N_0*lambda_1*(exp(lambda_1*t) - exp(lambda_2*t))*exp(-t*(lambda_1 + lambda_2))/(lambda_1 - lambda_2))

这个表达式等价于N_0*lambda_1*(exp(-lambda_2*t) - exp(-lambda_1*t))/(lambda_1 - lambda_2)虽然 SymPy 似乎不愿意在这里组合指数。

于 2018-01-16T23:09:40.483 回答