Алгоритм бройдена для решения систем нелинейных уравнений

Курсовая работа: Решение систем нелинейных уравнений методом Бройдена

Федеральное агентство по образованию

ГОСУДАРСТВЕННОЕ ОБРАЗОВАТЕЛЬНОЕ УЧРЕЖДЕНИЕ

ВЫСШЕГО ПРОФЕССИОНАЛЬНОГО ОБРАЗОВАНИЯ

«ВОРОНЕЖСКИЙ ГОСУДАРСТВЕННЫЙ ТЕХНИЧЕСКИЙ УНИВЕРСИТЕТ»

Факультет автоматики и электромеханики

Кафедра «Автоматизированные и вычислительные системы»

Специальность «Вычислительные машины, комплексы, системы и сети»

по дисциплине «Вычислительная математика»

Тема работы «Решение систем нелинейных уравнений методом Бройдена»

Пояснительная записка 26 с., 14 рисунка, 2 источника. Ключевые слова: МЕТОД БРОЙДЕНА, РЕШЕНИЕ СИСТЕМ МЕТОДОМ БРОЙДЕНА, РЕШЕНИЕ СИСТЕМ НЕЛИНЕЙНЫХ УРАВНЕНИЙ.

Объект исследования или разработки – решение систем нелинейных уравнений методом Бройдена.

Цель работы – создать программу, иллюстрирующую решение систем нелинейных уравнений методом Бройдена и исследовать результат ее работы.

Полученные результаты – листинг полученный программы, проверка соответствия найденных решений точным решениям заданной системы нелинейных уравнений.

Основные конструктивные, технологические и технико-эксплуатационные характеристики — персональная ЭВМ.

1. Алгоритм бройдена

1.1 Входные данные для алгоритма Бройдена

1.2 Содержание алгоритма Бройдена

1.3 Метод исключения Гаусса для решения СЛАУ

1.4 Вывод формулы пересчета Бройдена

2. Разработка программы и иследование результата ее работы

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

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

,

где — нелинейные функции от , в тождество.

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

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

В курсовой работе будет рассматриваться метод решения Бройдена для систем нелинейных уравнений.

1. АЛГОРИТМ БРОЙДЕНА.

1.1 Входные данные для алгоритма Бройдена

Входными данными для алгоритма Бройдена являются вектор начального решения, начальная матрица Якоби и заданная точность.

1. 2 Содержание алгоритма Бройдена

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

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

,

где . И весь процесс поиска решения повторяем по той же самой схеме до тех пор, пока не будет получено решение c заданной точностью [1].

Поскольку необходимо решить линейное уравнение, то рассмотрим метод решения Гаусса.

1.3 Метод исключения Гаусса для решения СЛАУ

Суть всех методов исключения состоит в приведении исходной системы уравнений к системе более простого вида, для которой легко найти решение. К этим методам можно отнести метод исключения Гаусса, который имеет много вычислительных схем и, как показали исследования, является идеальным алгоритмом для решения СЛАУ.

Рассмотрим сначала самую простую схему – схему единственного деления. Применение схемы единственного деления продемонстрируем на примере СЛАУ 4- го порядка

Разделив первое уравнение системы на , получим

Из второго уравнения системы вычтем первое, умноженное на коэффициент при , то есть на . В результате получаем:

=

Поступая таким же образом с третьим и последующими уравнениями системы, получим

;

;

.

К выделенной системе применим тот же алгоритм, что и к исходной. В результате получаем

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

, , .

1.4 Вывод формулы пересчета Бройдена

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

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

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

Равенство , переписанное в виде , называют соотношением секущих в Оно легко обобщается на n -мерный случай и лежит в основе вывода метода Бройдена. Опишем этот вывод.

В n-мерном векторном пространстве соотношение секущих представляется равенством

,

где — известные n-мерные векторы, — данное нелинейное отображение, а — некоторая матрица линейного преобразования в . С обозначениями , соотношение секущих в обретает более короткую запись . Аналогично одномерному случаю, а именно, по аналогии с формулой , будем искать приближения к решению векторного уравнения по формуле . Обратимую n x n-матрицу в ней нужно подобрать так, чтобы она удовлетворяла соотношению секущих . Но это соотношение не определяет однозначно матрицу : глядя на равенство , легко понять, что при n>1 существует множество матриц , преобразующих заданный n-мерный вектор в другой заданный вектор (отсюда — ясность в понимании того, что могут быть различные обобщения одномерного метода секущих).

При формировании матрицы будем рассуждать следующим образом. Переходя от имеющейся в точке аффинной модели функции F(x) к такой же модели в точке мы не имеем о матрице линейного преобразования никаких сведений, кроме соотношения секущих . Поэтому исходим из того, что при этом переходе изменения в модели должны быть минимальными. Эти изменения характеризует разность . Вычтем из равенства определяющее равенство и преобразуем результат, привлекая соотношение секущих . Имеем:

Представим вектор в виде линейной комбинации фиксированного вектора определенного в , и некоторого вектора t, ему ортогонального: ,

Подстановкой этого представления вектора в разность получаем другой ее вид

Анализируя выражение , замечаем, что первое слагаемое в нем не может быть изменено, поскольку — фиксированный вектор при фиксированном k . Поэтому минимальному изменению аффинной модели будет отвечать случай, когда второе слагаемое в будет нуль-вектором при всяких векторах t, ортогональных векторам , т.е. следует находить из условия

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

Таким образом, приходим к так называемой формуле пересчета С. Бройдена

2. РАЗРАБОТКА ПРОГРАММЫ И ИСЛЕДОВВАНИЕ РЕЗУЛЬТАТА ЕЕ РАБОТЫ

Задача. Разработать программу, реализующую метод Бройдена.

Структура программы. Программа была разработана в интегрированной среде разработке приложений Microsoft Visual Studio 2008 на языке программирования C#, проект программы Console Application. В ходные данные программы начальный вектор решения, начальная матрица Якоби и удовлетворяющая погрешность. Программа решает систему уравнений . Если программа не находит решения удовлетворяющего требуемой точности за 10 итераций, то поиск решения прекращается, а так же если процесс расходится (в соответствии с приложением А).

Введем матрицу Якоби , погрешность 0,3 начальное решение является точным решение. На 1 итерации получаем результат решения (рисунок 1).

Рисунок 1 – Первый пример работы программы

Результат точное решение на 1 шаге. Попробуем задать начальное решение отличное от точного (рисунок 2).

Рисунок 2 – второй пример работы программы

Получили близко решение к точному решению. Попробуем уменьшить погрешность (рисунок 3).

Рисунок 3 – третий пример работы программы

Получили точное решение. Попробуем сильнее отойти в начальном решении от точного (рисунок 4).

Рисунок 4 – Четвертый пример работы программы

Получаем точное решение. Уменьшим погрешность и сильнее отойдем от точного решения. Теперь начальное решение произвольное (рисунок 5).

Рисунок 5 – Пятый пример работы программы

Видим увеличение количества итераций. Решение получили точное. Немного изменим начальную матрицу Якоби (рисунок 6).

Рисунок 6 – Шестой пример работы программы

Увеличение количества итераций. Решение точное. Теперь возьмем другую матрицу Якоби (рисунок 7).

Рисунок 7 – Седьмой пример работы программы

Получили плохой результат решения. Попробуем выяснить из-за чего. Или матрица Якоби в начале исследования была близка к расчетной матрицы Якоби на основе конечно разностной аппроксимации производных или при таком начальном решении требуется слишком много итераций.

Попробуем с начальной матрицей Якоби. Процесс решения стал расходится. Делаем вывод, что не смогли найти решения из-за начального решения (рисунок 8).

Рисунок 8 – Восьмой пример работы программы

На основе рисунка 9, рисунка 10 и рисунка 11 видим, что наша первая матрица Якоби была удачно выбрана.

Матрица Якоби близка к первой матрице Якоби (рисунок 12).

Рисунок 9 – Девятый пример работы программы

Рисунок 10 – Десятый пример работы программы

Рисунок 11 – Одиннадцатый пример работы программы

Рисунок 12 – Двенадцатый пример работы программы

Попробуем изменить систему уравнений, решаемую программой и посмотрим на результаты работы программы (рисунок 13,14).

Рисунок 13 – Тринадцатый пример работы программы

Рисунок 14 – Четырнадцатый пример работы программы

ЗАКЛЮЧЕНИЕ

Если начальное приближение выбрано достаточно близко к решению и если начальная аппроксимация матрицы Якоби достаточно точна, то метод Бройдена обладает сверхлинейной сходимостью, но не квадратичной, как метод Ньютона.

Данная курсовая работа выполнена в полном объеме. В курсовой работе был рассмотрен метод Бройдена, написана программа реализующая его.

СПИСОК ЛИТЕРАТУРЫ

1. С.Л. Подвальный, Л.В. Холопкина. Вычислительная математика- учебное пособие ВГТУ, 2004 — 147 с.

2. Методы решения систем нелинейных уравнений. Метод Ньютона. Его реализации и модификации. — Электрон. дан. – Режим доступа: www.exponenta.ru/educat/referat/XVkonkurs/15/index.asp.

ПРИЛОЖЕНИЕ

/*программа предназначена для решения системы нелинейных уравнений.

Программа выполнена 1 ноября 2009 года. Обем необходимой памяти для работы составляет 124 КБ. Версия программы №1.Автор Харитонова Яна Андреевна.*/

static void Main(string[] args)

Console.WriteLine(«x+y-3» + «\n» + «x^2+y^2-9»);

double[,] yakob = new double[N, N];

Console.WriteLine(«введите элементы матрицы Якоби»);

for (int i = 0; i = 0; i—)

for (int j = i + 1; j = 0)

double[] Y = new double[N];

Y[0] = (XJ[0] + XJ[1] — 3) — Bnach[0];

Y[1] = (XJ[0] * XJ[0] + XJ[1] * XJ[1] — 9) — Bnach[1];

«Визуализация процесса нахождения корней нелинейных уравнений методом Бройдена»

ФЕДЕРАЛЬНОЕ АГЕНТСТВО ПО ОБРАЗОВАНИЮ

ГОУВПО «Самарский государственный архитектурно-строительный университет»

Кафедра прикладной математики и вычислительной техники

ПОЯСНИТЕЛЬНАЯ ЗАПИСКА К КУРСОВОЙ РАБОТЕ

ТЕХНОЛОГИЯ ПРОФЕССИОНАЛЬНОЙ ДЕЯТЕЛЬНОСТИ

«Визуализация процесса нахождения корней нелинейных уравнений методом Бройдена»

5 СЕМЕСТР 3 КУРС

Научный руководитель (ФИО) Бойцев Александр

Методический руководитель (ФИО)

студент ГИП-107 Давыдова Евгения

Оценка преподавателя _______________

Оценка комиссии по результатам защиты_______________

Наука в целом (информационные технологии — 004)Информационные технологии.

Компьютерные технологии. Теория вычислительных машин и систем

Методы решения задач

Алгоритмы составления программ

Алгоритм, Бройдена, метод, решение, уравнения, информационная, система

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

Экран оценки творческого уровня работы

Существует несколько методов решения нелинейных уравнений. В данной работе рассмотрены такие, как метод деления отрезка пополам, метод простой итерации, метод Ньютона (касательных), метод хорд, комбинированный метод и метод Бройдена. Метод секущих Бройдена является одним из оптимальных методов решения нелинейных уравнений. По-этому, возникла задача, целью которой является написание программы, визуализирующую процесс нахождения корней нелинейных уравнений методом Бройдена.

Наиболее простым методом, позволяющим найти корень нелинейного уравнения, является метод половинного деления. Сходимости мала, но к достоинствам метода относятся простота и безусловная сходимость итерационного процесса. Метод Ньютона (Касательных) имеет квадратичную скорость сходимости. Он является одним из оптимальных методов, так же как и комбинированный метод приближения корня уравнения по методу хорд и по методу касательных подходят к значению этого корня с противоположных сторон. Поэтому для быстроты нахождения корня удобно применять оба метода одновременно. Т. к. один метод даёт значение корня с недостатком, а другой – с избытком, то достаточно легко получить заданную степень точности корня.

Решение нелинейных уравнений комбинированным методом возможно только при выполнении следующих условий:

1) ,

2) и сохраняют знак на отрезке ,

Приближения корней находятся:

а) по методу касательных: ,

б) по методу хорд: .

При схождении к корню с разных сторон скорость схождения получается максимальной.

В процессе написания программы были реализованы два способа графического представления данного метода (пошаговый и стандартный).

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

При исследовании сходимости методов решения нелинейных уравнений метод Хорд и касательных оказался наиболее оптимальным из выбранных методов. На втором месте оказался метод Ньютона.

1. Методы решения нелинейных уравнений

1.1. Метод деления отрезка пополам

1.2. Метод простой итерации

1.3. Метод Ньютона (касательных)

1.5. Комбинированный метод

1.6. Метод простейших секущих

1.7. Метод секущих Бройдена

2. Выявление наиболее оптимального метода

3. Алгоритм решения уравнения комбинированным методом

4. Программа визуализации процесса нахождения корней нелинейных уравнений комбинированным методом и методом Бройдена.

5. Список литературы

1. Методы решения нелинейных уравнений

1.1 Метод деления отрезка пополам

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

Пусть на отрезке [a, b] задана непрерывная функция Если значения функции на концах отрезка имеют разные знаки, т. е. то это означает, что внутри данного отрезка находится нечетное число корней. Пусть для определенности корень один. Суть метода состоит в сокращении на каждой итерации вдвое длины отрезка. Находим середину отрезка [a, b] (см. рис. 1) Вычисляем значение функции и выбираем тот отрезок, на котором функция меняет свой знак. Новый отрезок вновь делим пополам. И этот процесс продолжаем до тех пор, пока длина отрезка не сравняется с наперед заданной погрешностью вычисления корня e. Построение нескольких последовательных приближений по формуле (3) приведено на рисунке 1.

Итак, алгоритм метода дихотомии:

1. Задать отрезок [a, b] и погрешность e.

2. Если f(a) и f(b) имеют одинаковые знаки, выдать сообщение о невозможности отыскания корня и остановиться.

Рис.1. Метод деления отрезка пополам для решения уравнения вида f(х)=0.

3. В противном случае вычислить c=(a+b)/2

4. Если f(a) и f(c) имеют разные знаки, положить b=c, в противном случае a=c.

5. Если длина нового отрезка , то вычислить значение корня c=(a+b)/2 и остановиться, в противном случае перейти к шагу 3.

Так как за N шагов длина отрезка [a, b] сокращается в 2N раз, то заданная погрешность отыскания корня e будет достигнута за итераций.

Как видно, скорость сходимости мала, но к достоинствам метода относятся простота и безусловная сходимость итерационного процесса. Если отрезок [a, b] содержит больше одного корня (но нечетное число), то всегда будет найден какой-нибудь один.

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

1.2 Метод простой итерации

При использовании этого метода исходное нелинейное уравнение (1) необходимо переписать в виде

(2)

Обозначим корень этого уравнения C*. Пусть известно начальное приближение корня . Подставляя это значение в правую часть уравнения (2), получаем новое приближение

и т. д. Для (n+1)- шага получим следующее приближение

(3)

Таким образом, по формуле (3) получаем последовательность С0, С1,…,Сn+1, которая стремиться к корню С* при n®¥. Итерационный процесс прекращается, если результаты двух последовательных итераций близки, т. е. выполняется условие

(4)

Исследуем условие и скорость сходимости числовой последовательности при n®¥. Напомним определение скорости сходимости. Последовательность , сходящаяся к пределу С*, имеет скорость сходимости порядка a, если при n®¥ выполняется условие

(5)

Допустим, что имеет непрерывную производную, тогда погрешность на (n+1)-м итерационном шаге en+1=Cn+1-C*=g(Cn)-g(C*) можно представить в виде ряда

en+1 » Cn+1 – C* = g¢(C*) (Cn-C*) +¼@ g¢(C*) en+¼

Таким образом, получаем, что при выполнении условия

Рис. 2. Графическая интерпретация метода простых итераций для решения уравнения вида x=g(х).

Построение нескольких последовательных приближений по формуле (3)

приведено на рисунке 2.

1.3 Метод Ньютона (Касательных)

В литературе этот метод часто называют методом линеаризации. Выбираем начальное приближение С0. Допустим, что отклонение С0 от истинного значения корня С* мало, тогда, разлагая f(C*) в ряд Тейлора в точке С0 , получим

f(C*) = f(C0) + f ¢(C0) (C*-C0) +¼ (8)

Если f ¢(C0) ¹ 0 , то в (8) можно ограничится линейными по DC =C-C0 членами. Учитывая, что f(C*)=0, из (9) можно найти следующее приближение для корня

C1 = C0 – f (C0) / f¢(C0)

или для (n+1)-го приближения

Cn+1= C n – f (C n) / f ¢(C n) (9)

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

çCn+1 – Cn ç 0 и ¦(a1) 0, ¦(x2) 0, то применяем эту формулу к отрезку [x1;a1]. Повторяя этот прием несколько раз, мы будем получать все более точные значения корня а2, а3 и т. д.

1.5 Комбинированный метод (метод хорд и касательных)

Если выполняются условия:

1) ,

2) и сохраняют знак на отрезке ,

то приближения корня уравнения по методу хорд и по методу касательных подходят к значению этого корня с противоположных сторон. Поэтому для быстроты нахождения корня удобно применять оба метода одновременно. Т. к. один метод даёт значение корня с недостатком, а другой – с избытком, то достаточно легко получить заданную степень точности корня.

1.6. Метод простейших секущих

Можно связать задание последовательности () с какой-либо сходящейся к нулю векторной последовательностью, например, с последовательность невязок () или поправок (). Так, полагая где j=1,…n, a k=1,2,…,приходим к простейшему методу секущих — обобщению скалярного метода секущих (5.32):

, (4.1.1)

Этот метод является. двухшаговым и требует задания двух начальных точек и . При п = 1 сходимость метода(4.1.1) имеет порядок . Можно рассчитывать на такую же скорость и в многомерном случае.

К методу секущих так же, как и к методу Ньютона, можно применить пошаговую аппроксимацию обратных матриц на основе метода Шульца. Расчетные формулы этой модификации легко выписать, заменив в совокупности формул ААМН (3.3.1) матрицу на матрицу из (4.1.1)

1.7. Метод секущих Бройдена

В процессе построения методов Ньютона и секущих решения нелинейного скалярного уравнения

(4.4.1)

функция f(x) в окрестности текущей точки подменяется линейной функцией

(4.4.1а)

Приравнивание к нулю последней, т. е. решение линейного уравнения

,

порождает итерационную формулу

(4.4.2)

для вычисления приближений к корню уравнения (4.4.1).

Если потребовать, чтобы заменяющая функцию f(x) вблизи точки аффинная модель имела в этой точке одинаковую с ней производную, то, дифференцируя (4.4.1а), получаем Значение коэффициента

,

подстановка которого в (4.4.2) приводит к известному методу Ньютона. Если же исходить из того, что наряду с равенством должно иметь место совпадение функций f(x) и в предшествующей точке т. е. из равенства

,

или, в соответствии с (4.4.1а)

, (4.4.3)

то получаем коэффициент

,

превращающий (4.4.2) в известную формулу секущих.

Равенство (4.4.3), переписанное в виде

,

называют соотношением секущих в Оно легко обобщается на n — мерный случай и лежит в основе вывода метода Бройдена. Опишем этот вывод.

В n-мерном векторном пространстве соотношение секущих представляется равенством

, (4.4.4)

где — известные n-мерные векторы, — данное нелинейное отображение, а — некоторая матрица линейного преобразования в . С обозначениями

, (4.4.5)

соотношение секущих в обретает более короткую запись:

(4.4.4а)

Аналогично одномерному случаю, а именно, по аналогии с формулой (4.4.2), будем искать приближения к решению векторного уравнения (2.1а) по формуле

(4.4.6)

Желая, чтобы эта формула обобщала метод секущих (5.32), обратимую n x n-матрицу в ней нужно подобрать так, чтобы она удовлетворяла соотношению секущих (4.4.4). Но это соотношение не определяет однозначно матрицу : глядя на равенство (4.4.4а), легко понять, что при n>1 существует множество матриц , преобразующих заданный n-мерный вектор в другой заданный вектор (отсюда — ясность в понимании того, что могут быть различные обобщения одномерного метода секущих).

При формировании матрицы будем рассуждать следующим образом.

Переходя от имеющейся в точке аффинной модели функции F(x)

(4.4.7)

к такой же модели в точке

(4.4.8)

мы не имеем о матрице линейного преобразования никаких сведений, кроме соотношения секущих (4.4.4). Поэтому исходим из того, что при этом переходе изменения в модели должны быть минимальными. Эти изменения характеризует разность . Вычтем из равенства (4.4.8) определяющее равенство (4.4.7) и преобразуем результат, привлекая соотношение секущих (4.4.4). Имеем:

Представим вектор в виде линейной комбинации фиксированного вектора определенного в (4.4.5), и некоторого вектора t, ему ортогонального:

,

Подстановкой этого представления вектора в разность получаем другой ее вид

(4.4.9)

Анализируя выражение (4.4.9), замечаем, что первое слагаемое в нем не может быть изменено, поскольку

— фиксированный вектор при фиксированном k. Поэтому минимальному изменению аффинной модели будет отвечать случай, когда второе слагаемое в (4.4.9) будет нуль-вектором при Iвсяких векторах t, ортогональных векторам , т. е. следует находить из условия

. (4.4.10)

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

.

Таким образом, приходим к так называемой формуле пересчета С. Бройдена (1965 г.)

(4.4.11)

которая позволяет простыми вычислениями перейти от старой матрицы к новой такой, чтобы выполнялось соотношение секущих (4.4.4а) в новой точке и при этом изменения в аффинной модели (4.4.7) были минимальны

Совокупность формул (4.4.6), (4.4.11) вместе с обозначениями (4.4.5) называют методом секущих Бройдена или просто методом Бройдена решения систем нелинейных числовых уравнений.

Хотя в методах секущих обычным является задание двух начальных векторов ( и ), для метода Бройдена характерно другое начало итерационного процесса. Здесь нужно задать один начальный вектор , начальную матрицу и далее в цикле по k = 0,1,2. последовательно выполнять следующие операции:

1. решить линейную систему

(4.4.12)

относительно вектора :

2. найти векторы и :

, ; (4.4.13)

3. сделать проверку на останов (например, с помощью проверки на малость величин и/или и если нужная точность не достигнута, вычислить новую матрицу по формуле пересчета (см. (4.4.11))

(4.4.14)

В качестве матрицы , требуемой равенством (4.4.12) для запуска итерационного процесса Бройдена, чаще всего берут матрицу Якоби или какую-нибудь ее аппроксимацию. При этом получаемые далее пересчетом (4.4.14) матрицы , . не всегда можно считать близкими к соответствующим матрицам Якоби , . (что может иногда сыграть полезную роль при вырождении матриц ). Но, в то же время, показывается, что при определенных требованиях к матрицам Якоби матрицы обладают «свойством ограниченного ухудшения», означающим, что если и происходит увеличение с увеличением номера итерации k, то достаточно медленно. С помощью этого свойства доказываются утверждения о линейной сходимости () к х*

при достаточной близости к х* и к а в тех предположениях, при которых можно доказать квадратичную сходимость метода Ньютона (3.1.2), — о сверхлинейной сходимости последовательности приближений по методу Бройдена.

Как и в случаях применения других методов решения нелинейных систем, проверка выполнимости каких-то условий сходимости итерационного процесса Бройдена весьма затруднительна.

Формуле пересчета (4.4.14) в итерационном процессе можно придать чуть более простой вид.

Так как, в силу (4.4.12) и (4.4.13),

,

,

то из формулы (4.4.14) получаем формально эквивалентную ей формулу пересчета

, (4.4.14а)

которую можно использовать вместо (4.4.14) в совокупности с формулой (4.4.6) или с (4.4.12), (4.4.13) (без вычисления вектора ). Такое преобразование итерационного процесса Бройдена несколько сокращает объем вычислений (на одно матрично-векторное умножение на каждой итерации). Не следует, правда, забывать, что при замене формулы (4.4.14) формулой (4.4.14а) может измениться ситуация с вычислительной устойчивостью метода; к счастью, это случается здесь крайне редко, а именно, в тех случаях, когда для получения решения с нужной точностью требуется много итераций по методу Бройдена, т. е. когда и применять его не стоит.

2. Выявление наиболее оптимального метода

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

В таблицу было выписано количество шагов сходимости для каждого метода:

Квазиньютоновские методы, или когда вторых производных для Атоса слишком много

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

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

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

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

Методы секущих

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

Но что делать, если матрицу Якоби мы по каким-либо причинам вычислить не можем? Первое, что в данном случае приходит в голову — если мы не можем вычислить частные производные аналитически, то вполне можем получить для них численное приближение. Самым простым (хотя и далеко не единственным) вариантом такого приближения может служить формула правых конечных разностей: , где – j-й базисный вектор. Матрицу, составленную из таких приближений, будем обозначать . Анализу того, насколько замена на в методе Ньютона влияет его сходимость, посвящено достаточно большое количество работ, но нас в данном случае интересует другой аспект. А именно, такое приближение требует вычисления функции в N дополнительных точках, и, кроме того, функция в этих точках интерполирует функцию , т.е.

Не каждая аппроксимация матрицы Якоби обладает таким свойством, но каждая матрица аффинной функции, обладающей таким свойством, является аппроксимацией матрицы Якоби. Действительно, если и , то при . Данное свойство, а именно — свойство интерполирования, дает нам конструктивный способ обобщить метода Ньютона.

Пусть – функция, удовлетворяющая требованию для некоторой системы линейно независимых векторов . Тогда такая функция называется секущей функции , а определяющее ее уравнение — уравнением секущей. Если при этом система векторов полна (то есть их ровно N и они все еще линейно независимы), и, кроме того, система векторов линейно независима, то определяется однозначно.

Любой метод, построенный на основе локальной замены уравнения уравнением вида , где удовлетворяет уравнению секущих, называют методом секущих.

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

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

Вообще, самым устойчивым оказывается 2-х точечный метод секущих. То есть метод, в котором на каждой итерации нам придется вычислять дополнительно N-1 значений функции. Это явно не подходит для наших практических целей.

Тогда вопрос — а к чему все это было?

Квазиньютоновские методы решения уравнений

Выход из положения прост, хотя и не очевиден. Если у нас нет технической возможности на основании уже вычисленных значений однозначно определить аффинную модель, удовлетворяющую уравнению секущих — то и не надо. Возьмем уравнение секущих за основу, но будем требовать, чтобы оно выполнялось только для некоторой неполной системы векторов . Иными словами, мы будем требовать выполнение условия интерполяции только для достаточно малого количества известных значений. Разумеется, в этом случае мы уже не можем гарантировать, что используемая в такой модели матрица будет стремиться к матрице Якоби, однако это нам и не потребуется. Добавив к этому, что аффинная модель должна интерполировать функцию в текущей точке, то есть , мы получаем следующую формулировку метода секущих:

Бройден был первым, кто рассмотрел методы такого вида при m=1, назвав их квазиньютоновскими. Понятно, что условие секущих в данном случае позволяет однозначно идентифицировать матрицу только если наложить на нее дополнительные условия, и каждое такое дополнительное условие порождает отдельный метод. Сам Бройден рассуждал следующим образом:

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

Идея Бройдена гениальна в своей простоте. Действительно, если мы не имеем новой информации о поведении функции, то лучшее, что мы можем сделать — постараться не профукать старую. Тогда дополнительное условие

для всех , таких, что

позволяет однозначно определить матрицу нового преобразования — она получается добавлением к старой матрице поправки ранга 1.

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

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

Квазиньютоновские методы оптимизации

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

Для исправления указанного недостатка добавим дополнительное ограничение к бройденовской задаче минимизации, явно потребовав, чтобы новая матрица была симметричной наряду со старой:

Решением данной задачи является

Здесь , а формула пересчета матрицы носит имя своих создателей — Пауэлла, Шанно и Бройдена (PSB). Получаемая в результате матрица симметрична, но явно не является положительно определенной, если только вдруг не окажется коллинеарным . А мы видели, что положительная определенность является крайне желательной в методах оптимизации.

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

Происхождение такой постановки вопроса — отдельная большая тема, однако интересно, что если матрица T такова, что (то есть G также является матрицей аффинного преобразования, удовлетворяющего уравнению секущих для направления p), то решение данной задачи оказывается независимым от выбора T и приводит к формуле обновления

известной как формула Давидона-Флетчера-Пауэлла. Данный способ обновления неплохо зарекомендовал себя на практике, поскольку обладает следующим свойством:

если 0″» data-src=»https://habrastorage.org/getpro/habr/post_images/e3e/c44/1c1/e3ec441c17e1b43df108a7d8e15d3dd6.gif»/> и положительно определена, то также положительно определена.

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

Если в задаче, приводящей к DFP-методу, взять в качестве меры расхождения аффинных моделей расстояние не между самими матрицами, а между обратными к ним, то получим задачу вида

Ее решение — широко известная формула, открытая почти одновременно Бройденом, Флетчером, Гольдфарбом и Шанно (BFGS).

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

Все описанные методы обновления матрицы предполагают внесение поправки ранга 2. Это позволяет легко и непринужденно обратить матрицу , используя формулу Шермана-Моррисона и значение .

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

Осталось рассмотреть вопрос о том, как должна выглядеть самая первая матрица в квазиньютоновском процессе. Здесь все очевидно — чем она ближе к матрице Гессе или к ее скорректированному варианту, если гессиан вдруг не окажется положительно определенным, тем будет лучше с точки зрения сходимости. Однако в принципе нам может подойти любая положительно определенная матрица. Самый простой вариант такой матрицы — единичная, и тогда первая итерация совпадает с итерацией градиентного спуска. Флетчером и Пауэллом было показано (естественно, для DFP-метода), что в случае минимизации квадратичной функции в независимости от того, какая (положительно определенная) матрица будет использована в качестве начальной DFP-итерации приведут к решению ровно за N итераций, где N – размерность задачи, причем квазиньютоновская матрица совпадет с матрицей Гессе в точке минимума. В общем нелинейном случае такого счастья мы, само собой, не дождемся, однако это по крайней мере дает повод не слишком переживать на счет неудачного выбора начальной матрицы.

Заключение

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

Следует также отметить, что методы переменной метрики, несмотря на свое название, не всегда приводят к построению матриц, которые собственно метриками являются, хотя и делают это каждый раз, когда это вообще возможно. Это обычно не составляет большой проблемы, однако желающие обезопасить себя от возможного конфуза вполне могут прибегнуть к тем же ухищрениям, которые делались для преодоления аналогичной проблемы у метода Ньютона — например, путем смены направления или применения схемы Левенберга-Марквардта. Правда, в этом случае снова приобретут актуальность вопросы выбора формы доверительного региона, но тут уж придется выбирать меньшее из зол. Другое решение проблемы заключается в использовании методов линейного поиска, гарантирующих выполнение необходимого условия сохранения положительной определенности. Правило Вульфа гарантирует выполнение данного условия, в то время как правила Армиджо и Гольдштейна — нет.

Теоретически определить, какой из огромного количества возможных квазиньютоновских методов окажется наиболее эффективным по отношению к некоторому классу задач практически невозможно. Обычно при формулировании метода ограничиваются тем, что показывают его эффективность при минимизации квадратичной функции (к слову, метод считается эффективным, если приводит к точному решению за N итераций, то есть не медленнее, чем прямые методы решения СЛАУ). В более редких случаях можно найти исследования порядка сходимости метода (который обычно сверхлинейный, то есть существенно лучше того, что нам дает градиентный спуск), устойчивости и других представляющих интерес характеристик. Но в целом, единственный разумный критерий, позволяющий судить об эффективности того или иного метода для конкретного класса задач — практика. Так что лопаты в руки — и успехов в применении.