ODE 求解器中的已知函数答案

【问题标题】：Known function in ODE solversODE 求解器中的已知函数
【发布时间】：2014-05-14 22:47:18
【问题描述】：

我是 Matlab 的新手，很抱歉问了一个平庸且可能已经回答的问题，但我无法解决它。

我要解决这个 ODE

g''(x)+f(x)g'(x)

g(0)=0 和 g'(0)=0.5

f(x) 是已知向量，所以我尝试了

xspan = linspace(0,10,100);
y0 = [0,0.5];
%f=f(xspan), known function, for example f=1./xspan;
[x,y] = ode45('odefun',xspan,y0);

对于odefun

function dy=odefun(x,y)

global f  

dy(1) = y(2);
dy(2) = -f.*y(2);

dy    = dy(:);
end

当然不行。如何将向量 f 的值传递给求解器？

编辑

感谢您的帮助，但这对我不起作用。也许我有一个旧的 Matlab 版本（7.11.0（R2010B））？我试着按照你说的做，但这是我在命令窗口中得到的

???  In an assignment  A(I) = B, the number of elements in B and
 I must be the same.

Error in ==> odefun at 3
dy(2) = -f.*y(2);

Error in ==> @(x,y)odefun(x,y,f)


Error in ==> odearguments at 109
f0 = feval(ode,t0,y0,args{:});   % ODE15I sets args{1} to yp0.

Error in ==> ode45 at 173
[neq, tspan, ntspan, next, t0, tfinal, tdir, y0, f0, odeArgs, odeFcn, ...

Error in ==> trying at 7
[x,y] = ode45(fun,xspan,y0);

这是正在尝试的代码。m

xspan=linspace(0,10,100);
y0=[0 0.5];
f=1./xspan;
[x,y] = ode45(@(x,y)odefun(x,y,f),xspan,y0);

还有奥德芬

function dy=odefun(x,y,f)
dy(1) = y(2);
dy(2) = -f.*y(2);
dy    = dy(:);
end

【问题讨论】：

Re：您的编辑（一个新问题）：您仍然必须满足编写 ODE 函数的所有常规规则。看起来您正在为f 传递一个向量，然后尝试将其设置为等于dy 的单个元素：dy(2) = -f.*y(2);。 dy(2) = -f(1).*y(2); 或 dy(2) = -f(end).*y(2); 会起作用。您是否可能尝试提取ode45 在每个时间步上使用的步长dt 以在odefun 中使用它？如果是这样，这种方法就行不通了。
我没想到，f是向量，我当然做不到。对不起这个错误。顺便说一句，我不想做 f(1)*y(2) 或 f(end)，我需要在 ode45 求解器中使用的每个步骤对 f 进行评估。那么，我应该按照建议找到步长，然后在这些步骤中插入 f 吗？我怎么能那样做？非常感谢您的帮助！
没有办法获得 ode45 在运行时在内部使用的步长（xspan 只是一组输出点，不一定是求解器所采取的真实步长 – @987654321 @ 可能会有所帮助）。你可以得到自变量x，但不能得到步长。 ODE 不应该是步长的函数，只有当前的自变量和因变量以及任何参数。您可能正在尝试一些数字和数学上不明智的事情。
其实我不需要知道ode45的步长，我需要解微分方程g''(x)+f(x)g'(x)=0 ，而且我不知道如何传递 f(x) 的值（我之前评估过的）。我知道 xspan 向量并不代表 ode45 使用的真实步骤，但我需要在这些步骤中评估 f(x) 以求解上面的方程，为此我考虑了插值。那么，我有没有办法做到这一点？如果我可以得到自变量 x，那么我可以使用插值来获取 ode 的每个“x”的 f(x) 值吗？如果没有，有没有办法解决？
odefun 的第一个输入是 x。

标签： matlab ode

【解决方案1】：

您似乎错误地使用了global 变量（f 也需要在您的主函数中声明为global），但您不应该一开始就使用全局变量。 You should rarely, if ever, use globals. 您还应该将积分函数 odefun 指定为 anonymous function 而不是现代 Matlab 中的字符串。以下是修改代码以解决这些问题的方法：

...
f = ...
[x,y] = ode45(@(x,y)odefun(x,y,f),xspan,y0);

而您的集成函数，可以是sub-function 或在单独的文件中，只需要为参数f 带一个额外的参数：

function dy=odefun(x,y,f)
dy(1) = y(2);
dy(2) = -f.*y(2);
dy    = dy(:);

这些方法也应该比使用字符串和全局更快。代码@(x,y)odefun(x,y,f) 创建了一个匿名函数x 和y（ode45 需要一个正好接受两个输入的函数），它从当前作用域捕获f 的值。这在计算机科学中也称为closure。您可以等效地创建一个函数句柄并将其传入：

...
f = ...
fun = @(x,y)odefun(x,y,f)
[x,y] = ode45(fun,xspan,y0);

另请参阅 The MathWorks 的 this article on parametrizing functions，了解有关此通用技术的更多信息。

【讨论】：

感谢您的关注，但我仍有问题。我编辑了这个问题，如果可以的话，我会很感激你的帮助

【解决方案2】：

我解决了，没那么难，我只需要在自变量x的每个值处插入f。这样代码就可以工作了

function dy=odefun(x,y,f,xspan)
f=interp1(xspan,f,x);
dy(1) = y(2);
dy(2) = -f.*y(2);
dy    = dy(:);
end

【讨论】：