Лабораторная работа N2.
Моделирование динамических систем. Решение систем дифференциальных уравнений.
Цель работы:
Разработка алгоритмов и программ моделирования систем дифференциальных уравнений методами:
1. Эйлера-Коши;
2. Эйлера-Коши c итерациями;
3. модифицированным методом Эйлера;
4. методом трапеций;
5. методом Рунге-Кутта 4_го порядка;
6. методом Рунге-Кутта с автоматическим изменением шага;
7. методом Рунге-Кутта-Мерсона с автоматическим изменением шага;
8. методом Рунге-Кутта-Фельберга с автоматическим изменением шага;
Порядок выполнения работы:
I. Разработать программную реализацию для приведенных методов решения систем дифференциальных уравнений.
2. Получить решения для заданных систем дифференциальных уравнений в соответствии с индивидуальными заданиями.
3. Представить графические результаты вычислений.
Отчет о работе должен содержать:
1. Содержание
2. Описание методов решения систем дифференциальных уравнений;
3. Блок-схема алгоритма решения поставленной задачи.
4. Листинг программы (с комментариями).
5. Скриншоты программы с графиками процессов.
6. Выводы
Требования к оформлению:
Размер лабораторной работы - около 10 стр. Выводы должны занимать не больше 1 стр;
Стандартный шрифт Times New Roman (14 пт.) с межстрочным единичным интервалом и отступом «первой строки» — 1,25 см;
Методические указания:
Задача Коши заключается в решений систем обыкновенных дифференциальных уравнений первого порядка, представляемых в виде
dy1/dx = F1 (x,y1,...,yJ,...,yN),
dyj/dx= FJ(x,y1,...,yJ,...,yN), (1)
dyn/dx= FN(x,y1,...,yJ,...,yN),
где j=1-N- номер каждой зависимой переменной уj,х- независимая переменная. Если задача Коши решается для анализа поведения системы или объекта во времени, то х является временем (x=t). Решение системы (I) при заданных начальных условиях х=х0, у1(x0)=y10, y2(x0)=y20,…,yN(x0)=yN0 сводится к наложению зависимостей (интегральных кривых) y1(x),y2(x),…,yN(х), проходящих через точки, заданные начальными условиями (x0,y10),(х0,y20),...,(x0,yj0),…,(x0,yN0). Задача Коши сводится к интегрирований дифференциальных уравнений. Порядок метода численного интегрирования при этом определяет и порядок метода решения (I).
Обобщенная форма записи каждого из уравнений системы (1) может быть представлена в общем виде
dyj /dx=Fj(x,y), (2)
где yj где в правой части уравнения - векторы переменных у1,y2,yj,…,yN, а Fj- правая часть каждого из уравнений (I). В частности, одно дифференциальное уравнение (y=Yj=Y1 и Fj=F=F1) записывается в виде
dy/dx=F1(x,y).
Дифференциальные уравнения высшего порядка
y(n)=F(x,y,y’,y’’,…y(n-1)), 3)
где (n) -порядок уравнения, могут быть сведены к системам вида (I) или (2) с помощью следующих преобразований:
dy/dx=y1,
dy1/dx=y2, (4)
dyn-2/dx=yn-1,
dyn-1/dx=F(x,y1,y2,…,yn-1),
Следовательно решение (з) сводится к решению системы обыкновенных дифференциальных уравнений первого порядка (4).
1. Метод Эйлера-Коши- простейший метод первого порядка для численного интегрирования дифференциальных уравнений. Он реализуется следующей рекурентной формулой:
yj(i+1)=yji+hFj(xi,yji).
где h- шаг интегрирования (приращение переменной x). Этот метод обладает большой погрешностью и имеет систематическое накопление ошибок. Погрешность метода R ~ (h2), т.е. пропорциональна h2.
2. Метод Эйлера-Коши с итерациями заключается в вычислении на каждом шаге начального значения
= Yji + hFj (xi,Yji).
Затем с помощью итерационной формулы
=yji+h/2 [ Fj(xi,yji)+Fj(xi+1, ) ]
решение уточняется. Итерации проводят до тех пор, пока не совпадет заданное число цифр результата на двух последних шагах итераций. Погрешность метода R~ (h2). Обычно число итераций не должно превышать 3-4, иначе нужно уменьшить шаг h.
3. Модифицированный метод Эйлера второго порядка реализуется следующими рекурентными формулами:
Yj(i+1)=yji+hFj(xi+h/2,y*j(i+1/2)),
где y*j(i+1/2)=yji+hFj(xi,yji)/2
Метод дает погрешность R ~ (h2) и имеет меньшее время вычислений, поскольку вместо нескольких итераций производится вычисление только одного значения.
4. Метод трапеций - одна из модификаций метода Эйлера второго порядка. Он реализуется применением на каждом шаге формулы
yj(i+1)=yji+1/2(Kj1+Kj2),
где Kj1 =hFj (xi,yji), Kj2=hFJ (xi+h, yji+Kj1), и дает погрешность
R ~ (h2). Этот метод относится к общим методам Рунге-Кутта. Не рекомендуется при сильной крутизне, т.к. возможно возникновение численной неустойчивости (уменьшить шаг).
5. Метод Рунге-Кутта четвертого порядка является наиболее
распространенным методом решения систем (2) при шаге h=const. Его достоинством является высокая точность - погрешность R ~ (h5) - и меньшая склонность к возникновению неустойчивости решения. Алгоритм реализации метода Рунге-Кутта заключается в циклических вычислениях yj(i+1) на каждом шаге i+1 шаге по следующим формулам:
K1j=hFj(xi,yji),
K2j=hFj(xi+h/2, yji+ K1j)
K3j=hFj(xi+h/2, yji+ K2j)
K4j=hFj(xi+h, yji+K3j)
Yj(i+1)=yji+1/6(K1j+2K2j+2K3j+K4j)
При переходе от одной формулы к другой задаются или вычисляются соответствующие значения х и Yj. и находятся по подпрограмме значения функций Fj(x,yj).
Решение одного дифференциального уравнения методом Рунге-Кутта производится по приведенным формулам, если в них опустить индекс j, а из алгоритма исключить циклы. Последнее резко упрощает программу и позволяет получить минимально возможное время счета.
Ввиду особого значения и широкого применения дифференциальных уравнений второго порядка полезно иметь специальную программу для их решения.
Метод Рунге-Кутта для дифференциального уравнения второго порядка вида
Y''= =F(x,y,y')
имеющий погрешность R ~ (h5), реализуется с помощью следующих формул:
К1 = hF(xi;yi;y'i);
K2=hF(xi+h/2; yi+ y'i+ K1; y'i+K1/2);
K3=hF(xi+h/2; yi+ y'i+ K1; y'i+K2/2);
K4=hF(xi+h; yi+hy'i+ K3; y'i+K3);
Y'i+1=y'i+1/6(K1+2K2+2K3+K4)
Перед началом вычислений надо задать шаг h и начальные значения х x0, y(x0)=y0 и y'(x0)=y'0.
Автоматическое изменение шага в ходе решения систем
дифференциальных уравнений необходимо, если решение требуется получить с заданной точностью. При высокой точности (погрешность e=Е=1*10-9) и решении в виде кривых с сильно различающейся крутизной автоматическое изменение шага обеспечивает уменьшение общего числа шагов в несколько раз, резко уменьшает вероятность возникновения числовой неустойчивости, дает более равномерное расположение точек графика кривых (решений) при их выводе на печать.
6. Метод Рунге-Кутта с автоматическим изменением шага заключается в том, что после вычисления Yj(i+1) c шагом h все вычисления проводятся повторно с шагом h/2. Полученный результат сравнивается с yj(i+1). Если |yj(i+1) - y*j(i+1)|<e, вычисления продолжается с шагом h, в противном случае шаг уменьшают. Если это неравенство слишком сильное, шаг, напротив, увеличивают. При той же погрешности R ~ (h5) лучшие результаты дает описанный ниже метод.