1

我对 Mathematica NDSolve 中的约束有疑问

如果我有两个线性微分方程

NDSolve[x'[t]= A x[t]+B y[t], 
        y'[t]= C x[t]+D y[t], 
        x[0]==0.002, y[0]==0.005, {x,y}, {t,0.10000}]

我怎样才能施加约束

x[t]+y[t] == 1 

对于任何时间 t

谢谢

4

2 回答 2

2

答案是一般情况下你不能。这不是 Mathematica 的限制,而是数学的限制。

您正在尝试同时求解三个方程:

x' = Ax+By
y' = Cx+Dy
1 = x + y

我将使用小写字母a,b,c,d来避免 Mathematica 出现问题。

正如贝利撒留建议的那样,您不需要手动进行替换。您可以简单地选择其中任何两个并要求 Mathematica 解决:

a = 1; b = 2; c = 3; d = 4; x0 = 0.2; y0 = 0.8;

NDSolve[{x'[t] == a x[t] + b y[t], y'[t] == c x[t] + d y[t], 
 x[0] == x0, y[0] == y0}, {x[t], y[t]}, {t, 0, 1000}];
Plot[Evaluate[ReplaceAll[{x[t], y[t]}, %]], {t, 0, 0.01}]

NDSolve[{x'[t] == a x[t] + b y[t], x[t]+y[t]==1, 
 x[0] == x0, y[0] == y0}, {x[t], y[t]}, {t, 0, 1000}];
Plot[Evaluate[ReplaceAll[{x[t], y[t]}, %]], {t, 0, 0.01}]

NDSolve[{y'[t] == c x[t] + d y[t], x[t]+y[t]==1,
 x[0] == x0, y[0] == y0}, {x[t], y[t]}, {t, 0, 1000}];
Plot[Evaluate[ReplaceAll[{x[t], y[t]}, %]], {t, 0, 0.01}]

请注意,要使其中任何一个起作用,您将需要指定常量a,b,c,d,x0,y0(并且您需要为第二个x0y0第三个求和为 1)。

正如你所问的,我正在使用 NDSolve,尽管其中任何一个也可以通过DSolve(有或没有常量的固定值)来解决:

Clear[a,b,c,d,x0,y0];

DSolve[{x'[t] == a x[t] + b y[t], y'[t] == c x[t] + d y[t], 
 x[0] == x0, y[0] == y0}, {x[t], y[t]}, t]

DSolve[{x'[t] == a x[t] + b y[t], x[t]+y[t]==1, 
 x[0] == x0, y[0] == 1-x0}, {x[t], y[t]}, t]

DSolve[{y'[t] == c x[t] + d y[t], x[t]+y[t]==1,
 x[0] == x0, y[0] == 1-x0}, {x[t], y[t]}, t]

但是,非常重要的是要注意,对于三个函数方程和两个函数,您不能同时求解所有这些方程。这就是为什么生成的三张图NDSolve和生成的三组解DSolve会有所不同。请注意,在后两种情况下,我必须通过给出1-x0初始值来明确强制 x[0]+y[0]=1。

因为你有三个方程,你要么有一个不一致的系统,要么是超定的。如果您尝试同时解决所有三个问题,Mathematica 将给出错误:

NDSolve::overdet:
"There are fewer dependent variables, {x[t],y[t]}, than equations, so the
system is overdetermined."

(类似的DSolve。)

我的选择a=1,b=2,c=3,d=4使它不一致(因此我的三个图表不同意)。你可以弄乱Manipulate一些东西来看看它们是如何关联的,但是我有两种基本的方法来确定如何(如果有的话)强制解决方案对所有三个方程保持一致:

  1. 只是区分1=x+y。这给了你0=x'+y'这意味着你有整个限制ax+by=cx+dy,换句话说y=(c-a)/(b-d)x。这实际上与您原来的约束不一致x+y=1

  2. 考虑 DEsx'=ax+by和的线性系统y'=cx+dy。为了使解成为一条直线,它必须在特征向量的方向上并且必须通过原点。您的约束x+y=1是一条不符合这些标准的线,因此两个 DE 的解都不会满足您的代数方程。

PS我认为您的原始代码中有错字。你可能想要{t,0,10000}是的?

于 2014-08-25T17:19:21.280 回答
-1

尝试例如:

DSolve[{x'[t] == A x[t] + B (1 - x[t]), x[0] == 2/10}, x[t], t]
于 2013-06-11T13:17:08.517 回答