有没有办法绘制一个带有方程（来自 ODE 系统）的函数 f？

Question

我正在尝试为我的微分方程课解决这个作业，但遇到了一些问题。

鉴于此 ODE 系统：

dy(1) = y(2) − y(3);
dy(2) = −3*y(2) − 5*sin(ω*y(1));
dy(3) = y(1)*y(2);

omega=3。初始值：y(1)=0; y(2)=4; y(3)=1;

在t=0 和t=20 之间绘制f：

f(t) = |2*cos(y1(t)) + y2(t)|

我首先尝试使用以下代码使用 ODE45 模拟方程组：

在编辑器上：

function dy=modelo10(t,y)
    global omega

    dy = zeros(3,1);
    dy(1)= y(2) - y(3);
    dy(2)= -3*y(2) - 5*sin(omega*y(1));
    dy(3)= y(1)*y(2);

end

在命令选项卡上：

>>global omega;
>>omega = 3;
>>[T,Y] = ode45(@modelo10, [0,20], [0,4,1]); 
%% I assign the function to variable f
>>f = @(t) 2.*cos(y(:,1)) + y(:,2);
%% And finally plot it with the values t
>>fplot(f, [0,20]);

我得到一个带有以下错误的空白图表：

Warning: Function behaves unexpectedly on array inputs. To improve 
performance, properly
vectorize your function to return an output with the same size and 
shape as the input
arguments. 
  In matlab.graphics.function.FunctionLine>getFunction
  In matlab.graphics.function.FunctionLine/updateFunction
  In matlab.graphics.function.FunctionLine/set.Function_I
  In matlab.graphics.function.FunctionLine/set.Function
  In matlab.graphics.function.FunctionLine
  In fplot>singleFplot (line 241)
  In fplot>@(f)singleFplot(cax,{f},limits,extraOpts,args) (line 196)
  In fplot>vectorizeFplot (line 196)
  In fplot (line 166) 
Warning: Error updating FunctionLine.

The following error was reported evaluating the function in 
FunctionLine update:
Non-scalar in Uniform output, at index 1, output 1.
Set 'UniformOutput' to false.

Warning: Error updating FunctionLine.

The following error was reported evaluating the function in 
FunctionLine update:
Non-scalar in Uniform output, at index 1, output 1.
Set 'UniformOutput' to false

所以我的问题是绘制f 的步骤（或命令）是什么？

我正在使用 MATLAB R2019a。

Answer 1

你的第一个问题是定义f：

[T,Y] = ode45(@modelo10, [0,20], [0,4,1]); 
f = @(t) 2.*cos(y(:,1)) + y(:,2);

请注意，第一行定义了一个变量Y，但在第二行中您使用了y。这些是不同的变量。这可能对您有用，因为您在工作区中定义了一个变量 y，但它做错了。

接下来，您的f 是一个函数f(t)，但它的定义不是使用输入t。它总是输出相同的东西。 fplot 需要一个函数输出与其输入大小相同的数组，这是您收到的错误消息。

但是，您不需要在此处定义函数，您可以直接计算 f 的所有值并使用 plot 绘制它们：

[t,y] = ode45(@modelo10, [0,20], [0,4,1]);
f = 2.*cos(y(:,1)) + y(:,2);
plot(t,f)

---

另外，我想建议你不要使用global。这里不需要，可以将omega定义为modelo10的输入参数：

function dy = modelo10(t,y,omega)
   dy = zeros(3,1);
   dy(1) = y(2) - y(3);
   dy(2) = -3*y(2) - 5*sin(omega*y(1));
   dy(3) = y(1)*y(2);
end

然后调用ode45如下：

[t,y] = ode45(@(t,y)modelo10(t,y,omega), [0,20], [0,4,1]);

这里，@(t,y)modelo10(t,y,omega) 定义了一个匿名函数，它带有 omega 的值。这个匿名函数有两个输入参数（t 和 y），这是 ode45 所要求的。

最后，您可以通过在一行中定义 modelo10 来简化代码：

function dy = modelo10(t,y,omega)
   dy = [ y(2) - y(3); -3*y(2) - 5*sin(omega*y(1)); y(1)*y(2) ];
end

您现在可以这样做：

modelo10 = @(t,y) [ y(2) - y(3); -3*y(2) - 5*sin(omega*y(1)); y(1)*y(2) ];
[t,y] = ode45(modelo10, [0,20], [0,4,1]);

Answer 2

global omega;
omega = 3;

[T,Y] = ode45(@modelo10, [0,20], [0,4,1]);

• 从这里 Y 是 T 的函数 --> 不需要定义新函数 Y = Y(T) for T = [0,20]，这不再是变量而是实数 Y = [Y1(0),..., Y1(20); Y2(0),..., Y2(20); Y3(0),..., Y3(20)]

• 我将函数分配给变量 f

  f = @(t) 2.*cos(y(:,1)) + y(:,2) --> f = abs( 2.*cos(Y(:,1)) + Y(:,2));
  

  小写 y 和大写 Y 是两个不同的变量。 matlab中的绝对值函数是abs

• 最后用值 t 绘制它

  Use plot instead of fplot
  

  实际上你不需要为函数 f 指定输入值，因为它有 在计算 ode45 时已经包含在内，而 fplot 需要一个键 参数，即函数表达式

  (function handle --> function defined as follow f = @(x)_____). 
  

  输入可以额外添加为第二个参数。

  这是可选的；默认情况下，fplot 计算 [0 5] 范围内的函数

• 这里的f是实数不是函数表达式，它已经被计算过了

      plot(T, f);
  

函数 f 绘制 T 范围从 0 到 20