Главная  Методы условной оптимизации 

[0] [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [ 31 ] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52] [53] [54] [55] [56] [57] [58] [59] [60] [61] [62] [63] [64] [65] [66] [67] [68] [69] [70] [71] [72] [73] [74] [75] [76] [77] [78] [79] [80] [81] [82] [83]

которые известны заранее и определяются природой минимизируемой функции.

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

4.8.1. ДИСКРЕТНЫЕ МЕТОДЫ НЬЮТОНА ДЛЯ ФУНКЦИЙ С РАЗРЕЖЕННЫМИ МАТРИЦА;ЧИ ГЕССЕ

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

При построении конечно-разностной аппроксимации плотной матрицы Гессе в качестве конечно-разностных векторов используются столбцы единичной матрицы, и расчет каждого столбца искомого приближения требует вычисления 1радиента в одной дополнительной точке (см. разд. 4.5.1). Если же матрица Гессе С(х) разрежена и структура ее заполнения известна, то, подобрав конечно-разностные векторы специальным образом, полную аппроксимацию для G(x) можно построить ценой вычисления градиентов в значительно меньшем чем п числе точек.

Пусть, например, G(x) является тридиагональной матрицей вида

и подсчитаны векторы 1

"л ii г,-и

(8(4 + 1!)-ВЫ), i=l,2.

где г, = (1,0, 1,0, ...), г, = (0, 1,0, 1, ...), а й-иекотпрый подходящий конечно-разностный интервал. Легко проверить, что й будет приближением суммы нечетных столбцов G(х), а у,- суммы четных. Таким образом, обозначив через у, ,- компоненту вектора у, с номером i и введя аналогичные обозначения для компонент у, получим

г/1.1:

dxi дХй

dxidxg дхдхц

и т. д. Отсюда, в частности, вндно, что

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

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

в силу симметрии для построения ее конечно-разностной аппроксимации G также достаточно двух дополнительных вычислений градиента Здесь в роли конечно-разностных векторов надо использовать г, = (1, 1. 1, 0) и г, = (0, О.....О, 1)г. Приращение градиента вдоль Zj даст первые п-1 диагональных элементов матрицы С; оставшиеся ненули G определятся по конечным разностям вдоль Zj.

Процедуры решения задач минимизации функций большого числа переменных с ра:1реженными матрицами Гессе, построенные на основе дискретного метода Ньютона, обычно работают в два этапа: сначала запускается алгоритм предварительной обработки, который, перенумеровывая переменные задачи, преобразует исходную структуру заполнения матрицы Гессе к виду, позволяющему уменьшить количество вычислений градиента при оценке вторых производных; затем начинается собственно минимизация. Поскольку сокращение количества вычислений градиента для расчета G,, проявится столько раз, сколько будет выполнено итераций спуска, экономить на предварительной обработке задачи и поиска подходящего набора конечно-разностных векторов не следует.

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

G,p, = -g. (4.73)

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



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

Итак, разреженность системы (4.73) еще не гарантирует возможности ее эффективного решения. Например, в модифицированной факторизации Холесского (см. разд. 4.4.2.2) столбцы треугольного фактора будут линейными комбинациями столбцов Сд, и, следовательно, матрица может оказаться довольно плотной даже при очень разреженной Gj. В данном случае воспользоваться факторизацией Холесского не удастся. Чтобы заполнение в процессе какой-то факторизации матрицы Gj было приемлемым, она должна иметь определенную структуру.

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

*4.8.2. КВАЗИНЬЮТОНОВСКИЕ МЕТОДЫ ДЛЯ ФУНКЦИЙ С РАЗРЕЖЕННЫМИ МАТРИЦАМИ ГЕССЕ

В разд. 4.5.2.3. говорилось о том, что некоторые кваэнньюта-новские формулы для расчета приближенной матрицы Гессе могут быть получены из решения задачи минимизации нормы уклонения последующего приближения от предыдущего при определенных ограничениях, среди которых всегда содержится квазиньютоиовское условие (4.30). В частности, была представлена задача (4.48), решение которой дает PSB-формулу. В этом разделе мы обсудим возможности применения оптимизационного подхода для вывода формул квазнньютоновской аппроксимации, обеспечиваюищх приближениям матрицы Гессе ту же структуру заполненности, что и у нее самой.

Когда матрица Гессе минимизируемой функции разрежена и структура ее заполненности известна, формулы расчета квазиньютоновских приближений, сохраняющих эту структуру, можно получить, введя Б задачу (4.48) требование равенства нулю соответ-

ствующнх элементов поправочной матрицы. Если обозначить через набор индексов {(i, /)0(>=0}, то обобщением (4.48) в указанном смысле будет задача вида

найти minl/(.= S

(4.74)

при ограничениях t/sj = уд-Bs,,,

Ее решением является матрица U, гарантирующая наличие у Bj+j той же структуры заполненности, что и у матрицы Гессе, если только этой структурой обладает Вд. Явные формулы расчета U удобно записываются через векторы которые определены

следующим правилом: чтобы получить вектор о*, надо взять вектор s* и обнулить все его компоненты, соответствующие нулевым позициям j-ro столбца матрицы Од. Довольно длинные в сложные выкладки показывают, что

(4.75)

Здесь ej есть j-u единичный вектор, а числа К/ представляют собой множители Лагранжа, связанные с оптимумом задачи (4.74). Составленный из них вектор Я является решением линейной системы уравнений

Q}.=y„--Bs„. (4.76)

= S(o;V/--0i>iie)e[.

Матрица Q симметрична и имеет ту же структуру заполненности, что и Вд; она положительно определена в том и только том случае, когда о">0 для всех /. Заметим, что ранг поправочной матрицы t/j (4.75), вообще говоря, равен п и что линейную систему (4.76) придется на каждом шаге решать «с нуля». Следует также иметь в виду, что сохранение положительной определенности квазиньютоновских приближений при использовании поправок (4.75) не гарантировано.

На каждой итерации рассматриваемой квазииьютоновской схемы кроме системы (4.76) придется решать систему (4.28) и тоже «с нуля», так как поправка (4.75) имеет высокий ранг (в отличие ог поправок в обычных квазиньютоновских методах). Таким образом, итерация схемы требует решения двух больших систем линейных уравнений с разреженными матрицами, причем для записи матрицы Q вспомогательной системы (4.76) в памяти машины придется выделить специальное место.



Исходя из полученного явного вида решения задачи (4.74), можно предложить универсальный прием модификации квазиньютов-скнх методов, позволяющий любой из них преобразовать в метод с сохранением структуры заполненности. Для того чтобы построить разреженную поправочную матрицу, ориентируясь на какую-нибудь обычную квазиньюггоновскую формулу, надо (i) вычислить матрицу и, структура заполненности которой та же, что и у матрицы Гессе, а ненулевые элементы определены по обычной формуле;

(ii) решть линейную систему уравнений

QX=Hj-BftSft-fs;

(iii) вычислить для найденного X матрицу (4.75) и сложить ее с (7 (результат и есть искомая разреженная поправка). В данном случае (4.75) будет решением задачи минимизации нормы корректирующей добавки к обычной квазиньютоновской матрице U, обнуляющей те элементы последней, которым отвечают нули матрицы Гессе.

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

4.8.3. МЕТОДЫ СОПРЯЖЕННЫХ ГРАДИЕНТОВ

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

4.8.3.1. Квадратичные функции. Рассмотрим задачу поиска минимума квадратичной функции:

Ф{х) = с-х+~х-Сх,

(4.77)

где G- симметричная положительно определенная матрица. Пусть Хд-текущее приближение и известны ft-f 1 линейно независимых

векторов Ро, pi.....Pi, порождающих некоторое подпространство

Sl,. Найдем выражение для точки минимума Ф на многообразии Xt+Sh- Для этого введем матрицу Рд, столбцами которой явля-

ются {pj}. Это позволит представить задачу минимизации Ф на многообразии Xi+Sk в виде

найти rain Ф(хд-(-РдШ).

в соответствии с определение.м Ф оптимальное значение аргумента W будет точкой минимума квадратичной функции

WPtgk + WPlCP,}. (4.78)

где g=VФ(xд)=c+Gxд. Эта точка вычисляется по формуле

ш=-(РдГ0Рд)-Р1гд.

и, следовательно, минимум Ф на многообразии Хн+Зн бУДет достигаться в

x„+t = Ч-PЛPlGPrPIek (4-79)

На равенство (4.79) можно смотреть как на формулу задания пошаговой процедуры минимизации Ф. Эта процедура имеет несколько интересных свойств. Прежде всего легко убедиться, что градиент gk+i функции Ф в точке Хд+х будет ортогонален векторам {pi) (столбцам матрицы Рд). В самом деле,

Plgk+t = PI (c + Gx,,) = Plg,-Pl GP, {PI GP,)- Pjg*=0,

T. e. gLiP,= 0. i = 0, k. Поскольку мы считаем, что точки Xj, /=1, k, тоже определены мини.мизацией Ф на соответствующих многообразиях, полученные равенства позволяют утверждать, что

gfp, = 0, i>i, (4.80)

для всех 1=1, . . ,k, k+l. Отсюда в свою очередь следует, что (4.79) можно переписать так:

x,+ix,+yPk{PlGP,)--e„ (4.81)

где Сд есть fe-й столбец единичной матрицы, а y = - glpi,.

Формула (4.81) сильно упростится, если матрица PlGP,, диагональна. Это будет, если векторы {pj\, 1 = 0, ...,k, взаимно сопряжены относительно G, т. е. если

prGpi = 0, 1ф;, i, ; = 0, 1.....к. (4.82)

в данном случае (4.81) преобразуется к виду

Хд« = х,-Ьадрд, (4.83)

где <4 = - gkPi,IPkGpi,. Легко проверить, что ад есть не что иное, как длина шага в точку минимума Ф вдоль рд. Таким образом, процедура последовательной минимизации Ф на многообразиях, определяемых сопряженными векторами, принимает стандартную форму метода спуска.



[0] [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [ 31 ] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52] [53] [54] [55] [56] [57] [58] [59] [60] [61] [62] [63] [64] [65] [66] [67] [68] [69] [70] [71] [72] [73] [74] [75] [76] [77] [78] [79] [80] [81] [82] [83]

0.001