Решение нежестких дифференциальных уравнений — метод среднего порядка точности
[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), вызвал функции события, являются нулем. В выводе t,y,te,ye,ie]
= ode45(odefun,tspan,y0,options)te является временем события, ye является решением во время события, и ie является индексом инициированного события.
Для каждой функции события задайте, должно ли интегрирование остановиться в нуле и имеет ли направление нулевого пересечения значение. Сделайте это путем установки свойства 'Events' на функцию, такую как myEventFcn или @myEventFcn, и создания соответствующей функции: [value, isterminal, direction] = myEventFcn (t, y). Для получения дополнительной информации смотрите Местоположение События ОДУ.
Простые ОДУ, которые имеют компонент единого решения, могут быть заданы как анонимная функция в вызове решателя. Анонимная функция должна принять два входных параметров (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. и М. В. Рейчелт, “Пакет ODE MATLAB”, SIAM Journal на Научных вычислениях, Издании 18, 1997, стр 1–22.