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

Тип работы:
Реферат
Предмет:
Физико-математические науки


Узнать стоимость

Детальная информация о работе

Выдержка из работы

ВЕСТНИК ТОМСКОГО ГОСУДАРСТВЕННОГО УНИВЕРСИТЕТА
2007 Управление, вычислительная техника и информатика № 1
УДК 519. 865
В. В. Поддубный, О.В. Романович
МОДИФИКАЦИЯ МЕТОДА ЭЙЛЕРА С УРАВНИВАНИЕМ ДЛЯ РЕШЕНИЯ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ С ПЕРЕМЕННЫМ ЗАПАЗДЫВАНИЕМ
Рассматривается задача численного интегрирования обыкновенных нелинейных дифференциальных уравнений, разрешенных относительно производной, с переменным, в том числе случайным запаздыванием аргумента. Получена интерполяционная модификация вычислительной схемы Эйлера с уравниванием, пригодная для решения дифференциальных уравнений с произвольными, в том числе случайно изменяющимися и как угодно малыми запаздываниями, когда становится неприменимым известный метод шагов. Показано, что сохраняется второй порядок точности этой модификации как для итерационной структуры алгоритма метода Эйлера с уравниванием, так и для ее первого приближения (метода Хьюна).
Проблема численного интегрирования обыкновенных дифференциальных уравнений с запаздывающим аргументом при постоянных или переменных, но строго положительных конечных запаздываниях решается сравнительно просто с использованием алгоритмов численного интегрирования обыкновенных дифференциальных уравнений без запаздывания на основе известного метода шагов [1]. Однако эта проблема резко усложняется в случае, если переменные запаздывания (в том числе случайные) могут принимать как угодно малые значения. В этом случае метод шагов неприменим. А поскольку любые методы численного интегрирование дифференциальных уравнений оперируют с конечным, отличным от нуля шагом интегрирования, не обязательно соизмеримым с произвольно изменяющимися запаздываниями, в процессе численного интегрирования приходится прибегать к интерполяции решения в точках, где значения решения не были вычислены.
В [4] предложен алгоритм такой интерполяции для случая произвольного, в том числе как угодно близкого к нулю, постоянного запаздывания, модифицирующий известную схему метода Эйлера с уравниванием второго порядка точности. В данной работе этот подход развивается и распространяется на случай переменных запаздываний, в том числе случайных и как угодно малых. При этом проводится оценка порядка точности модифицированного для случая переменного запаздывания итерационного метода Эйлера с уравниванием. Показывается, что модифицированный итерационный метод Эйлера с уравниванием, как и его первое приближение (метод Хьюна), сохраняют второй порядок точности.
1. Постановка задачи
Рассмотрим начальную задачу для обыкновенного (в общем случае нелинейного) дифференциального уравнения, разрешенного относительно производной, с переменным запаздыванием т (г):
= f ((У ()& gt-У (-т ()) *о ^ ^ T, (1)
dt
с заданными начальными условиями
У 0) = Ф О), t е [to — max т (t), t0) — y (t0) = y0,
(2)
где ф (г) — начальная функция, у0 — начальное значение решения, причем ф ((0) не обязательно равно у0.
Пусть функция ф (г) непрерывна, функция / (г, у, z) непрерывна по всем аргументам и удовлетворяет условиям Липшица по второму аргументу. Тогда, как известно [1], решение начальной задачи для дифференциального уравнения (1) существует и единственно. В случае возможного скачка начального условия (ф ((0) Ф у0), а также в случае, если начальная функция ф (г) не удовлетворяет дифференциальному уравнению (1), решение у (г) при г & gt- г0 испытывает скачки производной первого рода. Непрерывное и всюду дифференцируемое решение возможно лишь в случае, если ф (г) удовлетворяет условию
При постоянном, заранее известном т для интегрирования дифференциального уравнения с запаздывающем аргументом (1) при начальных условиях (2) можно использовать известный метод шагов [1]. Однако в случае переменного т = т (7), тем более случайно изменяющегося во времени, метод шагов не применим. В этом случае задача интегрирования нелинейного дифференциального уравнения
(1) существенно усложняется. Существует множество вычислительных схем решения начальной задачи для нелинейного дифференциального уравнения с переменным запаздыванием [2]. Однако по большей части эти методы либо носят частный характер, либо весьма сложны в реализации.
В данной работе предлагается метод численного решения начальной задачи для дифференциального уравнения (1), обобщающий метод Эйлера с уравниванием при отсутствии запаздывания [3] и при постоянном запаздывании [4] на случай переменного (в том числе случайного и как угодно малого) запаздывания т (/).
Для дальнейшего анализа необходимо конкретизировать свойства функции х (г).
Пусть х (?) случайно изменяется во времени. В простейшем случае такое случайное изменение можно описать кусочно-постоянной случайной функцией, порождаемой случайным временным рядом тк, к = 1,2,…, в котором все тк неотрицательны, статистически независимы и имеют одинаковую плотность распределения р (тк). Например, р (тк) = 1/ттах — равномерное в интервале [0, ттах & gt- 0] распределение. Тогда при любом разбиении г0 & lt-г1 & lt-••• & lt-ги = Т интервала интегрирования [, Т] на п примыкающих друг к другу интервалов гк, +1), к = 1, п -1, имеем т (г) в виде случайного «меандра»:
Поскольку любое запаздывание тк с некоторой отличной от нуля вероятностью может оказаться меньше любого заданного s & gt- 0, хорошо известный для решения дифференциального уравнения с запаздыванием метод шагов [1] оказывается, вообще говоря, неприменимым (для его применимости необходимо, чтобы min т (t) & gt- 0). В связи с этим для интегрирования дифференциального уравнения
dф (t0 — 0)/dt = f (t0, ф (t0), ф (t0 -т (t0))).
(3)
k=1
со случайно изменяющимся запаздыванием приходится применять вычислительные методы, использующие интерполяцию значений решения в точках, где не были вычислены значения решения при применении той или иной вычислительной схемы интегрирования.
Методы первого порядка точности (метод ломаных Эйлера и подобные ему) не используются при интегрировании дифференциальных уравнений, в том числе дифференциальных уравнений с запаздывающим аргументом, в силу его низкой точности (медленной сходимости к точному решению с уменьшением шага интегрирования). Поэтому рассмотрим схему второго порядка точности — метод Эйлера с уравниванием [3, 5], называемый также итерационным методом Эйлера [5], методом Эйлера с пересчетом [6], симметричной схемой Эйлера [7], и его упрощенный (неитерационный) вариант — метод Хьюна [8] (он же — усовершенствованный метод Эйлера — Коши [9]). Модифицируем эту схему для решения дифференциального уравнения с произвольным, в том числе и как угодно малым запаздыванием, введя в нее возможность интерполяции.
2. Обобщение метода Эйлера с уравниванием и метода Хьюна на случай переменного запаздывания
Рассмотрим дифференциальное уравнение (1) при произвольных изменениях запаздывания во времени: т = т (г). Разобьем интервал [г0,Т] на п примыкающих друг к другу подынтервалов равной длины к = (Т -г0)/п разбиением: г0 & lt- г1 & lt- - & lt- гк & lt- - & lt- гп = Т, где гк = г0 + кк, к = 1, п. Проинтегрируем формально дифференциальное уравнение (1) на интервале [г0, г ] с учетом начальных условий (2):
г
у () = Уо + | / (у (г), у (-т (г))).
Используя это выражение при г = гк и I = 1к+1, получаем рекуррентное соотношение
*к +1
У (гг+1) = У (гк)+ I /(у (г), у (-т (*))). (4)

Считая функцию f (г, у (г), у (г — т (г))) непрерывной функцией своих аргументов, представим приближенно при малом шаге интегрирования к значение интеграла в правой части точного соотношения (4) по формуле трапеции:
У (**+1) = У (**)+ -2(f (& gt-У (**). У -т (к))) +
+/ (гк+1& gt- У (гк+1)& gt- У (гк+1 — т (гк+1)))). (5)
Обозначая ук = у (?к), тк = т (гк), представим рекуррентное уравнение (5), определяющее приближенное решение уравнения (1) при начальных условиях (2), в виде
к —
Ук+1 = Ук + -(/ (к & gt- У к & gt- У (к -тк)) +/ (1к+1, У к+1, У (4+1 -тк+1))), к п, (6)
с начальным условием у0 в момент времени (0 и начальной функцией ф (гк — тк)
при к -тк & lt- Iо и Ф (?к+1 — тк+1) при гк+1 — тк+1 & lt- го.
В случае автономного уравнения (1) функция / (г, у (г), у (г-т (г))) не зависит явно от г: f = / (у (г), у (г-т (г))). В этом случае уравнение (6) принимает более простой вид:
к —
Ук+1 = Ук + ^(f (ук& gt-у (-тк^+/(ук+1 & gt-у (гк+1 -тк+1))), к п. (7)
Ограничимся в дальнейшем случаем автономного уравнения (1). Это не нарушает общности рассуждений, так как любое неавтономное уравнение вида (1) можно свести к автономной системе уравнений, расширив состояние у1 (() = у (() дополнительным состоянием у2 (/) = (, так что расширенное состояние примет вид двумерного вектора состояния у (г) = (у1 (г) у2 (г)). Тогда уравнение (1) станет системой уравнений автономного вида? у (г)/?1 = / (у (г), у (г -т (г)), где
/ = / (((), У1 (-т (/))) У2 (г)) — двумерная вектор-функция.
Для автономного дифференциального уравнения в качестве начальной функция ф (/) удобно взять функцию ф (?) = у*, тождественно равную постоянному равновесному значению решения (точке покоя) у*, обращающему в нуль правую часть автономного дифференциального уравнения ЖуУ /Ж = / (у*, у*) = 0.
Уравнение (7) обобщает известное соотношение метода Эйлера с уравниванием на случай переменного тк. Как и в схеме Эйлера с уравниванием, уравнение (6) будем решать при каждом фиксированном к относительно ук+1 методом последовательных приближений (простых итераций), используя в качестве начального приближения у^ решение по схеме метода Эйлера, использующего вычисление интеграла в правой части (3) по формуле прямоугольников:
•Ук+1 = У к +? {у к, у у-к)), (8)
Ук+10 = Ук + 2(/(Ук& gt-У (Ь-тк)) +/(У*+[, УМ (*+1 -т*+1))], v = 0,1,2,… (9)
2
Окончание процесса последовательных приближений происходит по достижении требуемой абсолютной точности решения е:
И--0 — уй| & lt-е. (10)
Упрощенным вариантом итерационной схемы (9−10) является ее первое приближение у^ (обобщение метода Хьюна [8]).
Поскольку тк может принимать любое неотрицательное значение, не обязательно кратное шагу к, моменты времени гк — тк также не обязательно будут кратными к. Для приближенного вычисления значения у {1к+1 -тк+1) через значения в точках, ближайших слева и справа к точке ^+1 — тк+1, но кратных к, проведем линейную интерполяцию решений внутри соответствующих отрезков длиной к:
у (+1- т+1) = у к +^1 --ту) у к+1, 0 *т к+1& lt- ь-
У (к+1 -Тк+1) = [1 -+1 ~%1+1 ~-т ] У т + к+1 +1 — ^ Ут+1 ,
к — тк+1 — гк+1 — г0, %+1 — тк+1 е [ш ' +1 ] -
У (к+1 -тк+1) = У, Тк+1 & lt- ?к — ?0. (11)
Как видим, интерполяция производится на том (первом при 0 & lt- т к+1& lt- к или т-м при к — тк+1 — гк+1 — г0) интервале, куда попадает гк+1 — тк+1.
Интерполяционная схема (11) обобщает на случай переменного запаздывания тк интерполяционную схему, предложенную в [4] для постоянного запаздывания т.
3. Порядок точности обобщенного метода Эйлера с уравниванием и интерполяцией
Найдем теперь порядок точности обобщенного на случай уравнений с переменным запаздыванием метода Эйлера с уравниванием и его упрощенного варианта — метода Хьюна. Воспользуемся определениями порядка точности и порядка аппроксимации вычислительных схем решения дифференциальных уравнений, приведенными, например, в [6]. Порядком аппроксимации (соответственно порядком точности) вычислительной схемы называют минимальную степень шага к в разложении в ряд Маклорена невязки вычислительной схемы метода численного решения (соответственно невязки численного решения). Невязкой вычислительной схемы называется разность между левой и правой частями схемы метода на точном решении, деленная на величину шага к на каждом к-м шаге решения. Невязкой решения называется разность между приближенным и точным решением на каждом к-м шаге решения. Можно показать, что порядок точности совпадает с порядком аппроксимации. Поэтому порядок аппроксимации часто называют также порядком точности. Конструктивный алгоритм вычисления порядка аппроксимации вытекает из определения невязки метода.
Запишем выражение для невязки метода Эйлера с уравниванием в соответствии с определением, обозначив точное решение дифференциального уравнения для к-го шага (в к-й момент времени) через ик = и (?к), где и (г) является точным решением уравнения (1). Очевидно, ик подчиняется точному интегральному аналогу (4) дифференциального уравнения (1). Для автономного уравнения (1) рекуррентное соотношение (4) принимает вид
*г +1 _
ик+1 = ик + | / (и (), и (-Т ()))& amp-, к = 1, п -1. (12)

Невязка метода Эйлера с уравниванием и интерполяцией при случайных запаздываниях вида (3) принимает вид
ч к (к) = К+1- ик- (к/2)[ / (ик& gt-и (ч -тк)) +/ (ик+1-и (гк+1 1 Ш/к, (13)
где от шага к на к-м шаге зависят ик+1 через верхний предел гк+1 = гк + к точного равенства (12) и и (гк+1 -тк+1) по формулам (11), куда вместо всех у подставляет-
ся и. При этом следует иметь в виду, что от к на к-м шаге зависят только ик+1 и гк+1 = гк + к. Величины к, стоящие в знаменателях формул (11), относятся не к к-му интервалу, а к более ранним интервалам, поэтому они в расчете порядка точности участия не принимают.
Разложив в ряд Маклорена правую часть равенства (12), получим с точностью до членов 4-го порядка малости
ик+1 = ик + /кк + (к (к + ЛкСт)) + [/11 к /к + 2/12к/кст + /22с1 +
к3
+/'-('-к/к + /кСт)] - + О (к4), (14)
3!
где обозначено:
/к = /(«к& gt-«(гк-Тк)), /1 '-к = д/(«к& gt-«(к-Тк))/5"к ,
/2к = д/ («к & gt-«(к — Тк))/(~хк), /1!* = д2/ (ик, и (1к — тк))/ди1 ,
/12к = д2/ («к & gt-«(к — Тк))//икд» (к -Тк), /22к = д2/(«к & gt-» (к — Тк))/д» (к -Тк)2 ,
Ст = (ит+1 — ит)/^т,к — Тк+1 ^ [т ' ^т+1 ]. (15)
В последнем равенстве через Нт обозначена длина т-го, более раннего, чем к-й, интервала [1т, 1т+1 ].
Разлагая аналогичным образом в ряд Маклорена функцию
/ («к+1, и (1к+1 тк+1)) во втором слагаемом выражения (13) как сложную функцию к, с учетом (12) и
(11) (с заменой у на и) и с учетом обозначений (15) получим
f (ик+1, и (к+1 -Хк +1)) = /к + (/1к/к + /2кст) к +
+ [/пк /к + 2/12к/кст + /22ск + Л (к (к + Лкст)]) + 0(). (16)
Подставляя выражения (14) и (16) с учетом обозначений (15) в (13), приведем ^?к (к) к виду
к (А) = 12[/к2 + 2/ш: /кС» + /^ + /'- (/?к/к + /1кст)] к1 + О (к3). (17)
Из разложения (17) видно, что порядок точности метода Эйлера с уравниванием и интерполяцией остался равным 2.
Аналогичным образом можно получить разложением в ряд Маклорена по степеням к на к-м шаге невязку метода Хьюна
к (А) = К-и -«к -(А/2)[/(«к& gt-«(?к -тк)) +
+/(ик + к/ (ик, и (1к -хк)), и (1к+1 -хк+1))]^к. (18)
Проводя разложение второго слагаемого в (18) как сложной функции к и используя разложение (14), с учетом обозначений (15) получим
к (А) = 12[-^к /2 — 2/2к/кс» — /?Ск2 +2// (/к/к + /2кСт)] А2 + О (к3). (19)
Отсюда видно, что порядок точности метода Хьюна с интерполяцией также остается равным двум. Заметим, что в случае линейных дифференциальных уравнений, для которых функция f линейна по своим аргументам, вторые производные в выражениях (17) и (19) обращаются в нуль, и мы видим, что функция? k (h)
для метода Хьюна оказывается при малых h ровно в 2 раза больше, чем для метода Эйлера с уравниванием, что свидетельствует все-таки о несколько меньшей скорости сходимости к точному решению вычислений по неитерационной схеме Хьюна по сравнению с итерационной схемой Эйлера с уравниванием.
ЛИТЕРАТУРА
1. Эльсгольц Л. Э. Введение в теорию дифференциальных уравнений с отклоняющимся аргументом. М.: Наука, 1964. 128 с.
2. Солодов А. В., Солодова Е. А. Системы с переменным запаздыванием. М.: Наука, 1980. 384 с.
3. Тихонов А. Н., Васильева А. Б., Свешников А. Г. Дифференциальные уравнения. М.: Наука, 1980. 233 с.
4. Поддубный В. В. Численное решение дифференциальных уравнений с постоянным запаздыванием методом Эйлера с уравниванием и интерполяцией // Обработка данных и управление в сложных системах. Вып. 7. Томск: Изд-во Том. ун-та, 2005. С. 165 — 174.
5. Эльсгольц Л. Э. Дифференциальные уравнения и вариационное исчисление. М.: Наука, 1965. 324 с.
6. Бабенко К. И. Основы численного анализа. М.: Наука, 1986. 744 с.
7. Самарский А. А., Гулин А. В. Численные методы. М.: Наука, 1989. 432 с.
8. Ортега Дж., Пул У. Введение в численные методы решения дифференциальных уравнений. М.: Наука, 1986. 288 с.
9. Коллатц Л. Численные методы решения дифференциальных уравнений. М.: ИЛ, 1953.
Статья представлена кафедрой прикладной информатики факультета информатики Томского государственного университета, поступила в научную редакцию 29 июня 2007 г.

ПоказатьСвернуть
Заполнить форму текущей работой