Этот пример показывает, как решить дифференциальное уравнение, представляющее модель хищника/добычи, использующую и ode23 и ode45. Эти функции для числового решения обыкновенных дифференциальных уравнений с помощью переменного размера шага методы интегрирования Рунге-Кутта. ode23 использует простую 2-ю и 3-ю пару порядка формул для средней точности, и ode45 использует 4-ю и 5-ю пару порядка для более высокой точности.
Считайте пару обыкновенных дифференциальных уравнений первого порядка известной как уравнения Лотки - Вольтерры или модель добычи хищника:
Переменные
и
мера размеры добычи и популяций хищников, соответственно. Квадратичный перекрестный термин составляет взаимодействия между разновидностями. Увеличение населения добычи, когда никакие хищники не присутствуют, и уменьшения популяции хищников, когда добыча недостаточны.
Чтобы моделировать систему, создайте функцию, которая возвращает вектор - столбец производных состояния, учитывая состояние и временные стоимости. Эти две переменные
и
могут быть представлены в MATLAB как первые два значения в векторном y. Точно так же производные являются первыми двумя значениями в векторном yp. Функция должна принять значения для t и y и возвратить значения, произведенные уравнениями в yp.
yp(1) = (1 - alpha*y(2))*y(1)
yp(2) = (-1 + beta*y(1))*y(2)
В этом примере уравнения содержатся в файле под названием lotka.m m. Этот файл использует значения параметров
и
.
type lotkafunction yp = lotka(t,y) %LOTKA Lotka-Volterra predator-prey model. % Copyright 1984-2014 The MathWorks, Inc. yp = diag([1 - .01*y(2), -1 + .02*y(1)])*y;
Используйте ode23, чтобы решить дифференциальное уравнение, заданное в lotka на интервале
. Используйте начальное условие того,
так, чтобы популяции хищников и добычи были равны.
t0 = 0; tfinal = 15; y0 = [20; 20]; [t,y] = ode23(@lotka,[t0 tfinal],y0);
Постройте график получившегося населения против времени.
plot(t,y) title('Predator/Prey Populations Over Time') xlabel('t') ylabel('Population') legend('Prey','Predators','Location','North')
Теперь постройте график населения друг против друга. Получившийся график плоскости фазы делает циклическое отношение между населением очень ясным.
plot(y(:,1),y(:,2)) title('Phase Plane Plot') xlabel('Prey Population') ylabel('Predator Population')
Решите систему во второй раз с помощью ode45 вместо ode23. Решатель ode45 занимает больше времени у каждого шага, но это также делает большие шаги. Тем не менее, вывод ode45 сглажен, потому что по умолчанию решатель использует непрерывную дополнительную формулу, чтобы произвести вывод в четырех равномерно распределенных моментах времени в промежутке каждого сделанного шага. (Можно настроить число точек с опцией 'Refine'.) Строят график обоих решений для сравнения.
[T,Y] = ode45(@lotka,[t0 tfinal],y0); plot(y(:,1),y(:,2),'-',Y(:,1),Y(:,2),'-'); title('Phase Plane Plot') legend('ode23','ode45')
Результаты показывают, что решение дифференциальных уравнений с помощью различных численных методов может произвести немного отличающиеся ответы.