Раздел: Документация
0 ... 107 108 109 110 111 112 113 ... 365 F = [Р*у<1> - Р*у(1)*у(2] -R*y<2> + г*у(1)*у(2)]; end % вложенная функция для вычисления матрицы Якоби, зависящей от параметров function J = JLotVolPar(t, у) J = (P - p*y<2) - p*y(l) r*y<2) - R + r*y(l)]; end; У0 = [1000; 1100]; % задание начальных условий % формирование управляющей структуры с матрицей Якоби options = odeset(1Jacobian1, SJLotVolPar); % вызов солвера [Т, Y] = ode23s(@LotVolPar, [0 100], У0, options); % построение графика решения figure plot(Y(:, 1), Y(:, 2)) end Системы, не разрешенные относительно производной, дифференциально-алгебраические уравнения До спх пор мы рассматривали задачу Коши для систем дифференциальных уравнений, разрешенных относительно производных Y = f(t,Y). MATLAB позволяет решать системы дифференциальных уравнений, заданных в неявной форме F(t, Y, г") = 0. Важным частным случаем таких систем являются системы с матрицей масс M(t, Y)Y=F(t,Y), причем М (г, К) может быть как невырожденной, так и вырожденной. Для решения систем вида Ft, Y, К) = 0 служит солвер odei5i, а системы с матрицей масс могут быть решены любым из солверов. Заметим, что область применения солвера ode23s ограничена только постоянными матри- цами, а в случае вырожденной матрицы необходимо прибегать к солверам odelSs иode23t. Вырожденная матрица масс соответствует системе дифференциально-алгебраических уравнений, которая наряду с дифференциальными у = f(t, у, и) содержит алгебраические связи g{t, у, н) = 0. Здесь fug — заданные вектор-функции, а у и и — искомые вектор-функции. Такая система может быть решена в MATLAB при условии, что ее индекс равен единице, т.е. матрица (Э,-/Эи;-) невырождена. При решении дифференциально-алгебраических уравнений применяется их сведение к системе М Y ~F{t, У) с вырожденной матрицей масс М = 1 о* о о где / — единичная матрица, размер которой совпадает с длиной вектора неизвестных у, а остальные блоки нулевые и Г = Если дифференциальные уравнения в системе дифференциально-алгебраических уравнений не разрешены относительно производной, то для решения такой системы следует применять солвер odel5i. Рассмотрим сначала решение систем с вырожденной матрицей масс на часто приводимом примере, опубликованном Робертсоном, для иллюстрации решения жестких уравнений — нелинейную задачу химической кинетики. Одно из дифференциальных уравнений системы ()з() = 3107у2(г)2) заменено уравнением баланса. Л#(/) = -0.04Л(/) + 104л(/)л(,); л(0=о-04л(0-ю4л(Ь(0-з-ю7й(/)2; (61) я()+л(0+л()-1=(> /е [0,1001, У,(0) = 1. у2(0) = 0. л(0) = 0. Примечание С некоторыми отличиями эта же задача приведена в демонстрационном примере из пакета MATLAB (функция ihbldae). В нашем случае матрица масс постоянна и имеет вид: г1 О о\ Г = О 1 О ООО Для иллюстрации применения солвера odel5s при решении задач указанного класса написана файл-функция (листинг 6.24) без аргументов, с использованием вложенной функции для вычисления вектора F(t, Y). График каждой функции выводится на свои оси, при этом автоматически масштабируется ось ординат для компоненты у2 • Листинг 6.24. Решение системы дифференциальных уравнений с матрицей масс function exaraple dae % Вложенная функция вычисления вектора f[x, у) function resid = rob(t, у) resid = zeros(3, 1); resid (1) = -0.04*у{1) + le4*y(2) *y(3) ,- resid (2) = +0.04*y(l) - le4*y(2)*y(3) - 3e7*y(2)"2; resid (3) = y(l) + y(2) + y(3) - 1; end % начальные условия. yO =0; 0]; % матрица масс m = diagf[1, 1, 0]); % интервал интегрирования. xtime = [0 100] ; % задание управляющей структуры (точность) options = odeset ( RelTol , 1.: 4, mass, m) ; % вызов солвера ft,у] = odel5s(@rob, xtime, yO, options); I вывод графиков решений subplot(3, 1, 1); plot(t, y(:, 1)); grid on legend (y {l>); subplot(3, 1, 2); plotft, y(:, 2)); 0 ... 107 108 109 110 111 112 113 ... 365
|
