Главная Методы условной оптимизации [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] жение Холесского матрицы Сд и в качестве Сд берется произведение LDL, где d,-= max{d,, 6}. Однако данный подход страдает двумя серьезными недостатками. Во-первых, для знако-неопределенной матрицы Сд факторизации Холесского может просто не существовать (см. соответствующий пример в разд. 2.2.5.2). Во-вторых, даже если она существует, процедура расчета факторов, вообще говоря, будет численно неустойчивой, поскольку их элементы могут быть сколь угодно большими независимо от обусловленности Сд. Кроме того, матрица Сд в данном случае может очень сильно отличаться от Сд, даже если Сд всего лишь «слегка знаконеопределена». Пример 4.7. Возьмем матрицу fi 1 Сь = I 1-1-10-" V2 3 Ij Ее собственные значения приблизительно равны 5.1131, -2.2019 и 0.0888 (слагаемое \0~" введено, чтобы можно было построить разложение Холесского; при вычислении собственных чисел в рамках принятой точности оно не сказывается). Точные факторы Холесского для матрицы из примера 4.7 выглядят так: (1 о о\ " " \ 1 I о ) , в = I о 10-» о 2 10» ij \0 О -{З + Ю")/ Поэго.чу при рассматриваемом подходе к построению Сд норма ЦСд-Сд[ будет величиной порядка 10». В то же время второй, более радикальный из методов предыдущего раздела дает матрицу Сд, для которой Сд-Сд« 4.4038. Итак, на основе обычной факторизации Холесского хорошего метода построения Сд не получить. Для этого нужно воспользоваться модифицированной факторизацией. Данный подход состоит в том. чтобы строить матрицы Холесского L н D, подчиняющиеся двум требованиям: все диагональные элементы D должны быть существенно положительными; модули всех элементов треугольного фактора L должны быть равномерно ограничены сверху. Точнее говоря, требуется, чтобы для всех А = 1, 2, .... п и некоторых заданных положительных 6 и р выполнялись неравенства б>6. к,* К Р. l>k, (4.22) где г,.д-введенные для удобства наложения вспомогательные величины, по определению равные /,-дКйд (выбор Р обсуждается ниже). При этом процедура расчета модифицированных L и D фактически представляет собой обычный алгоритм факторизации Холесского с попутным увеличением (по мере необходимости) диагонали исходной матрицы с целью добиться неравенств (4.22). Рассмотрим /-Й шаг предлагаемой процедуры. Пусть первые /-1 столбцов модифицированных факторов Холесского известны и для А=1, . . ., /-1 условия (4.22) выполнены. Вычислим величину Т/= ?,-2йл. (4.23) где через 5, обозначен диагональный элемент gjj матрицы Сд. В качестве пробного значения d для /-го диагонального элемента D возьмем d = raax{Yy, 6}, где 6-малое положительное число из (4.22). Если окажется, что при dj = d значения г, найденные в соответствии с формулой (4.21Ь). удовлетворяют (4.22), зафиксируем полученные dj и /-Й столбец L и перейдем к следуюо1ему шагу. Если же неравенства (4.22) при djd нарушаются, т. е. какие-то из r,j оказываются больше, чем р, окончате.пьное значение dj определим по формуле, получающейся из (4.23) заменой у на gjj + ejj, где число ejj подбирается так, чтобы максимальная из величин г,у совпала с р. Матрицы L и D, полученные по окончании описанной процедуры, будут факторами Холесского для положите.г1ьно определенной матрицы Сд, связанной с Сд следующим образом: Здесь £-неотрицательная диагональная матрица, /-й элемент которой равен еуу. Таким образом, положите.льно определенная матрица Сд может отличаться от исходной матрицы Гессе Сд только диагональными элементами. При заданной матрице Сд диагональная поправка Е с очевидностью зависит от величины р. Желательно, чтобы эта величина была достаточно большой, так как заниженное значение р приведет к неоправданной модификации. Для положительно определенной матрицы Сд иа (4.21а) следует, что при любом i=l, . . ., п и каждом /(/<() справедливо неравенство Ildign. Значит, чтобы обеспечить равенстю нулю поправки Е, если матрица G„ существенно положительно определена, в качестве р надо взять величину, удовлетворяющую неравенству Р>у. где у - значение максимального из диагональных элементов Сд. Когда у матрицы Сд есть отрицательные собственные значения, поправка Е будет ненулевой независимо от выбора р. При этом для п>1 справедлива оценка »£ (Р) 1. < (I-t-(«-1) Р)Ч2 (Y-I- (п-1) Р) -I-6. Здесь I - максимальный модуль недиагоиального элемента С. Правая часть неравенства достигает минимума при р = )/п-1 и неограниченно возрастает с увеличением р. Это говорит о том, что в общем случае слишком большие значения Р столь же нежелательны, как и слишком малые. Помимо того что завышенное значениеР чревато потерей численной устойчивости процедуры построения С,„ оно может еще и привести к неоправданно большому уклонению Gf. от С; выбор бесконечного р, грубо говоря, означает, что диагональными элементами модифицированного фактора D будут модули диагональных элементов обычного (если обычное разложение возможно). Учитывая приведенные соображения, величину р вычисляют по формуле p--maxiv, ем\. Здесь Bj,- машинная точность. Опа вводится в формулу расчета Р, чтобы обеспечить устойчивость вычислений, когда норма очень мала. Итак, что же такое представленная процедура модифицированной факторизации Холесского, если охарактеризовать ее двумя фразами? Это - численно устойчивый алгоритм, генерирующий по-JЮЖитeльнo определенную матрицу, которая может отличаться от исходной только диагональными элементами. Это - оптимизированный алгоритм в том смысле, что параметр р подбирается в нем путем минимизации априорной оценки нормы поправки Е при условии сохранения существенно положительно определенной матрицы неизменной. Следует отметить также, что реальная величина нормы Е почти всегда оказывается значительно меньше априорной оценки. Фактическое значение нормы Е можно дополнительно уменьшить, если использовать симметричные перестановки столбцов и строк Gf. На очередном, /-м шаге факторизации в качестве /-й строки и i-ro столбца надо брать ту из нетронутых ранее пар, для которой величина \j в (4.23) максимальна. Такая стратегия приведет к разложению вида PGtP-{-E = LDLT, где Р - некоторая перестановочная матрица. В методе с перестановками поправка Е инвариантна относительно нумерации переменных. Ниже приводится детальное описание всех операций, выполняемых по ходу построения модифицированного разложения Холесского с перестановками. Дана наиболее эффективная схема организации расчетов. Все фигурирующие в ией величины прн реализации вычислений на ЭВЛ\ могут размещаться в памяти, первоначально выделяемой для записи исходной матрицы С. При этом коэффициенты рассчитьшаелш1Х факторов занимают места ее использованных элементов. В процессе вычисления j-то столбца матрицы L участвуют вспомогательные величины Cis=tiids, s=l, . . ., /; i=j.....п. Их также естественно заносить в массив, выделенный для G. Точнее говоря, позиции использованных элементов Gj, сначала отводятся под запись чисел с„, а затем, по мере того как необходимость в каких-то Cij отпадает, вместо них записываются соответствующие коэффициенты фактора L. Алгоритм МС (метод разложения модифицированного разложения Холесского). МС1. [Расчет порога для элементов.] ВычислитьР=тах{у, где v=max j 1. Vtf-1}, а у и суть максимальные значения модулей диагонального и недиагонального элеметггов Gj. МС2. [Инициализация.] Присвоить индексу столбца / значение 1. Положить Cfi = g,-j, ( = 1, ...,п. МСЗ. [Перестановка строк и столбцов,] Найти индекс д, такой, что \с..\= птах \с„\. Поменять местами все данные, от-)<<" вечающие столбцам матрицы G с номерами q к i, а затем проделать то же самое с данными, отвечающими ее q-R и /-Й строкам, 1ЛС4. [Поиск максимальной по модулю величины/,.yd,..] Вычислить с,7 = Ви-Д ДЛЯ t"=/-f 1,..., 1 и найти 6= тах ]с,, (если / = п, взять ву = 0). МС5. [Расчет /-го диагонального элемента фактора С] Вычислить dy = max{6, Суу, е?/р2} и поправку £у ==dy-c,y. Если jn, вычисления прекратить. МС6. [Расчет /-й строки L.] Вычислить tjC/Jd, fui« s=l, ... / - 1. Для i - п пересчитать диагональные элементы с,-,, по формуле Сц = с,.,-1,,Си. Присвоить / значение /-1-1 и вернуться к шагу МСЗ. Для построения модифицированного разложения Холесского нужно выполнить около Ч,п арифметических операций (примерно столько же, сколько требуется для посгроения обычного разложения для положительно определенной матрицы). Применительно к матрице из примера 4.7 факторы модифицированной факторизации Холесского и поправка Е при Р=1.061 (с точностью до четырех значащих цифр) таковы: г 1 о o-v 0.2652 1 D \,0.5303 0.4295 1. 3.771 О f О О Е.760 О О О I и 1.121 У f 2.171 О О \ О 5,016 О I О О а.мз) в данном случае [С,-С=Яг» 6.154. Модифицированная факторизация Холесского помимо прочего позволяет определять и направления отрицательной кривизны. Пусть S-индекс, отвечающий минимальной из п-/ величин Cjj на шаге МСЗ представленной выше процедуры. Если матрица Сд не является знакоопределенной, при каком-то / окажется, что величина Cgg отрицательна. Тогда вектор р, найденный как решение уравнения Lp - e, будет направлением отрицательной Рис. 4к. Траектория поиска минимума фуикшш Розенброка модифицированным методом Ньютона. « - 1.861. кривизны. Так, для матрицы из примера 4.7 получим s = 3, и для соответствующего вектора р II р На рис. 4к изображена траектория метода минимизации с вычислением вторых производных, основанного на применении модифицированного разложения Холесского, в задаче с функцией Розенброка (пример 4.2). Как видно из этого рисунка, после первой итерации метод четко отслеживает дно оврага, почти «идеально» аппроксимируя его кусочно-линейной кривой. Замечания и избранная библиография к разделу 4.4 Известно много модификаций метода Ньютона. Большинство ранних численно неустойчивы (подробный разбор этих методов читатель найдет в работе Мюррея (1972а)). Подход с вычислением собственных векторов и значений был предложен Грииштадтом (1967). Метод с модифицированным разложением Холесского разработан Гиллом и Мюрреем (1974а). Бытует мнение, что модифицированные ньютоновские алгоритмы работают хорошо лишь вблизи решения. Оно породило массу «гибридных» схем, в которых на начальных итерациях, пока текущая точка далека от оптимальной, используется какой-нибудь относительно гр>€ый алгоритм (обычно метод наискорейшего спуска), а алгоритм ньютоновского типа включается на завершающей стадии счета. Мы на основании своего вычислительного опыта беремся утверждать, что объективных мотивов для такой гибридизации нет. Хорошо продуманный н тщательно реализованный метод ньюто.юв-ского типа будет эффективен и вдали от искомого решения. В подтверждение этого тезиса можно еще раз сослать:;я иа при.\1ер с функцией Розенброка. В то время как мгтод наискорейшего спуска с первых же и гераций работает из рук вон плохо (см. рис. 4j), ньютоновский метод с модифицированной факторизацией Холесского уже со второй итерации считает вполне удовлетворительно (см. рис. 4к). В основу модификации метода Ньютона можно положить знаконе-определеиное симметричное разложение Банча и Парлета (1971). Они показали, что для любой симметричной матрицы С/, существует представление prG„P = LBU, где Р - перестановочная, L - нижняя треугольная, а В - блтно-диагональная матрицы, причем блоки В имеют ра.змерности, ие превышающие двух. Возможность симметричных перестановок отражена в этой формуле потому, что такие перестановки необходимы для численной устойчивости процедуры расчета рассматриваемого разложения. Знаконеопределенная симметричная факторизация обладает одним существенным достоинством; она позволяет выяснить инерцию матрицы Сд (инерция - это тройка целых чисел, задающих количества положительных, отрицательных и нулевых собственных значений). Дело в том, что инерции Сд и фактора В ее разложения совпадают (хотя сами наборы собственных значений у этих матриц, вообще говоря, различны). При этом двумерные блоки фактора В строятся так, что у каждого из них всегда есть одно отрицательное и одно положительное собственные значения. Следовательно, число положительных собственных значений Сд равно сумме числа поло-житетьных одномерных блоков матрицы В и числа ее двумерных блоков. Идея использования знаконеопределенной симметричной факторизации для построения Сд принадлежит Море и Соренсену (1979). Состоит оиа в том, чтобы, получив спектральное разложение для В (BUAW) и cocTaBjiB матрицу В по правилу разд. 4.4.2.1, т. е. положив B = UAU, где Х,- = тах{Я,, 6}, [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 |