ВЕСТНИК САНКТ-ПЕТЕРБУРГСКОГО УНИВЕРСИТЕТА
Сер. 10. 2009. Вып. 4
УДК 517.97 Г. Ш. Тамасян
ГРАДИЕНТНЫЕ МЕТОДЫ РЕШЕНИЯ ЗАДАЧИ КОШИ*)
Введение. Для решения задачи Коши в настоящее время известно множество численных методов, например метод последовательных приближений Пикара, метод Эйлера, метод Рунге-Кутта [1]. В данной работе решение задачи Коши сводится к безусловной минимизации соответствующего функционала. С учетом специфики строения функционала для поиска минимизирующей последовательности применяются градиентные методы. Рассматриваемые ниже алгоритмы относятся к прямым методам вариационного исчисления [2, 3].
Постановка задачи Коши. Пусть Т > 0 - фиксированное число. Рассмотрим линейную неоднородную систему дифференциальных уравнений
х = Р (г)х + д(г), (1)
в которой х - п-мерный вектор искомых функций; Р(г) - (п х п)-матрица; д(г) - «.-мерная вектор-функция. Предположим, что элементы матрицы Р(г) и вектор-функции д(г) являются функциями, определенными и непрерывными при г € [0, Т].
Требуется среди всех решений системы (1) найти такое, которое будет удовлетворять начальному условию х(0) = хо, где хо € К" - заданный вектор. С учетом указанных требований к системе (1) решение задачи Коши существует и единственно [1].
Постановка вариационной задачи. Произведем следующую замену переменных:
г
т =хо+/4г > ¿т, (2)
о
где г (г) - непрерывная вектор-функция при г ^ 0. Отметим, что г(г) = х(г).
Итак, решение задачи Коши сводится к поиску такой вектор-функции г(г), которая удовлетворяет системе
г(г) = Р(г) ^хо + I г(т) ¿т^ + д(г). (3)
г
Далее для краткости записи будем писать х(г) вместо хо + J г(т) ¿т (см. (2)).
Тамасян Григорий Шаликович — кандидат физико-математических наук, доцент кафедры математической теории моделирования систем управления факультета прикладной математики—процессов управления Санкт-Петербургского государственного университета. Количество опубликованных работ: 17. Научные направления: негладкий анализ, вариационное исчисление, недифференцируемая оптимизация. E-mail: [email protected], [email protected].
*) Работа выполнена при финансовой поддержке Российского фонда фундаментальных исследований (грант № 09-01-00360).
© Г. Ш. Тамасян, 2009
Обозначим через Сп[0,Т] пространство непрерывных вещественных п-мерных вектор-функций на отрезке [0, Т]. В этом пространстве введем норму
(z(t),z(t)) dt.
Поиск решения системы (3) сведем к минимизации следующего функционала на всем пространстве Сп [0, Т].
Рассмотрим функционал и некоторые его свойства:
z
т
I (z) = \\z(t) - P (t)x(t) - g(t)\\2 = J (z(t) - P (t)x(t) - g(t), z(t) - P (t)x(t) - g(t)) dt,
0
t
где x(t) = xo + j z(t) dr.
0
Заметим, что I(z) ^ 0 при любых z G Cn[0, T]. Несложно также показать, что функционал I(z) достигает минимального значения, равного нулю I(z*) = 0, тогда и только тогда, когда z*(t) - решение задачи Коши (1) или (3).
Покажем, что I(z) является строго выпуклым функционалом [1, 3], доказав прежде следующие вспомогательные утверждения.
Обозначим через f (z) = z(t) - P(t)x(t) - g(t).
Утверждение 1. Пусть zi(t) и z2(t) - произвольные элементы из пространства Cn[0,T]. Для того чтобы f (zi) = f (z2) для всех t G [0,T], необходимо и достаточно, чтобы zi(t) = z2(t).
Доказательство. Достаточность очевидна. Докажем необходимость. Пусть f (zi) = f (z2), но zi(t) = z2(t). Тогда найдется такая функция ф(t) ф 0 из пространства Cn[0,T], что zi(t) = z2(t)+^(t).
Итак,
f (zi) = f (z2 + ф) = z2(t) + m - P(t)
t
xo + J(z2(r )+ф(т ))dr
- g(t) =
t
= f Ы + ф^) - P(t) J ф(т)dr.
Так как по условию f (zi) = f (z2), то
t
ф(t) - P(t) J ф(т)dr = 0.
Полученное выражение - однородное интегральное уравнение Вольтерра второго рода. При £ = 0 имеем ф(0) = 0, значит, единственным решением интегрального уравнения является ф(Ь) = 0. Получили противоречие.
Утверждение 2. Функционал I(z) строго выпуклый. Доказательство. Требуется показать, что
al(zi) + (1 - a)I(z2) - I(azi + (1 - a)z2) > 0
для всех a е (0,1), zi = z2, zi, z2 G Cn[0,T]. Действительно,
I(azi + (1- a)z2) = \\f (azi + (1 - a)z2)\\2 =
azi(T) + (1- a)z2(r) - P(t)
t
X0 + J(azi(t) + (1 - a)z^{T))¿t
- g(t)
= \\af (zi) + (1- a)f (z2)\\2 =
1 11
= a2 J(f (zi), f(zi)) dt + 2a(1 - a) j(f (zi), f(z2)) dt + (1 - a)2 j(f (z2), f(z2)) dt =
0 0 0 T T T
= a j(f (zi), f(zi)) dt + (1 - a) j(f (z2), f(z2)) dt + 2a(1 - a) j(f (zi), f(z2)) dt -
0 0 0 T T
- a(1 - a) j(f (zi), f(zi)) dt - a(1 - a) j(f Ы, f(z2)) dt =
00
T
= a\\f (zi)\\2 + (1- a)\\f (z2)\\2 - a(1 - a) j(f (zi) - f (z2), f(zi) - f Ы) dt.
0
Используя утверждение 1, имеем
aI (zi) + (1 - a)I (z2) - I (azi + (1 - a)z2) =
T
= a(1 - a) j(f (zi) - f (z2), f(zi) - f (z2)) dt > 0
0
для всех a G (0,1), zi = z2.
Итак, решение задачи Коши будем искать среди функций, доставляющих минимальное значение (равное нулю) функционалу
T
I(z) = J (z(t) - P(t)x(t) - g(t), z(t) - P(t)x(t) - g(tf) dt
min (4)
zecn[0,T ]
на всем пространстве Сп [0, Т].
Используя теорию точных штрафных функций, в [2] были получены необходимые условия минимума для поставленной проблемы.
Описание алгоритмов. В классическом вариационном исчислении наиболее распространенный поиск экстремалей функционала осуществляется либо из решений уравнений Эйлера [4], либо при помощи прямых методов [5, 6] типа Ритца, Галеркина,
2
Канторовича и т. п. Описанные ниже алгоритмы являются аналогами градиентных методов [2, 7], которые применяются в теории оптимизации конечномерных пространств, а именно методы наискорейшего спуска и сопряженных направлений. Указанные методы (в конечномерных задачах) обладают достаточно хорошей скоростью сходимости минимизирующей последовательности [7-9]. К примеру, в задачах минимизации выпуклых квадратичных функций метод сопряженных направлений сходится за конечное число шагов, которое не превосходит размерности пространства. Различные варианты градиентных и прямых методов решения задач вариационного исчисления можно найти в работах Б. Т. Поляка [3], Л. В. Канторовича [10], С. Г. Михлина [6], В. Ф. Демьянова [2] и др.
Рассмотрим классическую вариацию вектор-функции г(г):
ге(г) = г(г) + £у(г), (5)
где £ > 0, у(г) - произвольная вектор-функция пространства Сп[0,Т].
Применяя вариацию (5) к функционалу (4), получим классическую вариацию функционала I (г)
I(гЕ) = I(г) + £ I (С(г, г), у(г)) ¿г + о(£),
здесь ^ ^ -> 0, 6? (г, Ь) - градиент Гато [2, 3, 5] функционала 1{г) в точке г имеет вид
£ е|0
Т
С(г, г) = г(г) - Р(г)х(г) - д(г) - | Р* (т) [г(т) - Р(т)х(т) - д(т)] ¿т.
г
Через Р*(г) обозначается транспонированная матрица Р(г).
Далее рассмотрим несколько разновидностей градиентного метода. Изучаемые алгоритмы носят итерационный характер. Это значит, что строится некоторая последовательность ги(г), к = 0, 1,..., относительно которой можно утверждать, что она в той или иной мере сходится к решению задачи минимизации. Так как функционал (4) строго выпуклый, то равенство нулю нормы градиента в точке г* (стационарная точка) является необходимым и достаточным условием минимума функционала [2, 8, 9].
Градиентный метод. Выбираем начальное приближение г0(г). Строим последовательность {гь(г)} по правилу
гк+1(г) = гк(г) - ,г), 7к > 0, к = 0,1,...,
если \\С(гк,г)\\ =
Т
J (о(ги ,г), 0(ги ,г)^ ¿г = 0, то шаг 7к > 0 [8, 9] можно выбрать так,
чтобы I(ги+1) < I(ги).
Если \\С(ги,г)\\ = 0, то ги - стационарная точка, и процесс построения минимизирующей последовательности прекращается.
Метод наискорейшего спуска (МНС). Пользуясь тем, что градиент С(г,г) линейный относительно г, а подынтегральная функция функционала (4) квадратичная
относительно г, можно модифицировать предыдущий метод. А именно, шаг ^к > 0 будем выбирать из условия минимума функционала (4). Несложно показать, что таким шагом является величина
т г
У (хк(г) - Р(г)хк(г) - д(г),о(гк,г) - Р(г) ^о(гк,т) ¿^¿г
о
7й+1 — —
У (а(*к ,г) - Р (г)I С(гк,т) ¿т, о (г к ,г) - Р (г)! о(гк ,т) ¿^¿г о о о
Метод сопряженных градиентов (направлений) (МСГ). Пусть г0(г) - некоторое начальное приближение. Будем строить последовательность {гк(г)} по формулам
гк+1(г)= гк (г) - Чк Шк(г), к = 0, 1,..., (6)
где
Шо(г) = о(го ,г), Шк (г) = о(гк ,г)+вк Шк-1(г), к = 1,2,....
Величина 7к может определяться так же, как и в вышеуказанных методах, а вк - по одной из формул
т
У (о(гк,г),о(гк,г) - о(гк-1 ,г))¿г
13к = 2-?-, (7)
[ (а(гк-1 ,г),о(гк-1 ,г)) ¿г
о
т
I (о(гк ,г),о(гк ,г) )аь
о
¡Зк = -, (8)
о
т
У (о(гк-1 ,г),о(гк-1 ,г)) ¿г
У (а(гк,г),о(гк,г) - о(гк-1 ,г))¿г
13к = 2-?-, (9)
У (Шк-1 ,о(гк-1,г)) ¿г
о
т
У (а(гк ,г),о(гк ,г)) ¿г /?* = \-. (10)
Шк-1 ,о(гк-1,г)) ¿г
г
г
Как показывает практика, в каждом шаге алгоритма с неизбежностью накапливаются погрешности. Это может привести к тому, что векторы Шк перестают указывать направление убывания функционала, и релаксация метода может нарушиться. Для борьбы с таким явлением метод сопряженных направлений время от времени обновляют, полагая в (6) вк = 0, т. е. осуществляют градиентный спуск.
Заметим, что величины, вычисленные по формулам (7)—(10), будут одинаковыми, если градиенты попарно ортогональные.
Пример. Пусть Т =15. Решим задачу Коши
при начальном условии х(0) = 1.
Найдем решение данной задачи х*(£) = методом последовательного приближения Пикара. На к-м шаге получаем первые к членов разложенного в ряд Тейлора «искомого» решения, а именно
£2 к £в Х1(*) = 1, х2(*) = 1-*, хз(*) = 1-* + -, ••• = Х)(-1)в"г
2 8 = 1
Очевидно, что хк(£) -> х*(Ь) = .
к—
В табл. 1-3 приведены результаты решения задачи Коши алгоритмами Пикара, МНС и МСГ на отрезке [0,15]. Вычисления проводились в математическом пакете МаШса^ В качестве начального приближения в методе Пикара взято хо(£) = 1, а в алгоритмах МНС и МСГ выбрали ¿о(£) = 1. В табл. 1 представлены выборочные итерации алгоритма Пикара; в табл. 2 и 3 - первые 5 шагов соответствующих алгоритмов.
Таблица 1. Результаты алгоритма Пикара
к Р* - *к\ 1 \\Х || 1(*к)
7 1.168 • 104 2.224 • 104 3.38 • 104
14 1.204 • 10Б 1.205 • 10Б 2.4 • 10Б
21 3.4 • 104 2.3 • 104 5.7- 104
28 947 486.9 1440
35 4.58 1.897 6.48
40 0.043 0.016 0.058
Таблица 2. Результаты алгоритма МНС
к Р* - *к\ 1 II^к || 1(*к) ЦСЫ11
0 4.183 36.894 1635 824.327
1 1.328 2.597 11.527 20.506
2 0.732 1.043 1.911 7.66
3 0.678 1.011 1.594 8.05
4 0.613 0.828 1.34 6.09
5 0.6 0.823 1.14 6.48
Таблица 3. Результаты алгоритма МСГ
к - ¿к\\ \\х || /ы ЦСЫ11 №
0 4.183 36.894 1635 824.327 824.327
1 1.328 2.597 11.527 20.506 20.512
2 0.73 1.06 1.854 5.938 6.17
3 0.34 0.283 0.21 1.32 1.35
4 0.096 0.067 0.014 0.305 0.313
5 0.028 0.023 0.00132 0.105 -
Из примера видно, что рассматриваемые методы являются релаксационными.
Заключение. Результаты численных экспериментов показали эффективность использованных методов для решения задачи Коши. Метод сопряженных градиентов более трудоемкий по сравнению с методами наискорейшего спуска и градиентным, однако построенная минимизирующая последовательность имеет более высокую скорость сходимости. Указанные алгоритмы могут быть применены и для решения нелинейных систем дифференциальных уравнений. Остаются открытыми вопросы о сходимости и скорости сходимости построенных минимизирующих последовательностей.
Автор благодарит В. Ф. Демьянова и С. К. Мышкова за ценные указания и внимание к работе.
Литература
1. Матвеев Н. М. Методы интегрирования обыкновенных дифференциальных уравнений. Изд. 4-е, испр. и доп. Минск: Вышэйшая школа, 1974. 768 с.
2. Демьянов В. Ф. Условия экстремума и вариационные задачи. М.: Высшая школа, 2005. 335 с.
3. Поляк Б. Т. Градиентные методы минимизации функционалов // Журн. вычисл. математики и мат. физики. 1963. Т. 3, № 4. С. 643-653.
4. Гюнтер Н. М. Курс вариационного исчисления. М.: Гостехиздат, 1941. 308 с.
5. Вайнберг М. М. Вариационный метод и метод монотонных операторов в теории нелинейных уравнений. М.: Наука, 1972. 416 с.
6. Михлин С. Г. Численная реализация вариационных методов. М.: Наука, 1966. 432 с.
7. Карманов В. Г. Математическое программирование. М.: Наука, 1986. 228 с.
8. Пшеничный Б. Н., Данилин Ю. М. Численные методы в экстремальных задачах. М.: Наука, 1975. 320 с.
9. Васильев Ф. П. Численные методы решения экстремальных задач. М.: Наука, 1980. 550 с.
10. Канторович Л. В., Акилов Г. П. Функциональный анализ. 2-е изд., перераб. М.: Наука, 1977. 741 с.
Статья рекомендована к печати В. Ф. Демьяновым. Статья принята к печати 28 мая 2008 г.