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

Кафедра высшей математики

Отчет по лабораторной работе № 4

Тема: Решение дифференциальных уравнений в системе MathCAD.

Выполнил студент ИГЭУ Теплоэнергетического Факультета Кафедры тепловых электростанций

Группы I-2x

Атаманчук М. С.

Проверил: Крутов А.О.

Лабораторная работа № 4.

Тема: Решение дифференциальных уравнений

Цель: Научиться решать дифференциальные уравнения в системе MathCAD методами Рунге — Кутты (rkfixed, rkadapt), Булирша — Штера (Bulstoer) и Odesolve

Описание метода Рунге — Кутты Метод Рунге-Кутта четвертого порядка для решения уравнения первого порядка.

Идея Рунге-Кута состоит в том, чтобы использовать метод неопределённых коэффициентов. Наиболее употребительным методом Рунге-Кутта решения уравнения первого порядка y' = F (x, y) (2.1.1) является метод четвертого порядка, в котором вычисления производятся по формуле:

yk+1 = yk +(k1 +2k2 +2k3 +k4)/6, (2.1.2)

где

k1 = Fk h = F (xk, yk) h

k2 = F (xk +h/2, yk +k1 /2)h

k3 = F (xk +h/2, yk +k2 /2)h

k4 = F (xk +h, yk +k3)h,

k = 0, …, n-1

h = (xfx0)/n (2.1.3)

Описание метода Булирша — Штера (Bulstoer)

Метод Булирша-Штера — это метод решения системы обыкновенных дифференциальных уравнений первого порядка с гладкими правыми частями. Гладкость правых частей является необходимой для работы метода. Если правые части вашей системы не являются гладкими или содержат разрывы, то лучше использовать метод Рунге-Кутта. В случае же гладкой системы метод Булирша-Штера позволяет добиться существенно большей точности, чем метод Рунге-Кутта.

Принцип работы метода Основной идеей метода является вычисление состояния системы в точке x+h, как результата двух шагов длины h/2, четырех шагов длины h/4, восьми шагов длины h/8 и так далее с последующей экстраполяцией результатов. Метод строит рациональную интерполирующую функцию, которая в точке h/2 проходит через состояние системы после двух таких шагов, в точке h/4 проходит через состояние системы после четырех таких шагов, и т. д., а затем вычисляет значение этой функции в точке h = 0, проводя экстраполяцию.

Гладкость правых частей приводит к тому, что вычисленное при помощи экстраполяции состояние системы оказывается очень близко к действительному, а использование рациональной экстраполяции вместо полиномиальной позволяет ещё больше повысить точность.

Таким образом проводится один шаг метода, после чего принимается решение — следует ли изменять шаг, а если да — то в какую сторону. При этом используется оценка погрешности, которую мы получаем в качестве дополнительного результата при рациональной экстраполяции. Следует отметить, что алгоритм решает автономную систему, т. е. если уравнения системы содержат время, то необходимо ввести время в качестве переменной, производная от которой тождественно равна единице.

Задача № 1

Найдите два решения дифференциального уравнения на отрезке [1; 5] с начальными условиями соответственно,. Проконтролируйте достижение точности 0.001. Нарисуйте их графики (на одном чертеже). Определите значения этих решений в точке x = 4.85.

Метод Рунге — Кутты (rkfixed)

Создадим функцию Далее применим к ней функцию rkfixed

Для достижения точности 0.001 ниже полученных матриц Z, Z1 выполним оператор TOL :=0.001.

После этого снова найдем решение, но обозначим его иначе:

Возьмем два вектора и, являющиеся столбцами значений первого и второго решений. Найдем модуль разности этих векторов.

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

Построим графики этих функций.

Для нахождения значения этих решений в точке x = 4.85 применим к функции функцию rkfixed, но правое ограничение выстави в виде x = 4.85 и построим матрицу решений.

Решения дифференциального уравнения на отрезке [1; 5] с начальными условиями, в точке x = 4.85 будет y (4.85)=1.239 и y (4.85)=1.233 соответственно.

Метод Рунге — Кутты (rkadapt)

Создадим функцию

Далее применим к ней функцию rkadapt

Для достижения точности 0.001 ниже полученных матриц Z, Z1 выполним оператор TOL :=0.001.

После этого снова найдем решение, но обозначим его иначе:

Возьмем два вектора и, являющиеся столбцами значений первого и второго решений. Найдем модуль разности этих векторов.

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

Построим графики этих функций.

Для нахождения значения этих решений в точке x = 4.85 применим к функции функцию rkadapt, но правое ограничение выставим в виде x = 4.85 и построим матрицу решений.

Решения дифференциального уравнения на отрезке [1; 5] с начальными условиями, в точке x = 4.85 будет y (4.85)=-0.127 и y (4.85)=-0.142 соответственно.

Метод Булирша — Штера (Bulstoer)

Создадим функцию

Далее применим к ней функцию Bulstoer

Для достижения точности 0.001 ниже полученных матриц Z, Z1 выполним оператор TOL :=0.001.

После этого снова найдем решение, но обозначим его иначе:

Возьмем два вектора и, являющиеся столбцами значений первого и второго решений. Найдем модуль разности этих векторов.

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

Построим графики этих функций.

Для нахождения значения этих решений в точке x = 4.85 применим к функции функцию Bulstoer, но правое

ограничение выставим в виде x = 4.85 и построим матрицу решений.

Решения дифференциального уравнения на отрезке [1; 5] с начальными условиями, в точке x = 4.85 будет y (4.85)=-0.132 и y (4.85)=-0.146 соответственно.

Odesolve

Набираем в MathCad

Построим графики Решения дифференциального уравнения на отрезке [1; 5] с начальными условиями, в точке x = 4.85 будет y (4.85)=-0.129 и y (4.85)=-0.143 соответственно.

Общий график решения (Рунге — Кутты (rkfixed, rkadapt), Булирша

— Штера (Bulstoer) и Odesolve)

Таблица результатов с разными методами

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

Z — rkfixed, R — Rkadapt, B — Bulstoer

Задача № 2

Найдите решение дифференциального уравнения на отрезке [2; 4] с начальным условием. Нарисуйте график этого решения. Найдите значение этого решения в точке x = 3.75 с точностью .

Метод Рунге — Кутты (rkfixed)

Создадим функцию

Далее применим к ней функцию rkfixed.

Для достижения точности 10^-5 ниже полученной матрицы Z выполним оператор TOL := 10^-5.

После этого снова найдем решение, но обозначим его иначе:

Возьмем два вектора Z^<1> и Z1^<1>, являющиеся столбцами значений первого и второго решений. Найдем модуль разности этих векторов.

Видим что результаты совпадают в пределах заданной погрешности что и было необходимо для нашей задачи.

Построим графики этих функций.

Для нахождения значения этих решений в точке x = 3.75 применим к функции функцию rkfixed но правое ограничение выставим в виде x = 3.75 и построим матрицу решений.

Решением дифференциального уравнения на отрезке [2; 4] с начальными условиями, в точке x = 3.75 будет y (3.75)=0.45

Метод Рунге — Кутты (rkadapt)

Создадим функцию Далее применим к ней функцию rkadapt

Для достижения точности 10^-5 ниже полученной матрицы Z выполним оператор TOL := 10^-5.

После этого снова найдем решение, но обозначим его иначе:

Возьмем два вектора Z^<1> и Z1^<1>, являющиеся столбцами значений первого и второго решений. Найдем модуль разности этих векторов.

Видим что результаты совпадают в пределах заданной погрешности что и было необходимо для нашей задачи.

Построим графики этих функций.

Для нахождения значения этих решений в точке x = 3.75 применим к функции функцию rkfixed но правое ограничение выставим в виде x = 3.75 и построим матрицу решений.

Решением дифференциального уравнения на отрезке [2; 4] с начальными условиями, в точке x = 3.75 будет y (3.75)=0.45

Метод Булирша — Штера (Bulstoer)

Создадим функцию

Далее применим к ней функцию Bulstoer

Для достижения точности 10^-5 ниже полученной матрицы Z выполним оператор TOL := 10^-5.

После этого снова найдем решение, но обозначим его иначе Возьмем два вектора Z^<1> и Z1^<1>, являющиеся столбцами значений первого и второго решений. Найдем модуль разности этих векторов.

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

Построим графики этой функции.

Для нахождения значения этих решений в точке x = 3.75 применим к функции функцию Bulstoer но правое ограничение выставим в виде x = 3.75 и построим матрицу решений.

Решением дифференциального уравнения на отрезке [2; 4] с начальными условиями, в точке x = 3.75 будет y (3.75)=0.45

Odesolve

Набираем в MathCad

Построим график Решением дифференциального уравнения на отрезке [2; 4] с начальными условиями, в точке x = 3.75 будет y (3.75)=0.45

Общий график решения (Рунге — Кутты (rkfixed, rkadapt), Булирша — Штера (Bulstoer) и Odesolve)

Таблица результатов с разными методами

Результаты различных методов одинаковы но по общему графику видно что наиболее гладкая функция у метода Odesolve. Z — rkfixed, R — Rkadapt, B — Bulstoer

Задача № 3

Найдите решение системы дифференциальных уравнений на отрезке [0; 4] с начальными условиями,. Нарисуйте графики компонент и фазовый портрет решения.

Метод Рунге - Кутты rkfixed

Сформируем вектор начальных условий Сформируем вектор правых частей Найдем матрицу значений решения этой системы

Построим фазовый портрет решения графики компонент решения.

Фазовый портрет решения Графики компонент решения

Метод Рунге - Кутты Rkadapt

Сформируем вектор начальных условий сформируем вектор правых частей

Найдем матрицу значений решения этой системы

Построим фазовый портрет решения графики компонент решения

Фазовый портрет решения Графики компонент решения

Метод Булирша - Штёра Bulstoer

Сформируем вектор начальных условий Сформируем вектор правых частей Найдем матрицу значений решения этой системы Построим фазовый портрет решения графики компонент решения

Фазовый портрет решения Графики компонент решения

Odesolve

Введем в MathCAD

Построим фазовый портрет решения графики компонент решения Фазовый портрет решения Графики компонент решения

Общий график решения (Рунге — Кутты (rkfixed, rkadapt), Булирша — Штера (Bulstoer) и Odesolve)

Общий фазовый портрет решения

Общие графики компонент решения

Наиболее точным является метод Odesolve, так как из общих графиков видно что его функция более гладкая.

Z — rkfixed, R — Rkadapt, B — Bulstoer

Задание № 4

Найдите решение дифференциального уравнения с начальными условиями, , на отрезке [0; 4]. Нарисуйте графики

Метод Рунге - Кутты rkfixed

Сформируем вектор начальных условий Сформируем вектор правых частей Найдем матрицу значений решения этой системы

Построим график

Метод Рунге - Кутты Rkadapt

Сформируем вектор

начальных условий Сформируем вектор правых частей Найдем матрицу значений решения этой системы дифференциальный уравнение MathCAD

Построим график

Метод Булирша - Штёра Bulstoer

Сформируем вектор начальных условий Сформируем вектор

правых частей Найдем матрицу значений решения этой системы

Построим график

Odesolve

Построим график

Общий график

На графике видно что все методы совпадают по точности.Y — rkfixed, R — Rkadapt, B — Bulstoer

Вывод: мы научились решать дифференциальные уравнения в системе MathCAD методами Рунге —Кутты (rkfixed, rkadapt), Булирша — Штера (Bulstoer) и Odesolve. Наиболее точными являются методы Odesolve и Bulstoer, наихудшим rkfixed.