Решение нежестких дифференциальных уравнений — метод среднего порядка точности
[t,y] = ode45(odefun,tspan,y0)[t,y] = ode45(odefun,tspan,y0,options)[t,y,te,ye,ie] = ode45(odefun,tspan,y0,options)sol = ode45(___)[t,y] = ode45(odefun,tspan,y0)tspan = [t0 tf], интегрирует систему дифференциальных уравнений от t0 до tf с начальными условиями y0. Каждая строка в массиве решения y соответствует значению, возвращенному в векторе - столбце t.
Все решатели MATLAB® ODE могут решить системы уравнений формы , или проблемы, которые включают большую матрицу, . Решатели все использование подобные синтаксисы. Решатель ode23s только может решить проблемы с большой матрицей, если большая матрица является постоянной. ode15s и ode23t могут решить проблемы с большой матрицей, которая сингулярна, известна как дифференциально-алгебраические уравнения (ДАУ). Задайте большую матрицу с помощью опции Mass odeset.
ode45 является универсальным решателем ОДУ и является первым решателем, который необходимо попробовать за большинство проблем. Однако, если проблема жестка или требует высокой точности, то существуют другие решатели ОДУ, которые могут лучше подходить для проблемы. Смотрите Выбирают ODE Solver для получения дополнительной информации.
[t,y] = ode45(odefun,tspan,y0,options) options, который является созданным использованием аргумента функции odeset. Например, используйте AbsTol и опции RelTol, чтобы задать допуски абсолютной и относительной погрешности или опцию Mass, чтобы обеспечить большую матрицу.
[t,y,te,ye,ie] = ode45(odefun,tspan,y0,options) te является временем события, ye является решением во время события, и ie является индексом инициированного события.
Для каждой функции события задайте, должно ли интегрирование отключить в нуле и имеет ли направление нулевого пересечения значение. Сделайте это путем установки свойства 'Events' на функцию, такую как myEventFcn или @myEventFcn, и создания соответствующей функции: [value, isterminal, direction] = myEventFcn (t, y). Для получения дополнительной информации смотрите Местоположение События ОДУ.
sol = ode45(___) deval, чтобы оценить решение в любой точке на интервале [t0 tf]. Можно использовать любую из комбинаций входных аргументов в предыдущих синтаксисах.
Простые ОДУ, которые имеют компонент единого решения, могут быть заданы как анонимная функция в вызове решателя. Анонимная функция должна принять два входных параметра (t,y), даже если один из входных параметров не используется.
Решите ОДУ
Используйте временной интервал [0,5] и начального условия y0 = 0.
tspan = [0 5]; y0 = 0; [t,y] = ode45(@(t,y) 2*t, tspan, y0);
Постройте график решения.
plot(t,y,'-o')
Уравнение Ван дер Поля является ОДУ второго порядка
где
скалярный параметр. Перепишите это уравнение как систему ОДУ первого порядка путем создания замены
. Получившаяся система ОДУ первого порядка
Функциональный файл vdp1.m представляет использование уравнения Ван дер Поля
. Переменные
и
являются записями y(1) и y(2) двухэлементного вектора, dydt.
function dydt = vdp1(t,y) %VDP1 Evaluate the van der Pol ODEs for mu = 1 % % See also ODE113, ODE23, ODE45. % Jacek Kierzenka and Lawrence F. Shampine % Copyright 1984-2014 The MathWorks, Inc. dydt = [y(2); (1-y(1)^2)*y(2)-y(1)];
Решите ОДУ с помощью функции ode45 на временном интервале [0 20] с начальными значениями [2 0]. Получившийся вывод является вектором - столбцом моментов времени t и массив решения y. Каждая строка в y соответствует времени, возвращенному в соответствующей строке t. Первый столбец y соответствует, и второй столбец к.
[t,y] = ode45(@vdp1,[0 20],[2; 0]);
Постройте график решений для и против t.
plot(t,y(:,1),'-o',t,y(:,2),'-o') title('Solution of van der Pol Equation (\mu = 1) with ODE45'); xlabel('Time t'); ylabel('Solution y'); legend('y_1','y_2')
ode45 только работает с функциями, которые используют два входных параметра, t и y. Однако можно передать в дополнительных параметрах путем определения их вне функции и передачи их в том, когда вы задаете указатель на функцию.
Решите ОДУ
При перезаписи уравнения, когда уступает система первого порядка
odefcn.m представляет эту систему уравнений как функцию, которая принимает четыре входных параметра: t, y, A и B.
function dydt = odefcn(t,y,A,B)
dydt = zeros(2,1);
dydt(1) = y(2);
dydt(2) = (A/B)*t.*y(1);
Решить ОДУ используя решатель ode45. Задайте указатель на функцию, таким образом, что он передает в предопределенных значениях для A и B к odefcn.
A = 1; B = 2; tspan = [0 5]; y0 = [0 0.01]; [t,y] = ode45(@(t,y) odefcn(t,y,A,B), tspan, y0);
Постройте график результатов.
plot(t,y(:,1),'-o',t,y(:,2),'-.')
Рассмотрите следующее ОДУ с зависящими от времени параметрами
Начальное условие
. Функциональный f(t) задан n-1 векторным f, оцененным во времена ft. Функциональный g(t) задан m-1 векторным g, оцененным во времена gt.
Создайте векторы f и g.
ft = linspace(0,5,25); f = ft.^2 - ft - 3; gt = linspace(1,6,25); g = 3*sin(gt-0.25);
Запишите функции с именем myode, который интерполирует f и g, чтобы получить значение зависящих от времени условий в требуемое время. Сохраните функцию в своей текущей папке, чтобы запустить остальную часть примера.
Функция myode принимает, что дополнительные входные параметры оценивают ОДУ на каждом временном шаге, но ode45 только использует первые два входных параметра t и y.
function dydt = myode(t,y,ft,f,gt,g) f = interp1(ft,f,t); % Interpolate the data set (ft,f) at time t g = interp1(gt,g,t); % Interpolate the data set (gt,g) at time t dydt = -f.*y + g; % Evaluate ODE at time t
Решите уравнение по временному интервалу [1 5] с помощью ode45. Задайте функцию с помощью указателя на функцию так, чтобы ode45 использовал только первые два входных параметра myode. Также ослабьте пороги ошибок, используя odeset.
tspan = [1 5]; ic = 1; opts = odeset('RelTol',1e-2,'AbsTol',1e-4); [t,y] = ode45(@(t,y) myode(t,y,ft,f,gt,g), tspan, ic, opts);
Постройте график решения, y, как функция моментов времени, t.
plot(t,y)
Уравнение Ван дер Поля является ОДУ второго порядка
Решите уравнение Ван дер Поля с
использованием ode45. Функциональный vdp1.m поставляется с MATLAB® и кодирует уравнения. Задайте единственный вывод, чтобы возвратить структуру, содержащую информацию о решении, таком как точки оценки и решатель.
tspan = [0 20]; y0 = [2 0]; sol = ode45(@vdp1,tspan,y0)
sol = struct with fields:
solver: 'ode45'
extdata: [1x1 struct]
x: [1x60 double]
y: [2x60 double]
stats: [1x1 struct]
idata: [1x1 struct]
Используйте linspace, чтобы сгенерировать 250 точек в интервале [0 20]. Оцените решение в этих точках с помощью deval.
x = linspace(0,20,250); y = deval(sol,x);
Постройте график первого компонента решения.
plot(x,y(1,:))
Расширьте решение
использования odextend и добавьте результат к исходному графику.
sol_new = odextend(sol,@vdp1,35); x = linspace(20,35,350); y = deval(sol_new,x); hold on plot(x,y(1,:),'r')
ode45 основан на явной формуле Рунге-Кутта (4,5), модификации Дорманда-Принца. Это - одношаговый решатель – в вычислении y(tn), этому нужно только решение в непосредственно предшествующем моменте времени, y(tn-1) [1], [2].
[1] Dormand, J. R. and P. J. Prince, “A family of embedded Runge-Kutta formulae,” J. Comp. Appl. Math., Vol. 6, 1980, pp. 19–26.
[2] Шемпин, L. F. и М. W. Рейчелт, “Пакет ODE MATLAB”, SIAM Journal на Научных вычислениях, Издании 18, 1997, стр 1–22.