【问题标题】:Maple, maximum value of a numeric functionMaple,数值函数的最大值
【发布时间】:2016-03-19 11:23:35
【问题描述】:

我正在尝试使用 Maple 计算具有两个 ODE 的系统解的最大值。我首先解决了系统本身:

> with(DEtools):with(plots):
> a1:=0.00875;a2:=0.075;b1:=7.5;b2:=2.5;d1:=0.0001;d2:=0.0001;g:=4*10^(-8);K1:=5000;K2:=2500;n:=2;m:=2;

> dsol:= dsolve({
diff( x(t), t ) = a1+b1*x(t)^n/(K1^n+x(t)^n)-g*x(t)*y(t)-d1*x(t), 
diff( y(t), t ) = a2+b2*x(t)^m/(K2^m+x(t)^m)-d2*y(t), 
x(0) = 1000, y(0) = 1000}, numeric, output = listprocedure);

> xt:= eval( x(t), dsol );
yt:= eval( y(t), dsol );


> X:=plot(xt(t),t=0..50000,color=blue,legend="x(t)"):
Y:=plot(yt(t),t=0..50000,color=green,legend="y(t)"):
> display([X,Y]);

我在 xt 和 yt 上得到了系统的解,但它们是数值解。所以Maple的maximize()函数不起作用:

> maximize(xt);
> maximize(xt(t),t=0..20000);

是否可以用 Maple 计算数值函数的最大值?

【问题讨论】:

    标签: function numeric maple


    【解决方案1】:

    您的两条曲线xtyt 在您的t=0..50000 范围内都有一个局部最大值,因此您可以直接使用Optimization 包。

    restart;
    with(plots):
    a1:=0.00875: a2:=0.075: b1:=7.5: b2:=2.5: d1:=0.0001:
    d2:=0.0001: g:=4*10^(-8): K1:=5000: K2:=2500: n:=2: m:=2:
    dsol:= dsolve({diff(x(t),t)=a1+b1*x(t)^n/(K1^n+x(t)^n)-g*x(t)*y(t)-d1*x(t),
                   diff(y(t),t)=a2+b2*x(t)^m/(K2^m+x(t)^m)-d2*y(t),
                   x(0)=1000, y(0)=1000}, numeric, output=listprocedure):
    xt:= eval(x(t), dsol):
    yt:= eval(y(t), dsol):
    X:=plot(xt(t), t=0..50000, color=blue, legend="x(t)"):
    Y:=plot(yt(t), t=0..50000, color=green, legend="y(t)"):
    xmax:=Optimization:-Maximize(xt, 0..50000):
    [xmax[2][1],xmax[1]];
    
                   [9460.78688552799, 11193.0618953179]
    
    ymax:=Optimization:-Maximize(yt, 0..50000):
    [ymax[2][1],ymax[1]];
    
                   [21471.8648785947, 19006.6009784691]
    
    display( Y, pointplot([[ymax[2][1],ymax[1]]], symbolsize=20),
             X, pointplot([[xmax[2][1],xmax[1]]], symbolsize=20) );
    

    对于您的简单示例,这非常有效。

    如果您的xtyt 有许多局部最大值,您可以尝试使用method=branchandbound 选项调用Maximize

    然后还有其他方法,您可以使用新的因变量来增强您的 DE 系统,例如 xd(t)=diff(x(t),t)(以及适当的 IC),并让 dsolve/numeric 本身注意到它何时变为零(一个极值点)使用events 设施,或使用fsolve

    【讨论】:

    • 完美!很好的解释和非常有用的代码!!非常感谢!!
    猜你喜欢
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2023-03-20
    • 1970-01-01
    • 2016-01-27
    • 2015-10-29
    • 1970-01-01
    相关资源
    最近更新 更多