Вычислительная математика лекции

x(n+1) = x(n) + αn p(n) .

В случае квадратичной целевой функции её минимум в направлении p(n)

достигается при величине шага

При отсутствии вычислительных ошибок требуется не более n итераций.

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

F (x) 12 ( Ax, x) (b, x) эквивалентна решению системы линейных алгебраических уравнений Ax = b.

Метод сопряженных градиентов в общем случае.

В общем случае для расчета βn-1 можно использовать одну из следующих формул

Эти формулы не содержат явным образом матрицу А.

Минимизация идет в следующей последовательности

p(0)= -g(0), p(n) = -g(n) + βn-1p(n-1), n 1

n (n ) minn ( ),n ( ) f (x(n) p(n) ), 0.

x(n 1) x(n) n pn ,n 0.

При расчете βn-1 по третьей формуле приходится вычислять матрицу Гессе g'(n). Однако в этом случае возможно использование квазиньютоновских методов для аппроксимации этой матрицы.

Реализация такой процедуры использована в популярном BFGS –

алгоритме (Broyden, Fletcher, Goldfarb, Shanno).

191

Теперь получить точное решение за конечное число шагов невозможно даже при отсутствии вычислительных ошибок.

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

15.2.8. Метод Ньютона.

15.2.8.1. Простейший вариант метода Ньютона.

Пусть в заданной области ||x-x*||<ϭ выполнены необходимое и достаточное условия существования минимума функции f(x). Тогда задача

вычисления x*=arg min f(x) сводится к нахождении корней системы

x

нелинейных уравнений g(x)=0, где g(x)= grad(f(x)). Для решения задачи можно использовать метод Ньютона g(xn)+G(xn)(xn+1-xn)=0. Здесь xn n – ое приближение, G(xn) матрица Гессе.

Таким образом, чтобы попасть из точки xn в точку xn+1, нужно переместится вдоль вектора pn=xn+1-xn, который определяется из системы линейных алгебраических уравнений G(xn)pn=-g(xn).

Вектор pn называют ньютоновским направлением, а метод спуска xn+1=xn+αnpn – методом Ньютона. Ньютоновское направление

действительно является направлением спуска. В самом деле, pn=-G-1(xn)g(xn). Матрица G-1(xn) положительно определена (это следует из положительной определенности матрицы G(xn)). Поэтому

(g(xn),pn)=-(G-1(xn)g(xn),g(xn)<0.

При выборе αn=1, метод в точности совпадает с методом Ньютона решения системы нелинейных уравнений g(x)=0.

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

192

Метод Ньютона может быть интерпретирован как последовательная квадратичная аппроксимация. Действительно, при выполнении условий непрерывности любую функцию f(x) можно аппроксимировать в окрестности точки x0 функцией

(x) f (x0 ) (x x0 ) f (x0 ) 12 (x x0 )T G(x0 )(x x0 ) , где

G(x0) – матрица Гессе, вычисленная в точке x0.

Разумной аппроксимацией минимума функции f(x) может быть минимум функции φ(x). Если последний находится в точке xm, то

'(xm ) 0 или f (x0 ) G(x0 )(xm x0 ) 0,

откуда

xm x0 -G 1 (x0 ) f (x0 ),

или

xm x0 -G 1 (x0 )g(x0 ).

Полагая x0=xn, xm=xn+1, получаем формулу метода Ньютона. xn 1 xn -G 1(xn )g(xn ).

15.2.8.2. Метод Ньютона с дроблением шага.

Позволяет облегчить проблему выбора начального приближения,

однако не обеспечивает достаточную скорость вычисления.

Аналогично методу наискорейшего спуска очередной шаг αn выбирается из условия f (xn 1) f (xn ) (gn , pn ), 0< <1.

Теорема о сходимости.

Пусть выполнены следующие условия:

1.f(x) ограничена снизу.

2 Матрица Гессе удовлетворяет условиям

|| y ||2 (G(x) y, y) R || y ||2 ; , R 0

xn+1 = xn+αnpn .

193

Тогда ||gn||→0 при n→∞. Целевая функция f(x) монотонно уменьшается с ростом n.

В процессе доказательства теоремы формулируются требования к выбору коэффициента 2 2 (1 ).

Алгоритм метода, учитывающий эти условия, выглядит следующим образом.

1. Выбираем произвольную точку x0 – начальное приближение.

Полагаем xn=x0.

2.Выбираем произвольное значение α.

3.Определяем новую точку xn+1=xn+αpn=xn-αG-1(xn)g(xn). .

4.Вычисляем f(xn+1)=f(xn+αpn) и f(xn).

5.Если f(xn+1)-f(xn) γα(g(xn), pn ), 0<γ<1, то значение α берем в качестве искомого αn и переходим к шагу 6, иначе уменьшаем α ( путем

умножения на множитель γ, 0<γ<1) и переходим к шагу 4.

6.Вычисляем новое приближение xn+1=xn-αG-1(xn)g(xn). .

7.Проверяем выполнение условия стационарности

||gn+1||= (gn+1,gn+1) ε – заданная погрешность. Если условие выполнено переходим к шагу 8, иначе к шагу 2, положив xn=xn+1.

8. Останов.

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

15.2.8.3.Понятие о квазиньютоновских методах.

Вметодах, которые называются квазиньютоновскими, также как и в методах Ньютона используется квадратичная аппроксимация целевой функции. Однако, здесь не производится вычисление матрицы Гессе.

Информация о кривизнах функции здесь накапливается постепенно в

результате наблюдений за изменением градиента.

194

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

B(n) p(n) g(n) ,

в которой B(n) - текущее приближение к матрице Гессе. Начальное приближение обычно полагают равным единичной матрице. В таком случае p(0) = -g(0) и первая итерация совпадает с одним шагом метода наискорейшего спуска.

После того, как найдено приближение x(n+1) вычисляют матрицу

B(n+1) = B(n) + ΔB(n),

где ΔB(n) некоторая поправочная матрица. В пределе матрица B(n)

становится матрицей Гессе.

Среди наиболее популярных квазиньютоновских алгоритмов известен

BFGS – алгоритм (Broyden, Fletcher, Goldfarb, Shanno), в котором аппроксимация B(n 1) производится итерационно по формуле

где

s(n) = x(n+1) – x(n), q(n) = g(n+1) - g(n)

В некоторых вариантах квазиньютоновских методов аппроксимируется матрица H (n) , обратная B(n) (метод Давидона - Флетчера – Пауэла).

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

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

так и свойство сопряженных направлений.

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

1. Гилл Ф., Мюррей У., Райт М. Практическая оптимизация. – М.; Мир, 1985. 195

2. Банди Б. Методы оптимизации. Вводный курс. – М. : Радио и

связь. 198

15.2.9. Метод Левенберга – Марквардта.

Рассматриваемый метод является модификацией метода Ньютона,

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

Направление поиска Левенберга — Марквардта определяется из системы: (G(xn)+λnE) pn=-g(xn).

где λn — некоторая неотрицательная константа, своя для каждого шага, E — единичная матрица, xn+1=xn + pn.

Нетрудно заметить, что при λn = 0 алгоритм вырождается в метод Ньютона. При достаточно большом λn метод незначительно отличается от метода наискорейшего спуска, сходимость которого гарантирована при произвольном выборе начального приближения. Таким образом, при правильном подборе параметра λn можно добиться сочетания достоинств методов наискорейшего спуска и Ньютона . Величину λn можно увеличивать, пока выполняется условие f(xn+1) < f(xn) . Однако при чрезмерно большом λn проявляются все недостатки метода градиентного спуска, связанные с проблемой овражности. Дополнительное решение этой проблемы предложил Марквардт. Он заметил, что если заменить единичную матрицу на диагональ матрицы Гессе, то можно достичь увеличения шага вдоль пологих участков и уменьшения вдоль крутых спусков:

(G(xn)+λn diag(G(xn)) pn=-g(xn).

196

16.Методы решения краевой задачи для линейных

обыкновенных дифференциальных уравнений.

16.1. Постановка задачи.

Ранее уже рассматривались методы решения обыкновенных

дифференциальных уравнений – методы решения задачи Коши:

(t) = f(x,t), x(t0) = x0 , t [t0,tk] , x(t) = ?

Краевая задача отличается от данной заданием краевых вместо

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

(t) = f(x,t), t [a,b] ,

αTx(a) = γa ,

βTx(b) = γb ,

где α, β – заданные числовые векторы одной размерности с x, γa , γb

– заданные числа. Граничные условия могут быть заданы и в более сложной форме.

В отличие от имеющей единственное решение задачи Коши краевая

задача может иметь или одно решение, или бесконечное множество решений, или, наконец, может совсем не иметь решений.

Рассмотрим наиболее часто встречающуюся на практике краевую

Теорема. Для того, чтобы существовало единственное решение

неоднородной краевой задачи (1), (2) необходимо и достаточно, чтобы

197

≡ ,t [a,b] имела только

Рассматриваемый метод весьма трудоемок, так как требует

многократного интегрирования системы (3).

16.3. Разностный метод (метод конечных разностей).

заданные числа.

Теорема. При выполнении указанных условий задача (1), (2) имеет единственное решение x(t) C4[0,1] .

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

постоянным шагом {tk}Nk=0 , где tk =kh, h=1/N, N≥2. Значения сеточной функции в узлах будем обозначать xk .

Аппроксимируем значения производных соответствующими разностными операторами

x '(tk ) xk 1 xk 1 , 2h

x ''(tk ) xk 1 2xk xk 1 . h2

Тогда вместо дифференциальной задачи (1), (2) получаем разностную

Разностное уравнение (2) при всех значениях k эквивалентно системе линейных алгебраических уравнений

199

b1x0 + c1x1 +d1x2 = f1 , b2x1 + c2x2 +d2x3 = f2 ,

………………

bN-1xN-2 + cN-1xN-1 +dN-1xN = fN-1 .

Если дополнить эти уравнения сверху соотношением x0 = x0 и снизу xN = x1, то получим линейную алгебраическую систему

Ax = r ,

где x = (x0 ,x1 ,…,xN )T, r = (x0 ,f1 ,f2 ,…,fN-1, x1 )T и А – трёхдиагональная матрица

Если выполняется условие |ck|≥|bk| + |dk|, то говорят, что матрица А является матрицей с доминирующей главной диагональю. Указанное

условие гарантирет существование и единственность решения. Оно всегда выполняется, если |pk|h 2.

Наиболее экономично решать задачу с трехдиагональной матрицей

методом прогонки.

17.Сингулярное разложение и метод наименьших квадратов.

17.1. Решение переопределенной системы уравнений.

Вернемся к задаче среднеквадратичного приближения сеточной функции. Заданную табличную функцию {fj , xj } размера m (j=0,1,…,m-1)

требуется заменить обобщенным полиномом степени (n-1), имеющим набор базисных функций φi (x), i=0,1,…,n-1.

200