Matlab：微分方程（ode45）：我可以反转 tspan 以获得更好的初始条件吗？

Question

我正在使用ode45 在 Matlab 中求解/绘制二阶微分方程。我的tspan 是从 0 到 0.25。但是接近零的初始条件定义不明确（斜率趋于无穷大，复数值）。 0.25 附近的条件定义明确（斜率和值均为零）。

问题：

1. 我可以反转tspan，并使用“最终条件”作为初始条件吗？

2. 嗯，我知道我能做到（见下面的代码），我得到的情节看起来像我所期望的那样，但一般来说这是一个有效的事情吗？在这种情况下我很幸运吗？

3. ode45 提供数值解，并不精确。反转tspan 后可能会出现更大的错误吗？

这是我的代码，应该独立运行：

function ReverseTspan()
% solve diff-eq backward from tspan end to tspan start using ode45()
%  - Good initial conditions at the end, but not start.
%  - Is this a valid thing to do?

% clean slate
clc; clear all; close all;

% tspan - reversed!
R = 0.25:-0.001:0;

% initial values
hinits=[0.0000001;0]; 

% solve
[R,v] = ode45(@equ7,R,hinits);

% strip imaginary values (can't plot 'em)
v(find(real(v)~=v)) = NaN;

% plot first column 
plot(R,v(:,1)); 


function vprime = equ7(R,v);
% Solve second order non-linear differential equation 7:
%   v''(R) + 2/R*v'(R) = K_plus/(R^2)*( v^(-1/2) - lamda_plus*(1-v)^(-1/2)
%
% Matlab ode45 only likes first derivatives, so let:  
%    v_1(R) = v(R)
%    v_2(R) = v'(R)
%
% And create a system of first order diff eqs:
%    v_1'(R) = v_2(R)
%    v_2'(R) = -2/R*v_2(R) + K_plus/(R^2)*( v_1(R)^(-1/2) - lamda_plus*(1-v_1(R))^(-1/2)
%
% Constant Parameters:
K_plus = 0.0859;
lambda_plus = 3.7;

% Build result in pieces for easier debugging of problematic terms
int1 = 1 - v(1);
int2 = int1^(-1/2);
int3 = v(1)^(-1/2);
int4 = K_plus/(R^2);
vprime2 = -2/R*v(2);
vprime2 = vprime2 + int4*( int3 - lambda_plus*(int2) );
vprime=[v(2); vprime2 ];

Answer 1

从数字上讲，您永远无法使用该 ODE 从（或到达）R=0 开始。 NaN 将在您的初始条件之后的第一步返回，从而破坏任何未来的结果。可能可以重新调整您的系统以摆脱奇点，但这是一个数学问题。如果您指定有效的初始条件（满足您的 ODE），您也可以尝试在比零大一点的时间开始积分。

1. 我可以反转tspan，并使用“最终条件”作为初始条件吗？

是的。这是一个标准程序。特别是，它经常被用来寻找不稳定的流形和不稳定的平衡。

1. 好吧，我知道我可以做到（见下面的代码），我得到的情节看起来像我所期望的那样，但一般来说这是一个有效的事情吗？在这种情况下我会走运吗？

很难概括所有可能的系统。在某些情况下，由于稳定​​性或数值错误，这可能不起作用。

1. ode45 提供数值解，并不精确。反转tspan后会不会有更大的错误？

我不明白为什么它一定会更大，但这又取决于所使用的系统和求解器。尝试使用其他求解器（例如，ode15s 或 ode113）和/或通过 odeset 更改 absolute and relative tolerances。