89129

Численное решение систем дифференциальных уравнений

Контрольная

Информатика, кибернетика и программирование

В процессе выполнения данной контрольной работе была написана программа Matlab, основной задачей которой является решение системы ОДУ методом Рунге-Кутты 4-5 порядка. Система ОДУ была решена по уравнениям и данным, заданным согласно варианту.

Русский

2015-05-09

392.63 KB

4 чел.

Министерство образования и науки РФ

ПЕНЗЕНСКИЙ ГОСУДАРСТВЕННЫЙ ТЕХНОЛОГИЧЕСКИЙ

УНИВЕРСИТЕТ

Кафедра «Прикладная математика и исследование операций в экономике»

КОНТРОЛЬНАЯ РАБОТА № 1 (13ИС2б)

по дисциплине «Дискретная математика. Методы оптимизации. Численные методы»

на тему «Численное решение систем дифференциальных уравнений»

Вариант № __

Автор работы:              ____________________________

Специальность:                                       _____________________________

Группа:       __________________________

Руководитель:      _____________________________

Работа защищена «__»_____20__г.  Оценка _______________

2014 г.

Содержание

Введение 3

1. ПОСТАНОВКА ЗАДАЧИ 3

2. Теоретические сведения 5

2.1 Метод Рунге–Кутты 5

2.2 Априорный выбор шага интегрирования 5

3. Решение задания 7

Вывод 10

Список литературы 10

Приложение А – Текст программы 11

Приложение Б – Таблица погрешности 12


Введение

В процессе выполнения данной контрольной работе была написана программа Matlab, основной задачей которой является решение системы ОДУ методом Рунге-Кутты 4-5 порядка. Система ОДУ была решена по уравнениям и данным, заданным согласно варианту. Также было реализовано решение системы ОДУ стандартным решателем MATLAB – функцией ode45. Результат ее решения сравнен с результатом написанной программы в точке T/2. Относительные погрешности, полученные в результате сравнения, занесены в таблицу погрешностей. Также был создан видеофайл формата avi с помощью функции VideoWritter, в котором показано движение точки в декартовой системе координат. Данный видеофайл был записан на диск формата DVD-RW.

1. ПОСТАНОВКА ЗАДАЧИ

Цель работы.

Научиться решать системы обыкновенных дифференциальных уравнений в программе Matlab, научиться реализовывать алгоритмы шага интегрирования, строить трехмерный график движения точки в декартовой системе координат и реализовывать график движения в видеофайле.

Задание на контрольную работу.

  1.  Решить заданную систему обыкновенных дифференциальных уравнений (ОДУ) методом Рунге - Кутты 4-5-го порядка, разработав собственную программу в Matlab в виде m-файла, а также решить задачу с помощью решателя Matlab (использовать как эталонное решение).
  2.  В разработанной программе реализовать выбор шага интегрирования по алгоритмам, приведенным в соответствии с заданным вариантом. При решении стандартным решателем Matlab, использовать автоматический шаг.
  3.  Решение, полученное с помощью разработанной программы, сравнить с эталонным решением в точке. Результаты сравнения представить в виде таблицы относительных погрешностей решения. Сделать выводы о точности решения.
  4.  Построить отдельно графики , , , а также трехмерный график движения точки в декартовой системе координат средствами Matlab.
  5.  Создать видеофайл решения задачи, используя функцию VideoWriter: движение точки в трехмерной декартовой системе координат (представить на CD).

Индивидуальное задание.

№ п/п

Система ОДУ

Начальные условия

Граничные условия

Метод выбора шага интегрирования

24

0.0

0.1

0.0

6.0

Априорный


2. Теоретические сведения

2.1 Метод Рунге–Кутты

Методы Рунге — Кутты — важное семейство численных алгоритмов решения обыкновенных дифференциальных уравнений и их систем. Данные итеративные методы явного и неявного приближённого вычисления были разработаны около 1900 года немецкими математиками К. Рунге и М. В. Куттой.

Формально, методом Рунге — Кутты является модифицированный и исправленный метод Эйлера, они представляют собой схемы второго порядка точности. Существуют стандартные схемы третьего порядка, не получившие широкого распространения. Наиболее часто используется и реализована в различных математических пакетах стандартная схема четвёртого порядка. Иногда при выполнении расчётов с повышенной точностью применяются схемы пятого и шестого порядков. Построение схем более высокого порядка сопряжено с большими вычислительными трудностями. Методы седьмого порядка должны иметь по меньшей мере девять стадий, в схему восьмого порядка входит 11 стадий. Хотя схемы девятого порядка не имеют большой практической значимости, неизвестно, сколько стадий необходимо для достижения этого порядка точности. Аналогичная задача существует для схем десятого и более высоких порядков.

Метод Рунге — Кутты четвёртого порядка столь широко распространён, что его часто называют просто методом Рунге — Кутты. Рассмотрим задачу Коши:

Тогда приближенное значение в последующих точках вычисляется по итерационной формуле:

Вычисление нового значения проходит в четыре стадии:

где  — величина шага сетки по .

Этот метод имеет четвёртый порядок точности, то есть суммарная ошибка на конечном интервале интегрирования имеет порядок  (ошибка на каждом шаге порядка ).

2.2 Априорный выбор шага интегрирования

Величина x0 (значение решения в узле, для которого выбирается шаг h), ε, ∆ (допустимые относительная и абсолютная погрешности) считаются заданными. Вначале вычисляются масштабирующие множители αi, i [1 : n], если она задана пользователем, иначе:

Затем вычисляется величина τ:

1/u, 1/v вычисляются интерполированием, ρ = ρ(α) и h = τρ.


3. Решение задания

Программа запускается в файле b2.m. Сначала была объявлена функция a2 с добавлением трех дифференциальных уравнений из файла a2, начальные условия, граничные условия и точность интегрирования. Далее функция запускается с учетом объявленных значений. Также в файле b2.m выполняется графическая часть работы. Сначала строятся графики с приближенными значениями функций f(x), f(y), f(z) по отдельности.

Рисунок 1 – сравнение метода Рунге-Кутты с функцией ode45 по уравнению X(T).

Рисунок 2 – сравнение метода Рунге-Кутты с функцией ode45 по уравнению Y(T).

Рисунок 3 – сравнение метода Рунге-Кутты с функцией ode45 по уравнению Z(T).

Затем была объявлена функция реализующая решение ОДУ стандартным решателем MATLAB (функция ode45). После вычислений был построен график по полученным значениям.

Далее был создан видеофайл формата avi с помощью функции VideoWriter и открывается для того, чтобы вписать в него информацию. В цикле for для всех приближенных значений функций f(x), f(y), f(z) создаем трехмерный график и добавляем точку, которая должна двигаться по полученной трехмерной декартовой системе координат. В цикле были заданы координаты точки относительно каждой оси, по которым она должна двигаться. В конце файла b2.m файл закрывается. Видеофайл с движением точки в декартовой системе координат был записан на диск формата DVD-RW.

Рисунок 4 – декартова система координат движения точки.

В файле a2.m (функции) были заданы три функции, составляющие систему дифференциальных уравнений и алгоритм решения ОДУ стандартным решателем MATLAB (функция ode45). и основной цикл программы численного решения системы ОДУ. Как начальное условие были введены три заданные функции для решения ОДУ, границы диапазона уравнений, точность решения, начальный шаг интегрирования и начальные условия x(0), y(0), x(0), которые записаны в запускающем файле b2.m. Основной цикл while выполняется до тех пор, пока диапазон не сдвинется к верхней границе. Далее идет вложенный цикл for, который работает до тех пор, пока логическое значение истинно. Сначала берется исходный шаг h и по нему высчитывается значение производных функций до производной 4 порядка по формулам:

Далее высчитывалось приближенное значение трех дифференциальных уравнений в последующих точках по формуле: . Далее идет вычисление априорного шага интегрирования и новых значений для следующей итерации.


Вывод

В процессе выполнения контрольной работы была разработана реализация системы обыкновенных дифференциальных уравнений в программе Matlab, алгоритма шага интегрирования, трехмерного графика движения точки в декартовой системе координат и реализация графика движения в видеофайле.

Список литературы

1) Шампайн Л. Ф., Гладвел И., Томпсон С. Решение обыкновенных дифференциальных уравнений с использованием MATLAB: Учебное пособие. 1-е изд. – СПб.: Лань, 2009, 304 с.

2) Чарльз Генри Эдвардс, Дэвид Э. Пенни. Дифференциальные уравнения и краевые задачи: моделирование и вычисление с помощью Mathematica, Maple и MATLAB. 3-е издание. – Диалектика-Вильямс, 2007, 450 с.

3) Е. Р. Алексеев, О. В. Чеснокова Решение задач вычислительной математики в пакетах Mathcad 12, MATLAB 7, Maple 9. Серия: Самоучитель. – М.: НТ Пресс, 2006,496 стр.


Приложение А – Текст программы

Файл a2.m

function [x, y, te, ye] = a2(f,x,h0,y0,e)

function с2 = div(t,Y)

% далее запишем систему ДУ

dy = ones(length(y),1);

dy(1) = y(1)*cos(y(2))-y(3);

dy(2) = sqrt(abs(y(2)^2-y(1)^3))+3*y(3);

dy(3) = y(2)^2*sin(y(1))+0.1;

end

% f - имя функции системы,y0 - начальные условия

% x - интервал времени

% ho - начальный шаг интегрирования

% e - точность

jmax = 1000;

h    = h0;

j    = 1;

while j < jmax

   N    = length(y0);

   H    = max(x);

   n  = ceil((H-x(1))/h);

   y  = zeros(N,n+1);

   x(1) = min(x);

   for k=1:1:N

       y(k,1) = y0(k);

   end

   for i=1:1:n

       x(i+1) = x(i) + h;

   end

   for i=2:1:n+1

       k1 = h.*f(x(i-1),y(:,i-1));

       k2 = h.*f(x(i-1)+h/2,y(:,i-1)+k1/2);

       k3 = h.*f(x(i-1)+h/2,y(:,i-1)+k2/2);

       k4 = h.*f(x(i),y(:,i-1)+k3);

       dy = (k1 + 2*k2 + 2*k3 + k4)/6;

       y(:,i) = y(:,i-1) + dy;

   end

   L(j) = n+1;

   if j ~= 1

      for i3=1:1:N

          M(i3) = max(abs(Y2(i3,:) - Y1(i3,:)));

          M11(i3) = max(abs(y(i3,:)));

          M22(i3) = max(abs(y(i3,:)));

      end   

      M1 = max(M11);

      M2 = max(M22);

      for i1=1:1:N

          for i2=1:1:L(j-1)

              Y2(i1,i2) = y(i1,2*i2-1)/M2;

          end

          for i2=1:1:L(j)

           Y1(i1,i2) = y(i1,i2)/M1;

       end

      end   

   end

   if (j > 1) && (max(M) <= e)

       break;

   end

   h = h/2;

   j = j+1;

end

x = x';

y = y';

te = x(n+1);

ye = y(n+1,:);

end

Файл b2.m

[t,Y,te1,ye1]  = a2(@div,[0 6],0.2,[0, 0.1, 0],0.01);

figure('NumberTitle', 'off', 'Name', 'Метод Рунге-Кутта (X(T))')

title('Метод Рунге-Кутта (X(T))')

legend('X(T)');

plot(t,Y(:,1),'k');         

grid on;

xlabel('T');

ylabel('X');

figure('NumberTitle', 'off', 'Name', 'Метод Рунге-Кутта (Y(T))')

title('Метод Рунге-Кутта (Y(T))')

legend('Y(T)');

plot(t,Y(:,2),'b');         

grid on;

xlabel('T');

ylabel('Y');

figure('NumberTitle', 'off', 'Name', 'Метод Рунге-Кутта (Z(T))')

title('Метод Рунге-Кутта (Z(T))')

legend('Z(T)');

plot(t,Y(:,3))         

grid on;

xlabel('T');

ylabel('Z');

 

[x, y, te2, ye2] = ode45(@div,[ [0 6],[0, 0.1, 0]); % решение через ode45

figure('NumberTitle', 'off', 'Name', 'ode45 (X(T))')

title('ode45 (X(T))')

plot(x, y(:,1))

legend('X(T) ode45');

grid on;

xlabel('T');

ylabel('X');

figure('NumberTitle', 'off', 'Name', 'ode45 (Y(T))')

title('ode45 (Y(T))')

plot(x, y(:,2))

legend('Y(T) ode45');

grid on;

figure('NumberTitle', 'off', 'Name', 'ode45 (Z(T))')

title('ode45 (Z(T))')

plot(x, y(:,3))

legend('Z(T) ode45');

grid on;

 

t=[0 6];                              % границы

[x, y, te, ye] = ode45(@div,t,[0, 0.1, 0]);

figure('NumberTitle', 'off', 'Name', 'Метод Рунге-Кутта (движение точки)')

mov = VideoWriter('rgr.avi');       % создание видеофайла

mov.FrameRate = 25;

open(mov);

for i = 1:length(T)

   plot3(Y(:,1),Y(:,2),Y(:,3),'-k',...  % построение графика

       Y(i,1),Y(i,2),Y(i,3),'*r')

   view(-38+0.1*i,26+0.1*i)

   view( [ 15 , 25 ] )

   xlim([min(Y(:,1)), max(Y(:,1))])

   ylim([min(Y(:,2)), max(Y(:,2))])

   zlim([min(Y(:,3)), max(Y(:,3))])

   grid on

   title(['T=',num2str(T(i),'%1.3f'),])

   % заголовок графика со значением Т

   xlabel('X(T)')

   ylabel('Y(T)')

   zlabel('Z(T)')

   F = getframe(gcf);

   writeVideo(mov,F);    % запись видео

end

close(mov);

Приложение Б – Таблица погрешности

Функция ode45

Собственная функция

Относительная погрешность (%)

X(t)

-1.763379798047849

-1.763118495173552

0.999852

Y(t)

0.004430856890985

0.004528856919895

1.022118

Z(t)

-0.006494568945723

-0.006351795448609

0.978017

Вывод: погрешность собственной функции, реализующей метод Рунге-Кутты, не превышает 1,02 %.


 

А также другие работы, которые могут Вас заинтересовать

24062. Витамины. Этапы нарушений обмена витаминов 81.5 KB
  Витамины не синтезируются в организме или синтезируются в таких количествах которые не достаточны для выполнения функций и поэтому должны поступать в составе пищевых продуктов при резкой недостаточности витаминов в организме развивается характерный симптомокомплекс. Функции витаминов. Нарушение функций витаминов: Нарушение обмена витаминов может быть связано с нарушением всасывания витаминов или их транспорта с кровью. Нарушение образования активной формы кофермента или нарушение синтеза апофермента может привести к нарушению функций...
24063. Тиамин – В1 113.5 KB
  Патология: При недостаточности тиамина наблюдается неврологическое заболевание берибери я не могу. Для берибери характерны мышечная слабость истощение плохая координация периферический неврит спутанность сознания снижение частоты сердечных сокращений и увеличение размеров сердца. Биохимическая диагностика берибери свидетельствует о повышении концентрации пирувата что свидетельствует об участии ТПФ в качестве кофермента в пируватдегидрогеназном комплексе.
24064. Витамин В5(РР) 68.5 KB
  Никотиновая кислота синтезируется из триптофана через кинуренин и оксихинолиновую кислоту. окислении SH2 НАД НАДНН ФАД ФАДН2 КоQ КоQН2 цит b цит с цит а цит а3 О2 Никотинамид синтезируется из триптофана Триптофан кинурениназа Кинуреновая кислота В6 Кинуренин 1 В6 Антраниловая кислота 2 Ксантуреновая кислота Оксикинуренин Оксиантраниловая кислота Никотинамид Хинолиновая кислота Патология обмена витамина В5.
24065. Витамин В2 – рибофлавин 41 KB
  ФАД – участвует в следующих реакциях: Окислительное декарбоксилирование пирувата – входит в состав пируватдегидрогеназного комплекса: СН3СОСООН СН3СОSКоА Окислительное декарбоксилирование кетоглутарата – входит в состав кетоглутаратдегидрогеназного комплекса: НООССН2СН2СОСООН НООССН2СН2СОSКоА В окислении сукцината при СДГ В окислении жирных кислот в митохондриях: RСН2СН2СОSКоА RСН=СНСОSКоА Участие в работе дыхательной цепи Недостаточность рибофлавина проявляется в снижении содержания коферментных форм в тканях. КоА участвует...
24066. Витамин В6 99 KB
  Триптофан кинурениназа Кинуреновая кислота В6 Кинуренин 1 В6 Антраниловая кислота 2 Ксантуреновая кислота Оксикинуренин Оксиантраниловая кислота Никотинамид Хинолиновая кислота В6 входит в состав кинурениназы которая обеспечивает превращение кинуренина в антраниловую и оксикинуринина в оксиантраниловую кислоту реакция 2.
24067. Обмен витамина Н (биотин) 43 KB
  Карбоксилирование ацетилКоА с образованием малонилКоА СН3СОSКоА НООССН2СОSКоА Подготовительным этапом биосинтеза жирных кислот. Карбоксилирование пропионилКоА с образованием метилмалонилКоА: СН3СН2СОSКоА НООССНСН3СОSКоА 4. В основе – дефект метилкротонилКоАкарбоксилазы. ПропионилКоА образуется при расщеплении изолейцина метионина треонина жирных кислот с нечетным числом атомов углерода.
24068. Фолиевая кислота – витамин В9, Вс 32.5 KB
  Всасывание фолатов осуществляются с помощью специфического механизма активного транспорта требует затраты энергии и обеспечивает поступление фолиевой кислоты в кровоток против концентрационного градиента. Недостаток биотина нарушает образование активной формы витамина – тетрагидрофолиевой кислоты. Первая стадия образования коферментных форм – это восстановление фолиевой кислоты в тетрагидрофолиевую кислоту при участии дегидрофолатредуктазы. Наиболее важной функцией коферментных форм фолиевой кислоты является их участие в биосинтезе пуриновых...
24069. Витамин В12-кобаламин 40.5 KB
  Коферментная форма витамина В12дезоксиаденозилкобаламин необходима для функционирования метилмалонилКоАмутазы которая обеспечивает изомеризацию метилмалонилКоА в сукцинилКоА: С разветвленной цепью Жирные кислоты С нечетным числом атомов С Холестерин Изолейцин Метионин Треонин Нарушения обмена витамина В12. Это нарушение приводит к накоплению метилмалонилКоА. МетилмалонилКоА ингибирует пируваткарбоксилазу и это нарушает превращение пирувата в оксалоацетат и в результате тормозится глюконеогенез развивается гипогликемия...
24070. Аскорбиновая кислота (витамин С) 98 KB
  Аскорбиновая кислота являясь донором водорода участвует в окислительновосстановительных реакциях и превращается при этом в дегидроаскорбиновую кислоту: Аскорбиновая кислота участвует в следующих биохимических процессах: Гидроксилирование триптофана в 5гидрокситриптофан синтез серотонина. Аскорбиновая кислота метгемоглобин ДАК гемоглобин ДАК глутатион АК окисленный глутатион Аскорбиновая кислота восстанавливает метгемоглобин в гемоглобин сама окисляется в дегидроксиаскорбиновую кислоту. Дегидроксиаскарбиновая кислота...