Этот пример показывает, как использовать MATLAB®, чтобы сформулировать и решить несколько различных типов дифференциальных уравнений. MATLAB предлагает несколько числовых алгоритмов, чтобы решить большое разнообразие дифференциальных уравнений:
Задачи с начальными значениями
Краевые задачи
Дифференциальные уравнения с запаздывающим аргументом
Дифференциальные уравнения с частными производными
vanderpoldemo является функцией, которая задает уравнение Ван дер Поля
type vanderpoldemofunction dydt = vanderpoldemo(t,y,Mu) %VANDERPOLDEMO Defines the van der Pol equation for ODEDEMO. % Copyright 1984-2014 The MathWorks, Inc. dydt = [y(2); Mu*(1-y(1)^2)*y(2)-y(1)];
Уравнение записано как система двух обыкновенных дифференциальных уравнений первого порядка (ОДУ). Эти уравнения оценены для различных значений параметра . Для более быстрого интегрирования необходимо выбрать соответствующий решатель на основе значения .
Для , любой из решателей ОДУ MATLAB может решить уравнение Ван дер Поля эффективно. Решатель ode45 является одним таким примером. Уравнение решено в области с начальными условиями и .
tspan = [0 20]; y0 = [2; 0]; Mu = 1; ode = @(t,y) vanderpoldemo(t,y,Mu); [t,y] = ode45(ode, tspan, y0); % Plot solution plot(t,y(:,1)) xlabel('t') ylabel('solution y') title('van der Pol Equation, \mu = 1')

Для больших значений , проблема становится жесткой. Эта метка для проблем, которые сопротивляются попыткам, которые будут оценены с обычными методами. Вместо этого специальные численные методы необходимы для интеграции FAST. ode15s, ode23s, ode23t и функции ode23tb могут решить жесткие проблемы эффективно.
Это решение уравнения Ван дер Поля для использование ode15s с теми же начальными условиями. Необходимо протянуть отрезок времени решительно к смочь видеть периодическое перемещение решения.
tspan = [0, 3000]; y0 = [2; 0]; Mu = 1000; ode = @(t,y) vanderpoldemo(t,y,Mu); [t,y] = ode15s(ode, tspan, y0); plot(t,y(:,1)) title('van der Pol Equation, \mu = 1000') axis([0 3000 -3 3]) xlabel('t') ylabel('solution y')

bvp4c и bvp5c решают краевые задачи для обыкновенных дифференциальных уравнений.
Функции, взятой в качестве примера, twoode записали дифференциальное уравнение как систему двух ОДУ первого порядка. Дифференциальное уравнение
type twoodefunction dydx = twoode(x,y) %TWOODE Evaluate the differential equations for TWOBVP. % % See also TWOBC, TWOBVP. % Lawrence F. Shampine and Jacek Kierzenka % Copyright 1984-2014 The MathWorks, Inc. dydx = [ y(2); -abs(y(1)) ];
Функциональный twobc имеет граничные условия для проблемы: и .
type twobcfunction res = twobc(ya,yb) %TWOBC Evaluate the residual in the boundary conditions for TWOBVP. % % See also TWOODE, TWOBVP. % Lawrence F. Shampine and Jacek Kierzenka % Copyright 1984-2014 The MathWorks, Inc. res = [ ya(1); yb(1) + 2 ];
До вызова bvp4c необходимо обеспечить предположение для решения, которое вы хотите представленный в mesh. Решатель затем адаптирует mesh, когда это совершенствовало решение.
Функция bvpinit собирает исходное предположение в форме, можно передать решателю bvp4c. Для сетки [0 1 2 3 4] и постоянных предположений и , вызов bvpinit:
solinit = bvpinit([0 1 2 3 4],[1; 0]);
С этим исходным предположением можно решить проблему с bvp4c. Оцените решение, возвращенное bvp4c в некоторых точках с помощью deval, и затем постройте получившиеся значения.
sol = bvp4c(@twoode, @twobc, solinit); xint = linspace(0, 4, 50); yint = deval(sol, xint); plot(xint, yint(1,:)); xlabel('x') ylabel('solution y') hold on

Эта конкретная краевая задача имеет точно два решения. Можно получить другое решение путем изменения исходных предположений на и .
solinit = bvpinit([0 1 2 3 4],[-1; 0]); sol = bvp4c(@twoode,@twobc,solinit); xint = linspace(0,4,50); yint = deval(sol,xint); plot(xint,yint(1,:)); legend('Solution 1','Solution 2') hold off

dde23, ddesd и ddensd решают дифференциальные уравнения с запаздывающим аргументом с различными задержками. Примеры ddex1, ddex2, ddex3, ddex4 и ddex5 формируют мини-пример при использовании этих решателей.
Пример ddex1 показывает, как решить систему дифференциальных уравнений
Можно представлять эти уравнения с анонимной функцией
ddex1fun = @(t,y,Z) [Z(1,1); Z(1,1)+Z(2,2); y(2)];
История проблемы (для ) является постоянным:
Можно представлять историю как вектор из единиц.
ddex1hist = ones(3,1);
Двухэлементный вектор представляет задержки системы уравнений.
lags = [1 0.2];
Передайте функцию, задержки, историю решения и интервал интегрирования к решателю как входные параметры. Решатель производит непрерывное решение на целом интервале интегрирования, которое подходит для графического вывода.
sol = dde23(ddex1fun, lags, ddex1hist, [0 5]);
plot(sol.x,sol.y);
title({'An example of Wille and Baker', 'DDE with Constant Delays'});
xlabel('time t');
ylabel('solution y');
legend('y_1','y_2','y_3','Location','NorthWest');
pdepe решает дифференциальные уравнения с частными производными в одной пространственной переменной и время. Примеры pdex1, pdex2, pdex3, pdex4 и pdex5 формируют мини-пример при использовании pdepe.
Эта проблема в качестве примера использует функции pdex1pde, pdex1ic и pdex1bc.
pdex1pde определяет дифференциальное уравнение
type pdex1pdefunction [c,f,s] = pdex1pde(x,t,u,DuDx) %PDEX1PDE Evaluate the differential equations components for the PDEX1 problem. % % See also PDEPE, PDEX1. % Lawrence F. Shampine and Jacek Kierzenka % Copyright 1984-2014 The MathWorks, Inc. c = pi^2; f = DuDx; s = 0;
pdex1ic настраивает начальное условие
type pdex1icfunction u0 = pdex1ic(x) %PDEX1IC Evaluate the initial conditions for the problem coded in PDEX1. % % See also PDEPE, PDEX1. % Lawrence F. Shampine and Jacek Kierzenka % Copyright 1984-2014 The MathWorks, Inc. u0 = sin(pi*x);
pdex1bc настраивает граничные условия
type pdex1bcfunction [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t) %PDEX1BC Evaluate the boundary conditions for the problem coded in PDEX1. % % See also PDEPE, PDEX1. % Lawrence F. Shampine and Jacek Kierzenka % Copyright 1984-2014 The MathWorks, Inc. pl = ul; ql = 0; pr = pi * exp(-t); qr = 1;
pdepe требует пространственной дискретизации x и вектор времен t (в котором вы хотите снимок состояния решения). Решите проблему с помощью сетки 20 узлов и запросите решение в пяти значениях t. Извлеките и постройте первый компонент решения.
x = linspace(0,1,20); t = [0 0.5 1 1.5 2]; sol = pdepe(0,@pdex1pde,@pdex1ic,@pdex1bc,x,t); u1 = sol(:,:,1); surf(x,t,u1); xlabel('x'); ylabel('t'); zlabel('u');
