无法使用 Octave 求解简单的 ODE

Question

我是 Octave 的新手，所以在转向更复杂的项目之前，我试图让一些简单的示例发挥作用。

我正在尝试解决 ODE 问题dy/dx = a*x+b，但没有成功。这是代码：

%Funzione retta y = a*x + b. Ingressi: vettore valori t; coefficienti a,b
clear all;
%Inizializza argomenti
b = 1;
a = 1;
x = ones(1,20);
function y = retta(a, x, b) %Definisce funzione
y = ones(1,20);
y = a .* x .+ b;
endfunction
%Calcola retta
x = [-10:10];
a = 2;
b = 2;
r = retta(a, x, b)
c = b;
p1 = (a/2)*x.^2+b.*x+c  %Sol. analitica di dy/dx = retta %
plot(x, r, x, p1);
% Risolve eq. differenziale dy/dx = retta %
y0 = b; x0 = 0;
p2 = lsode(@retta, y0, x)

输出是：

retta3code
r =

 -18  -16  -14  -12  -10   -8   -6   -4   -2    0    2    4    6    8   10 12   14   16   18   20   22

p1 =

Columns 1 through 18:

82    65    50    37    26    17    10     5     2     1     2     5    10    17    26    37    50    65

Columns 19 through 21:

82   101   122

error: 'b' undefined near line 9 column 16
error: called from:
error:   retta at line 9, column 4
error: lsode: evaluation of user-supplied function failed
error: lsode: inconsistent sizes for state and derivative vectors
error:   /home/fabio/octave_file/retta3code.m at line 21, column 4

所以，该函数retta第一次运行正常，但在lsode. 为什么会发生这种情况？需要更改哪些内容才能使代码正常工作？

Answer 1

不知怎的，你仍然错过了故事的一些重要部分。要求解 ODE，y'=f(y,x)您需要定义一个函数

function ydot = f(y,x)

其中ydot与 具有相同的维度y，两者都必须是向量，即使 f 它们的维度为 1。x是标量。由于某些传统原因，lsode（在多个求解器包中使用的 FORTRAN 代码）更喜欢较少使用的顺序(y,x)，在大多数教科书和其他求解器中您都可以找到该顺序(x,y)。

然后通过调用的ylist样本点获取解决方案样本xlist

ylist = lsode("f", y0, xlist)

哪里xlist(1)是初始时间。

的内部结构f与样本列表list及其大小无关。这是一个单独的问题，您可以使用多重评估来计算精确的解决方案，例如

yexact = solexact(xlist)

要传递参数，请使用匿名函数，例如

function ydot = f(y,x,a,b)
    ydot = [ a*x+b ]
end

a_val = ...
b_val = ...
lsode(@(y,x) f(y,x,a_val, b_val), y0, xlist)