Математика/1. Дифференциальные и интегральные уравнения

 

К.т.н. Фомина Е.Е.

Тверской государственный технический университет, Россия

Особенности решения дифференциальных уравнений разных типов в MatLab

 

Введение

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

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

Решение задачи Коши для дифференциального уравнения

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

  , (1)

где Y, Y0 – скаляры или векторы-столбцы, F – функции правых частей уравнений.

Библиотека MatLab содержит следующие основные функции, позволяющие решать задачу (1):

1) нежесткие задачи:

  а) с небольшой точностью - ode45 (использует явный метод Рунге-Кутты 4-го и 5-го порядков); ode23 (использует явный метод Рунге-Кутты 2-го и 3-го порядков);

  б) с большой точностью - ode113 (использует метод Адамса типа предиктор-корректор переменного порядка);

2) жесткие задачи:

  а) с небольшой точностью - ode15s (использует численные формулы «дифференцирования назад»); ode23s (использует модифицированный метод Розенброка); ode23tb (использует неявный метод Рунге-Кутты с первым шагом по методу трапеций и вторым шагом по методу «дифференцирования назад»
2-го порядка);

  б) с большой точностью ode23t (использует неявный метод трапеций).

Наиболее частое обращение к функции имеет вид:

[T,Y]=ode*(odefun,tspan,y0),

где odefun – указатель на функцию вычисления правых частей уравнения (1); tspan – отрезок интегрирования дифференциального уравнения или системы; y0 – начальное значение зависимой переменно (скаляр или вектор); T – выходной вектор-столбец значений независимой переменной; Y – выходной вектор-столбец зависимой переменной, в котором каждая строка соответствует одному значению в столбце T.

Пример: найти решение задачи Коши , , где S, S0, h0 – произвольные константы большие нуля, g=9,8.

Решение:

1. Создается функция vt для расчета правой части уравнения (рис. 1).

2. Оформляется форма, предназначенная для ввода исходных данных и вывода результата (рис. 2).

3. Для кнопки Расчет записывается текст программы, позволяющий производить решение уравнения и выводить результаты в виде графика (рис. 3).

 

Рис. 1. Функция vt                                    Рис. 2. Форма

 

Рис. 3. Текст программы

 

Решение краевой задачи для дифференциального уравнения

Система позволяет решать краевую задачу, сводящуюся к следующему виду:

  , (2)

где Y, Y0 – скаляры или векторы-столбцы, F – функции правых частей уравнений.

Граничные условия могут иметь вид: , , , , ,  и др.

Алгоритм решения краевой задачи следующий:

1. Задается начальное приближение (xi,yi) с помощью функции

solinit=bvpinit(x,yinit),

где x – вектор-строка a=x1<x2<…<xn=b; yinit – гипотетическое значение для yi; solinit – структура, в которой заполняются два поля solinit.x=xi и solinit.y=yi.

2. Краевая задача решается, используя функцию bvp4c, которая имеет следующий синтаксис:

sol = bvp4c(odefun,bcfun,solinit),

где odefun – указатель на функцию, вычисляющую вектор правых частей уравнения (2); bcfun – указатель на функцию, вычисляющую вектор граничных условий; sol – выходная структура функции bvpinit, содержащая решения краевой задачи и включающая поля sol.x, sol.y, sol.yp.

Пример: решить краевую задачу , y(0)=1, y(1)=1, где p=const>0

Решение:

1. Уравнение путем замены сводится к системе уравнений следующего вида: , .

2. Создается функция sistema для расчета правых частей уравнений системы (рис. 4).

3. Создается функция border1 для расчета ограничений (рис. 5).

  

Рис. 4. Функция sistema                             Рис. 5. Функция border1

4. Задаются начальные значения переменной x и начальное приближение вектора y в этих точках. Формируется сетка (рис. 6).

5. Осуществляется поиск решения (рис. 7).

 

         Рис. 6. Формирование сетки                         Рис. 7. Поиск решения

 

Решение дифференциального уравнения неявного типа

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

[T,Y]=ode15i(odefun,tspan,y0,yp0),

где odefun – указатель на функцию вычисления правых частей уравнения (3); tspan – отрезок интегрирования дифференциального уравнения или системы; y0 – начальное значение зависимой переменно (скаляр или вектор); yp0 – начальное значение производной, такое, что F(t0,y0,yp0)=0.

Для подбора значение y0 и yp0 служит функция decic, которая имеет следующий синтаксис:

[y0mod,x0mod]=decic(odefun,t0,y0,fixed_y0,yp0,fixed_yp0),

где y0, yp0 – задаваемые пользователем гипотетические начальные значения, которые не обязательно должны удовлетворять уравнению (3); y0mod, yp0mod – вычисленные начальные значения, которые удовлетворяют уравнению (3).

Параметры fixed_y0, fixed_yp0 позволяют управлять ходом решения. Если fixed_y0=1, то y0mod=y0, если fixed_y0=0, то y0mod может быть отличным от y0 (аналогично для yp0mod).

Пример: найти решение уравнения .

Решение:

1. Создается функция ne_y для расчета левой части уравнения (рис. 8).

Рис. 8. Функция border1

2. Осуществляется подбор значений y0 и yp0, удовлетворяющих уравнению (3) (рис. 9).

3. Осуществляется поиск решения (рис. 9, 10).

 

            Рис. 9. Формирование сетки                   Рис. 10. График функции

 

Литература

1. Алексеев, Е.Р. MATLAB 7 [Текст] / Е.Р. Алексеев, О.В. Чеснокова. М.: НТ Пресс, 2006. 464 с.

2. Кетков, Ю.Л. MATLAB 7: программирование, численные методы [Текст] / Ю.Л. Кетков, А.Ю. Кетков, М.М. Шульц. СПб.: БХВ-Петербург, 2005. 752 с.

3. MATLAB Version 7.5.0.342 (R2007b). August 15, 2007. License Number 161052.