Эта тема показывает, как ограничить решение ОДУ быть неотрицательной. Внушительная неотрицательность не всегда тривиальна, но иногда это необходимо из-за физической интерпретации уравнений или из-за природы решения. Необходимо только наложить это ограничение на решение при необходимости, такой как в случаях, где интегрирование перестало работать без него, или где решение было бы неподходящим.
Если определенные компоненты решения должны быть неотрицательными, то используйте odeset
установить NonNegative
опция для индексов этих компонентов. Эта опция не доступна для ode23s
ode15i
, или для неявных решателей (ode15s
ode23t
ode23tb
) примененный проблемы с большой матрицей. В частности, вы не можете наложить ограничения неотрицательности на проблему ДАУ, которая обязательно имеет сингулярную большую матрицу.
Рассмотрите задачу с начальными значениями
решенный на интервале с начальным условием . Решение этого ОДУ затухает, чтобы обнулить. Если решатель производит отрицательное значение решения, то он начинает отслеживать решение ОДУ через это значение, и расчет в конечном счете перестал работать, когда расчетное решение отличается к . Используя NonNegative
опция предотвращает этот отказ интегрирования.
Сравните аналитическое решение к решению ОДУ с помощью ode45
без дополнительных опций, и одной с NonNegative
опция установлена.
ode = @(t,y) -abs(y); % Standard solution with |ode45| options1 = odeset('Refine',1); [t0,y0] = ode45(ode,[0 40],1,options1); % Solution with nonnegative constraint options2 = odeset(options1,'NonNegative',1); [t1,y1] = ode45(ode,[0 40],1,options2); % Analytic solution t = linspace(0,40,1000); y = exp(-t);
Постройте эти три решения для сравнения. Внушительная неотрицательность крайне важна, чтобы помешать решению повернуть прочь к .
plot(t,y,'b-',t0,y0,'ro',t1,y1,'k*'); legend('Exact solution','No constraints','Nonnegativity', ... 'Location','SouthWest')
Другим примером проблемы, которая требует неотрицательного решения, является проблема колена, закодированная в файле в качестве примера kneeode
. Уравнение
решенный на интервале с начальным условием . Параметр обычно берется, чтобы удовлетворить , и эта проблема использование . Решение этого ОДУ подходы для и для . Однако вычисление числового решения с допусками по умолчанию показывает, что решение следует изоклина для целого интервала интегрирования. Внушительные ограничения неотрицательности приводят к правильному решению.
Решите задачу колена с и без ограничений неотрицательности.
epsilon = 1e-6; y0 = 1; xspan = [0 2]; odefcn = @(x,y,epsilon) ((1-x)*y - y^2)/epsilon; % Solve without imposing constraints [x1,y1] = ode15s(@(x,y) odefcn(x,y,epsilon), xspan, y0); % Impose a nonnegativity constraint options = odeset('NonNegative',1); [x2,y2] = ode15s(@(x,y) odefcn(x,y,epsilon), xspan, y0, options);
Постройте решения для сравнения.
plot(x1,y1,'ro',x2,y2,'b*') axis([0,2,-1,1]) title('The "knee problem"') legend('No constraints','Non-negativity') xlabel('x') ylabel('y')
[1] Шемпин, L.F., С. Томпсон, Дж.А. Кирженка и Г.Д. Бирн, "Неотрицательные решения ОДУ", Прикладная математика и Издание 170, 2005 Расчета, стр 556-569.