Количество событий ODE изменяется при разных условиях простого ODE - PullRequest
0 голосов
/ 04 марта 2019

У меня есть следующий простой ODE:

dx/dt=-1

При начальном условии x (0) = 5 меня интересует, когда x (t) == 1.Итак, у меня есть следующая функция событий:

function [value,isterminal,direction] = test_events(t,x)
    value = x-1;
    isterminal = 0;
    direction = 0;
end

Это должно создать событие при t = 4.Однако, если я запускаю следующий код, я получаю два события, одно в момент времени t = 4 и одно в соседнем месте t = 4 + 5.7e-14:

options = odeset('Events',@test_events);
sol = ode45(@(t,x)-1,[0 10],5,options);
fprintf('%.16f\n',sol.xe)
% 4.0000000000000000
% 4.0000000000000568

Если я запускаю похожие коды, чтобы найтикогда x (t) == 0 или x (t) == - 1 (значение = x; или значение = x + 1; соответственно), у меня происходит только одно событие.Почему это генерирует два события?

ОБНОВЛЕНИЕ: Если структура параметров изменяется на следующее:

options = odeset('Events',@test_events,'RelTol',1e-4);

... тогда ODE возвращает только одно событие при t = 4 + 5,7.е-14.Если RelTol установлен на 1e-5, он возвращает одно событие при t = 4.Если для RelTol установлено значение 1e-8, оно возвращает те же два события, что и по умолчанию (RelTol = 1e-3).Кроме того, изменение начального условия с x (0) = 5 на x (0) = 4 приводит к одному событию, но установка x (0) = 4 и «RelTol» = 1e-8 приводит к двум событиям.

ОБНОВЛЕНИЕ 2: Наблюдая выходные данные sol.x и sol.y (t и x соответственно), время прогрессирует как целые числа [0 1 2 3 4 5 6 7 ...], а x прогрессирует как целые числа вплоть до x (t = 5) примерно так: [5 4 3 2 1 1.11e-16 -1,000 -2,000 ...].Это указывает на то, что между t = 4 и t = 5 возникает нечто, что создает «удар» в решении ODE.Почему?

1 Ответ

0 голосов
/ 07 марта 2019

Одно предположение, которое может объяснить, как ошибки округления могут возникнуть в этой простой задаче: решение интерполируется между внутренними шагами с использованием оценок k_n функции производных ODE, также называемой "плотным выходом".Теоретическая форма:

b_1(u)k_1 + b_2(u)k_2 + ...b_s(u)k_s

, где 0 <= u<= 1 это параметр в интервале между внутренними точками, то есть t = (1-u)*t_k+u*t_{k+1}.

Полиномы коэффициентов нетривиальны.Хотя в этом примере все k_i=1 являются постоянными, при оценке суммы b_1(u)+...+b_s(u) могут накапливаться ошибки округления, которые становятся видимыми в значении решения, близком к корню, даже если y_k и y_{k+1} являются точными.В этом диапазоне накопленного шума с плавающей запятой значение может колебаться вокруг корня, что приводит к обнаружению нескольких пересечений нуля.

...