Решение обыкновенных дифференциальных уравнений первого порядка

Рассмотрим обыкновенное дифференциальное уравнение (ОДУ) n-го порядка

, (1)

где х ­­ ­­­ ­­­– независимая переменная;

– неизвестная функция;

­­­­ – ее производные.

Функция рассматривается, как функция от n+2 независимых переменных и всегда будет предполагаться непрерывной.

Порядком уравнения называется максимальный порядок входящей в него производной (или дифференциала)­­.

Запишем ОДУ (1) в нормальной форме

(2)

На практике часто приходится решать задачу Коши или краевую задачу.

Сформулируем задачу Коши для ОДУ первого порядка.

Требуется найти такое решение уравнения

,

которое удовлетворяет начальному условию .

В пакете MathCAD можно найти аналитическое решение этого уравнения, используя команду Переменная ­– Интегрировать из меню Символика. А затем найти значение константы С для заданных начальных условий и протабулировать искомую функцию, если это необходимо.

Рассмотрим пример.

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

Но не всегда удается получить аналитическое решение ОДУ и тогда используют численные методы. В пакете Mathcad имеется 13 встроенных функций для решения дифференциальных уравнений и систем различных порядков, использующих численные методы[5,6].

Почти всегда справляется с поставленной задачей функция rkfixed. В нее заложен широко распространенный метод решения дифференциальных уравнений – метод Рунге-Кутта четвертого порядка (сам метод будет рассмотрен позже). Однако есть случаи, когда лучше использовать более сложные методы. Нередко приходится пробовать на одном дифференциальном уравнении несколько методов, чтобы определить, какой метод точнее, быстрее.

Рассмотрим функцию rkfixed и ее использование для ОДУ первого и более высокого порядка. Эта функция ищет решение на отрезке в равноотстоящих точках.

Функция rkfixed имеет следующие аргументы:

– вектор начальных условий размерности , где – порядок дифференциального уравнения (ДУ) или число уравнений в системе ( если решается система уравнений). Для ДУ первого порядка вектор начальных условий вырождается в одну точку ;

– граничные точки интервала, на котором ищется решение ДУ. Начальные условия, заданные в векторе , – это значение решения в точке ;

– число точек(не считая начальной), в которых ищется приближенное решение. Этот аргумент определяет число строк в матрице, возвращаемой функцией rkfixed;

– функция, возвращающая значение в виде вектора из элементов, содержащих первые производные неизвестных функций.

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

При использовании функции rkfixed наибольшие трудности вызывает определение функции .

Для дифференциального уравнения первого порядка достаточно разрешить это уравнение относительно первой производной . В случае нелинейных уравнений это может вызывать затруднения. Тогда можно воспользоваться возможностями символьной математики пакета Mathcad. Для этого используем команду Переменная–Решить из меню Символика.

Рассмотрим пример решения ДУ первого порядка с использованием функции rkfixed.

Решить задачу Коши . Будем искать решение на интервале .

Сначала найдем аналитическое решение этого уравнения:

.

Мы вывели в таблице первые 15 значений решения. Если необходимо, то таблицу можно увеличить, растянув вниз поле для таблицы.

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

Функция rkfixed почти всегда справляется с поставленной задачей, т.к. в нее заложен широко распространенный метод решения дифференциальных уравнений – метод Рунге-Кутта. Однако есть случаи, когда лучше использовать более сложные методы. Нередко на одном дифференциальном уравнении приходится пробовать несколько методов, чтобы определить, какой метод точнее, быстрее.

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

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

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

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

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

Список аргументов и матрица, возвращаемая функцией Bulstoer, те же, что и для функции rkfixed : .

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

Рассмотрим пример решения ДУ с помощью функции rkfixed и функции Bulstoer и сравним результаты.

Решить задачу Коши . Будем искать решение на интервале .

Сначала найдем аналитическое решение этого уравнения:

.

Для каждой функции в таблицах 1 и 2 вывели первые 15 значений решения. Чтобы показать, что функция Bulstoer дает более точные результаты, чем функция rkfixed мы вывели в таблицу 3 разность между полученными решениями для последних значений аргумента.

Можно решить задачу более точно (более быстро), если уменьшать шаг там, где производная меняется быстро, и увеличивать шаг там, где она ведет себя более спокойно. Для этого можно использовать функцию Rkadapt, основанную на методе Рунге-Кутты. Этот адаптивный контроль величины шага дает возможность повысить точность и сократит время, требуемое для решения уравнения. Эта функция обеспечивает также решение для точек, находящихся на одинаковом расстоянии (интервале), заданном пользователем.

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

Список аргументов и матрица, возвращаемая функцией Rkadapt, те же, что и для функции rkfixed :

.

Используем для вышеприведенного примера функцию Rkadapt.

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

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

Если надо найти приближенное решение задачи Коши только в конечной точке отрезка , то тогда лучше использовать следующие функции: rkadapt и bulstoer.

Функции rkadapt и bulstoer имеют следующие аргументы:

­­– вектор начальных условий размерности , где – порядок дифференциального уравнения (ДУ) или число уравнений в системе ( если решается система уравнений);

– граничные точки интервала, на котором ищется решение ДУ. Начальные условия, заданные в векторе , – это значение решения в точке ;

– параметр, контролирующий точность решения. Малое значение заставляет алгоритм брать меньшие шаги вдоль траектории, что увеличивает точность решения. Значения, близкие к 0.001, дают хорошую точность решения;

– функция, возвращающая значение в виде вектора из элементов, содержащих первые производные неизвестных функций;

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

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

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

Количество строк в матрице решения определяется параметром . Чем больше это значение, тем меньше строк будет в матрице решения.

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

Таким образом, .

Задание

1. Найти аналитическое решение задачи Коши вручную и с помощью пакета MathCad - с использованием символьной математики.

2. Решить задачу Коши, используя функции rkfixed, Rkadapt и Bulstoer.

3. Построить интегральную кривую, вывести таблицы полученных результатов и сравнить их (вывести таблицу разностей).

4. Найти значение функции в конечной точке отрезка, на котором ищется решение, используя функции rkadapt или bulstoer

Отчет оформить как файл пакета MathCad.

Варианты работы № 1

Работа № 2