Метод проекции градиента
Назначение сервиса. В онлайн-калькуляторе для нахождения условного экстремума функции используется основной алгоритм проекции градиента (случай линейных ограничений).Рассмотрим задачу оптимизации при единственном линейном ограничении в виде равенства
f(x) → min, (1)
. (2)
В заданной точке xk, в которой ▽f(xk)≠0, делается попытка найти направление поиска, которое бы лежало на поверхности ограничения и являлось направлением спуска. Такое направление можно получить геометрически ортогонально проектируя вектор, противоположный ▽f(xk) на поверхность ограничения (см. рис.1).
Здесь ▽fck - проекция антиградиента на поверхность ограничений, которая приводит в допустимые точки. Действительно, для ∀α ≥ 0 точки, заданные соотношением
x = xk-α▽fck (3)
удовлетворяют линейному ограничению
.
Это направление задает спуск, т. к. угол между ▽fk и ▽fck больше 90°. Процесс ортогонального проектирования состоит в разложении вектора на две ортогональные компоненты: параллельную поверхности, заданной ограничением, и перпендикулярную к ней. Параллельная компонента является искомой проекцией градиента.
Пусть вектор a - нормаль к поверхности ограничения. Отметим, что из выражения aTS=0 следует допустимость направления, задаваемого вектором S (S параллельно поверхности).
Т. о. все векторы, перпендикулярные к поверхности ограничения должны быть параллельны к a. Следовательно, для любого вектора S его компонента S’, перпендикулярная к поверхности ограничения, равняется значению a, умноженному на константу.
Обозначим через S” компоненту S, параллельную поверхности ограничения. Тогда S” удовлетворяет соотношению
aTS”=0. (4)
Т.о. любой вектор можно представить в виде векторной суммы:
S=S′+S″, (5)
где S′=λa, а S”удовлетворяет уравнению aTS”=0.
Найдем λ. Рассмотрим скалярное произведение aTS. В силу (5) и (4) имеем
aTS=aTλa+aTS″= λaTa, (6)
откуда λ=(aTa)-1aTS. (7)
Из (5) найдем S″=S-S′=S- λa. Подставим сюда (7), получим
S″=S-a(aTa)-1aTS=(I- a(aTa)-1aT)S , (8)
где I- единичная матрица, порядок которой согласован с S. Матрица P=I- a(aTa)-1aT - проекционная матрица. Она проектирует вектор S на плоскость, задаваемую ограничением h(x).
Отметим, что P является симметрической и положительно полуопределенной. Симметричность P очевидна. Для доказательства положительной полуопределенности рассмотрим произведение yTPy для произвольного y≠0. Тогда
.
Используя неравенство Шварца (yTa)2≤(yTy) (aTa), убеждаемся, что числитель неотрицателен.
Свойства проекций
Пусть S”=-P▽f (xk) - направление спуска.Если S”=0, то точка xk удовлетворяет необходимым условиям Лагранжа. Вектор множителей Лагранжа задается выражением
λ=(AAT)-1A▽f. (9)
Первое утверждение следует из того, что матрица P симметрическая и положительно полуопределена. Действительно, имеем
xk+1=xk-αP▽f (xk); α≥0;
f(xk+1)=f(xk- αP▽f (xk))= f(xk)- α▽fT(xk+1)P▽f (xk)
Так как
, то fk+1≤fk, что и требовалось доказать.
Имеем S”=-P▽f, т. е. S“ - это проекция вектора ▽f на поверхность ограничений. Поэтому, если S”=0, значит градиент ▽f перпендикулярен поверхности ограничений.
Рассмотрим второе свойство. Т.к. S”=-P▽f, то произведение ▽fTS”=-▽fTP▽f≤0. Если S”=-P▽f=0, то следовательно ▽f∟S”, т. е. ▽f перпендикулярен поверхности ограничений. Тогда из формулы
S=S’+S”=ATλ+S”=0 следует
▽f =ATλ. (10)
Т.к. строками матрицы A являются векторы коэффициентов в линейных ограничениях, то (10) представляет собой другую форму записи необходимого условия Лагранжа:
, (11)
где ak - k-й вектор ограничений. Отсюда мы видим, что λk- это множители Лагранжа, т. е. множители функции
.
Необходимое условие экстремума:
;
▽hk=ak.
Нетрудно получить выражение для λ: S=S’+S”;
где (AAT)-1 существует только в том случае, если значения ak линейно независимы. Если имеются линейно зависимые ограничения, их следует исключить из рассмотрения.
Основной алгоритм проекций градиента
Вычислить матрицу P=I-AT(AAT)-1A в предположении, что векторы ak линейно независимы. Задать ε>0- погрешность сходимости. Пусть найдена допустимая точка xk.Шаг 1. Вычислить Sk=-P▽f (xk).
Шаг 2. Если |Sk|≤ ε, то вычислить λ по формуле (9) и закончить вычисления. В противном случае - продолжить вычисления. Переход на Шаг 3.
Шаг 3. Определить максимальную длину шага:
f(xk+αSk)→min; 0≤α≤ αmax.
Шаг 5. Положить xk+1=xk+α Sk → Шаг 1.
Замечание. Мы рассмотрели алгоритм, в котором используются линейные равенства. Но его можно легко распространить на неравенства, используя либо дополнительные переменные, либо активные ограничения, что предпочтительнее, т. к. второй способ позволяет уменьшить размерность.
При использовании второго способа с помощью получаемой по формуле (12) оценки множителей Лагранжа осуществляется поочередное исключение ограничений из множества активных ограничений. Модификация выглядит следующим образом. В заданной точке xk для определения активного множества проверяются ограничения в виде неравенств: ajTx≥bj, j=1..m.
Составляется матрица ограничений из строк, соответствующих активным ограничениям. Вычисляется проекционный оператор P и проекция Sk. Если |Sk|≤ ε, то вычисляем множители Лагранжа:
λ=(AAT)-1A▽f (xk).
Если все λi≥0, то решение найдено. В противном случае ограничение с наибольшим по модулю множителем Лагранжа исключается из множества активных ограничений, вычисляется заново P и делается переход на Шаг 1.
Пример №1. f(x)=2*x12+2*x22-2*x1*x2-4*x1-6*x2
X0=(0;0)
g1: x1+x2≤2
g2: x1+5*x2≤5
g3: x1 ≥ 0
g4: x2 ≥ 0
Определим градиент целевой функции:
| ▽ f(X) = |
|
Тогда A0:
| A0= |
|
| ▽ f(X0) = |
|
| A0= |
|
| P=I-AT(AAT)-1A = |
| - |
| = |
|
| S0=-P ▽ f(X0) = - |
|
| = |
|
Следовательно, необходимо проверить соответствующие ограничения и множители:
| λ = (AAT)-1A ▽ f(X0) = |
|
A0=(1;0)
Шаг 1.
В точке x0 нет активных ограничений. P=I
| S0=-P ▽ f(X1) = - |
|
| = |
|
Шаг 3.
λmin = (0.2; 0.147; 0; 0) = 0.147
Шаг 4. Ищем минимум вдоль прямой:
| X = X1 + hS1 = |
| + h |
| = |
|
Найдем такой шаг h, чтобы целевая функция достигала минимума вдоль этого направления. Из необходимого условия существования экстремума функции (f'(x1)=0):
112.0*h-52.0=0
0 ≤ h ≤ 0.147
h=0.464, лежит вне интервала.
Шаг 5. X1 = (0.588;0.882)
Итерация 1.
| ▽ f(X2) = |
|
| A1= |
|
| P=I-AT(AAT)-1A = |
| - |
| = |
|
| S1=-P ▽ f(X2) = - |
|
| = |
|
Шаг 3.
g2: λmax = max(0, ∞)
λmin = (0.969; ∞; 0; 0) = 0.969
Шаг 4. Ищем минимум вдоль прямой:
| X = X2 + hS2 = |
| + h |
| = |
|
Найдем такой шаг h, чтобы целевая функция достигала минимума вдоль этого направления. Из необходимого условия существования экстремума функции (f'(x2)=0):
32.995*h-6.9183=0
0 ≤ h ≤ 0.969
h=0.21; x=(1.129;0.774)
Шаг 5. X2 = (1.129;0.774)
Итерация 2.
| ▽ f(X3) = |
|
| A2= |
|
| P=I-AT(AAT)-1A = |
| - |
| = |
|
| S2=-P ▽ f(X3) = - |
|
| = |
|
Следовательно, необходимо проверить соответствующие ограничения и множители:
λ = (AAT)-1A ▽f(X3) = 1.032
λ≥0. Таким образом, точка x2 удовлетворяет условиям Куна-Таккера.
Ответ: x = (1.129;0.774)
Пример №2. f(x)=(x1-1)2+(x2-2)2→min
g1=x1-2x2≥-2
g2=-x1-x2≥-4
g3=x1≥0
g4=x2≥0
Рис. 3
В точке x0 ограничения g1 и g2 неактивные, а ограничения g3=0 и g4=0, т. е. активные. Поэтому a3=(1,0); a4=(0,1).
Имеем 
; AAT=I, P=(I-I) =0.
Шаг 1. 
. Следовательно, необходимо проверить соответствующие ограничения и множители
.
Так как λ2 - наибольший по модулю множитель, второе активное ограничение g4≥0 исключается. При этом активное множество сокращается до одного единственного ограничения g3≥0.
Таким образом, A=(1,0); AAT=1 - скаляр.
Следовательно 
.
Шаг 1. 
.
Шаг 2. Т. к. |S0|=4>ε, то переход на Шаг 3.
.
, 
;
.
Шаг 5. x1=(0;1).
В этой новой точке опять проверяем ограничения на активность. Ясно, что ограничения g1≥0 и g3≥0 активные. Таким образом
, 
.
.
Вектор ▽f(x1)=(-2;-2);
=(1;-3).
Т. к. λ2<0, то второе активное ограничение g3≥0 исключаем из множества активных ограничений. Теперь активное множество состоит только из g1. Следовательно A=(1,-2);
.
Шаг 1.
.
Шаг 2. |S1|>0, следовательно, продолжаем.
Шаг 3. Длина шага λmax=5/6 определяется ограничением g2≥0.
Шаг 4. Ищем минимум функции f(x) вдоль прямой:
;
f(x)=(2.4α-1)2+(1+1.2α-2)2, откуда α=1/2.
Шаг 5. x2=(1.2;1.6). В этой точке активно только первое ограничение. Множество активных ограничений остается без изменений и следовательно нет необходимости пересчитывать P. В точке x2 имеем:
;
.
