Линейное приближение точки пересечения кривой решения системы обыкновенных дифференциальных уравнений с поверхностью разрыва правой части Текст научной статьи по специальности «Математика»

УДК 517.91

ЛИНЕЙНОЕ ПРИБЛИЖЕНИЕ ТОЧКИ ПЕРЕСЕЧЕНИЯ КРИВОЙ РЕШЕНИЯ СИСТЕМЫ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ С ПОВЕРХНОСТЬЮ РАЗРЫВА

ПРАВОЙ ЧАСТИ

Р.Г. Идрисов, В.В. Коробицын

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

Введение

Многие математические модели, описывающие процессы в механике, электротехнике, биологии и других областях, задаются системами обыкновенных дифференциальных уравнений (ОДУ) с разрывной правой частью [6], [12]. Траектория решения таких систем теряет гладкость на поверхности разрыва, что может приводить к росту погрешности численного решения при пересечении границы непрерывности.

Традиционно для нахождения численного решения обыкновенных дифференциальных уравнений используют методы Рун го Купы или Адамса, предполагая выполнение условия Липшица для функции правой части [1], [7], [8]. Однако для систем с разрывами применение этих методов приводит к неуправляемому росту погрешности. Даже если метод снабжен процедурой изменения шага интегрирования, все равно погрешность может увеличиваться в точке разрыва, и при этом может происходить рост числа итераций численного метода [2].

В настоящее время много работ посвящено численному решению систем ОДУ с разрывом в правой части [5], [9], [10], [11]. Однако, проведя обзор этих

Copyright (с

Омский государственный университет им. Ф.М. Достоевского. E-mail: idrisov@omskportal.ru

работ, нельзя сказать, что существуют универсальные алгоритмы решения систем ОДУ с разрывом в правой части,

В данной статье предлагается подход для дополнения традиционных методов Рун го Купы процедурой нахождения точки пересечения траектории решения с поверхностью разрыва. Дополнительные аспекты решения систем ОДУ с разрывами освещены в работах [3], [4].

1. Описание подхода вычисления точки разрыва

Рассмотрим систему обыкновенных дифференциальных уравнений вида,

Зх

— = /(Ь,х), £ > 0, ж(0) = х0,

или в покомпонентной записи

^ = /*(*,жь ...,хп), г > О, аЬ

Хг(0) = х0, г = 1, 2,п,

где х0 — координаты начальной точки траектории решения х^(Ь). Некоторые функции ^(Ь,х 1,...,хп) терпят разрыв па поверхностях, задаваемых уравнениями фу (Ь, х1, ...,хп) = 0, ] = 1, 2, ...к.

Существует два подхода для численного решения систем ОДУ с разрывом в правой части: 1) игнорировать разрыв; 2) на поверхности разрыва остановить вычисление и продолжить с новыми значениями функций,

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

Для реализации второго подхода предлагаем следующий алгоритм:

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

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

фу (Ь,х1 ,...,хп) 0,

ф(Ь,х1 ,...,хп) = 0,

где ф(Ь,х1, ...,хп) = 0 - уравнение прямой, составленной из точки решения до пересечения с поверхностью разрыва и точки решения после пересечения поверхности разрыва;

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

Для построения функции ф(Ь,х1, ...,хп) необходимо взять точку решения М1 (Ь(1), х1\ ..., х^), полученную на последнем до пересечения с поверхностью

разрыва шаге, и точку М2(г(2),х12),...,хП)) — первую точку решения, полученную после пересечения с поверхностью разрыва. Тогда уравнение прямой,

г — г(1)

проходящей через эти две точки, примет вил

Ь(2) — г«

(1)

/у» _ /у» 4 '

хп хП

~(2) (ТУ

/у V ' _ /у V '

ХП ХП

(1)

-•у* _ -г*4 '

х1 х1

(2) (Ту

,"у»4 ' _ ,"у»4 '

11

Получаем

(1)

пг>__ -г*

11

(2)

(1)

(1)

/у» __ /у» 4 '

хп хП

~(2) оУ

^у» 4 ' __ ^у»4 '

V Хп Хп

Решение этой системы можно выразить в явном виде

хДЬ) = Х(1) +

(2)

/у»4 ' _ /у»

хх

(1)

(г — Ь(1)).

г(2) — г«

Тогда решение системы (1) преобразуется к виду

Фз(ь,Х1 (г),... ,Хп(г)) = о,

где все Хг(г) выражены явно через г. Задача вычисления точки пересечения кривой решения с поверхностью разрыва сводится к нахождению значения свободной переменной г, значение которой соответствует искомой точке. Однако в нашем случае эта точка будет показывать точку пересечения отрезка прямой линии, концы которого лежат в точках до и после пересечения поверхности разрыва.

2. Вычислительный эксперимент

Для тестирования описанного подхода для численного решения систем ОДУ с разрывом в правой части было исследовано поведение численных методов на решении следующей системы ОДУ:

“гг = У~& ь

Х — с1,

~тг =

Х — С2,

при х < а,

(2)

при х > а,

с заданными параметрами й1 = ^2 = 0.5, с1 = 0.2, с2 = 0.8, а = 0.5 и начальными данными г0 = 0,х0 = 0.4999, у0 = 0.3,

Дня этой системы известно точное решение:

где

х(г) = Л1 ■ в- + Л2 ■ в-(1-Ьо) + ах,

у (г) = Л1 ■ в1-10 — Л2 ■ в-(1-10) + ау,

(х0 — ах) + (Уо — ау) . (х0 — ах) — (Уо — ау)

= ------------^-----------> 2 = ---------------------------1

ах

ау

с1,

d1,

с2,

d2,

х0 < а, х0 > а.

у2

х=а

периодическим с периодом

Т = 21п

- х01 + л/(ах - Хо)2 - 4А1А2 2Ал

Траектория точного решения на фазовой плоскости изображена на рисунке 1,

Рис. 1. Траектория решения системы (2) на фазовой плоскости

а

х

Дня проведения эксперимента была написана программа па языке программирования С---, реализующая следующие методы решения систем ОДУ: мето-

ды Рупге—Кутты 3 и 4 порядка, метод Руиге—Кутты—Фельберга 4(5) порядка,

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

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

Егг(гк) = у/(хк - х*(гк))2 + (ук - у*{ь))2,

где хк ,ук — координаты точки численного решения, соответствующего моменту времени гк, Х*(гк), У*(гк) — координаты точки точного решения в момент време-гк

Рунге—Кутты 3 и 4 порядков и соответствующих им модификаций представлены на рисунке 2, Как видно из рисунка, локальная погрешность стандартных методов растет с каждым периодом, в то время как погрешность модифицированного метода остается в пределах заданной точности.

б) модифицированный метод Р-КЗ

г) модифицированный метод Р-К4

Рис. 2. Графики локальных погрешностей методов Рунге Кутты 3 и 4 порядков до модификации и после (заданная точность 10-4)

Графики изменения локальной погрешности дня методов, соответствующих схеме Рунге—Кутты—Фельберга 4 порядка с модификацией, представлены на рисунке 3, Из рисунка видно, что локальная погрешность метода до модификации растет, а после модификации метода остается в пределах заданной точности.

Графики изменения локальной погрешности метода Розенброка (неявный метод Рунге—Кутты) представлены на рисунке 4, Этот метод показан себя луч-

а) стандартный метод Р-К-Ф4 (заданная точность 10-4)

в) стандартный метод Р-К-Ф4 (заданная точность 10-8)

б) модифицированный метод

10-4

г) модифицированный метод

10-8

Рис. 3. Графики локальных погрешностей методов Рунге Кутты Фельберга до

модификации и после иее

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

И

V

а) стандартный метод Розенброка Рис. 4. Графики локальных погрешностей метода Розенброка до модификации и после нее

б) модифицированный метод Розенброка

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

через несколько периодов.

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

Литература

1. Деккер К., Вервер Я. Устойчивость методов Рунге-Кутты для жестких нелинейных дифференциальных уравнений. М.: Мир, 1988. 334 с.

2. Коробицын В.В., Маренич В.Б., Фролова Ю.В. Исследование поведения явных методов Рунге-Кутты при решении систем обыкновенных дифференциальных уравнений с разрывной правой частью // Математические структуры и моделирование. 2007. Вып. 17. С. 19-25.

3. Коробицын В.В., Фролова Ю.В. Алгоритм численного решения дифференциальных уравнений с разрывной правой частью // Математические структуры и моделирование. 2005. Вып. 15. С. 46-54.

4. Коробицын В.В., Фролова Ю.В., Маренич В.Б. Алгоритм численного решения кусочно-сшитых систем // Вычислительные технологии. 2008. Т. 13, N. 2. С. 70-81.

5. Новиков Е.А., Шорников Ю.В. Численное моделирование гибридных систем методом Рунге—Кутты второго порядка точности // Вычислительные технологии. 2008. Т. 13, N. 2. С. 99-105.

6. Филиппов А. Ф. Дифференциальные уравнения с разрывной правой частью. М.: Наука, 1985. 225 с.

7. Хайрер Э., Ваннер Г. Решение обыкновенных дифференциальных уравнений. Жесткие и дифференциально-алгебраические задачи. М.: Мир, 1999. 685 с.

8. Хайрер Э., Нерсетт С., Ваннер Г. Решение обыкновенных дифференциальных уравнений: Нежесткие задачи. М.: Мир, 1990. 512 с.

9. Gear C.W., Osterbv О. Solving ordinary differential equations with discontinuities // ACM Trans. Math. Software. 1984. Vol. 10. P. 23-44.

10. Piiroinen P.T., Kuznetsov Yu.A. An event-driven method to simulate Filippov systems with accurate computing of sliding motions // ACM Trans. Math. Software. 2008. Vol. 34, N. 13. P. 1-24.

11. Shampine L.F., Thompson S. Event location for ordinary differential equations // Computer and Mathematics with Application. 2000. Vol. 39. P. 43-54.

12. Zhusubalivev Z., Mosekilde E. Bifurcations and chaos in piecewise-smooth dynamical systems. Singapore: World Scientific, 2003.