Федоренко Р.П. Введение в вычислительную физику (1185915), страница 101
Текст из файла (страница 101)
о' о, (9) Не станем выводить этой формулы, укажем лишь основные опера- ции ее вывода. Подставляем в (б) вместо т и з соответственно т + Ьт и з+ Ьж Пользуясь тем, что []з[[ > е в гх'хб,, разлагаем подынтегральное вы- ражение в ряд Тейлора с учетом первого порядка малых величин. Заменяя Ьх на — ~ ~ В Ьт (в соответствии с (8)), получаем выраже- ние для ЬР [Ьт( )] в виде четырехкратного интеграла. В этом ин- теграле меняем очередность интегрирования (переобозначая х', у' через х, у и наоборот), «Внутренний» интеграл и есть,0(х, у). Что касается приращения недифференцируемой части, то тут никакик особых упрощений нет.
Итак, приращение функционала при вариа- ции т может быть записано в форме Р[т( ) + Ьт(.)] — Р[т( )] = ~ ~ (Ю(х, у), Ьт(х, у)) с1х х(у + а + ~ ~ У ]1] 1х(х, у) — ~ ~ В(х — х', у — у') Ьт(х', у') 4(х' 4(у' ~~ 41х 41у— с, с — ~~ У(х, у) [[з(х, у)]] 4(х с(у+ О([[ЬХ[]х). (1О) Пренебрегая О([[ЬХ[[з), определим процедуру вычисления вариации Ьт, обеспечиваквцей убывание Р. Обычно решение таких сложнык задач начинают с относительно простых алгоритмов. Совсем не обязательно в конкретных расчетах должны появиться все неприятности, которые в принципе возможны.
Начнем поиск минимума с начального приближении т(х, у) ш О. Пока Ор Я и функционал «зо пгиьлижзнньш мзтоды вычислительной»итаки ~ч. и дифференцируем, работает метод спуска по градиенту. На каждом шаге т пересчитывается по формуле Т:=Р(т — ЬР), где Р— оператор проецирования, работающий локально, в каждой точке (х, у) независимо от других, Я вЂ” шаг процесса. Оператор Р введен с целью учета поточечного ограничения [) т(х, у) )) < /(х, у) и реализуется просто. Определяем т* = т — ЯР и, если [) т" (х, у) [[ ~ Г, полагаем т(х, у) = т'(х, у)/(х, у)/)) т'(х, у))).
Вычисляем фактическую вариацию Ьт(х, у) = т(х, у) — т(х, у) н предсказанную вариацию функционала ЬР= ~~ (В, Ьт) Ых Ыу. с Выбор шага процесса 5 играет большую роль, если нас интересует не только факт сходнмости, но и скорость процесса минимизации. Затем вычисляем новое значение функционала Р[т( ) ) и фактическое приращение ЛР= Р[т( )) — Г[т( )). Если ЬР > О, итерация считается неудачной, шаг 5 уменьшается вдвое и с той же производной повтори- ется вариация т.
Если ЛР с О, итерация выполняется, т.е. т заменяется 'на т, а шаг Б корректируется в зависимости от точности линейною приближенна — величины и = 2[ба — ЛР[/[ЬР+ ЬР[. Если эта величина мала, шаг Я увеличивается, если велика, — уменьшается. Заметим, что трудоемкость итерации велика: она определяется необходимостью вычисления функционала Р и его производной Р. Обе операции стоят О(л 4) операций (в каждой точке двумерной сетки нужно вычислить двумерный интеграл). Поэтому здесь не применяется надежный способ выбора шага, связанный с решением задачи ш1п Р[Р(т — БР)) по 5. Даже не очень точное ее решение требует нескольких вычислений Р.
Расчеты по этой простой схеме показали, что сначала функционал достаточно быстро убывает, затем образуется небольшая область сцепления б,, которая растет. По мере ее роста все большую роль в ЛР начинает играть недифференцируемая составляющая, шаг Б катастрофически уменьшается и алюрнтм «застревает» в заведомо неоптимальный точке т( ). Метод регуляризации. Наиболее простой и дешевый способ продвинуться дальше, почти не усложняя алгоритма, состоит в регуляризации задачи, т.е. в данном случае в аппроксимации недифференцируемой функции [[ [[ дифференцируемой.
Практически это „ях'Р ~я+,+* <. 0>.в, у ° но затем, по мере достижения минимума (для данною з), постепенно уменьшать, используя найденное ранее решение как начальное приближение при новом значении з. Методы регуляризации задач весьма популярны, но, к сожалению, не очень эффективны. Дело в том, что нужно не только ап- $29] влгилционныз зьдьчи мах«ники проксимировать недифференцируемую функцию дифференцируемой (эта цель легко достигается в данном случае), ио и получить функцию, с хорошей точностью аппроксимируемую своей касательной на таких расстояниях от точки линеаризации, которые следует использовать в эффективном алгоритме построения минимизирующей последовательности.
При замене ] з ] иа г Р + е возникает конфликт между точностью аппроксимации и гладкостью регуляризоваиной функции. В описываемом алгоритме разумный компромисс между точностью и гладкостью аппроксимации достигался следующим образом. Наряду с основными счетными массивами т«, зк, П«использовался массив е«, и вместо ]]з„]] в формулы входила величина (]]з«]]з+ е„)пз. После осуществления вариации (т- т+ Ьт, з- г+ Ьз) величины з пересчитывались.
При этом предполагалось, что на следующей итерации значение Ьз будет примерно таким же. Для хорошей линеаризации нужно, чтобы значение ]]з«]]з+ е«было раз в пять больше ожидаемой вариации Ьз Таким образом, величины е„автоматически убывали в процессе расчета при уменьшении шага 5. Таблица 21 дает представление о том, как протекал процесс минимизации. В ней представлены: номер шага ч (звездочкой отмечены неудачные итерации с йР > О), г, ЬР, АР и шаг спуска Я. Видно, что неудачные шаги сравнительно редки. Обратим внимание на то, что величина ЬР (в принципе пропорциональная 5]]Щ]з) убывает намного быстрее, чем Я. Это связано с убыванием производной 2У, т.е. с приближением к минимуму. Возникает вопрос: насколько рациональна вырабатывающаяся в ходе расчета величина Я? Свидетельством в пользу этого алгоритма служит хорошо видный из таблицы факт: обычно резкое уменьшение шага Я сопровождается ростом фактического убывания функционала Лг.
Экспериментальные попытки волевым образом увеличить 5 не приводили к успеху: возникала ситуация пг" > О, шаг последовательно дробился несколько раз подряд и приходил к старому значению. Вышеприведенный расчет носил методический характер, поэтому число итераций относительно велико (впрочем, размерность конечномерного пространства, в котором решалась дискретная задача на минимум, здесь был» около 1500), Видно, что после 15-й итерации расчет практически «стоит на месте» и продолжение его бесполезно. Каковы же полученные при этом результаты? Анализ показал, что почти всюду в области гг (это был эллипс) значение функции 1 = ]Дз]] — (т, з) ~ было очень мало. Более точно это означает следующее.
В начальном приближении среднее значение 1, 2 максимальное значение 1„,„, 5. В конце расчета (при Рж0.18 в большей части области для 1 0.0002 . 0.002) пгиклижкннык мктолы вычислительной бианки !ч. н 482 1,р 0.027, У„а„, = 2.2, Имеется небольшая подобласть в б, в которой значение функции 1!х, у) достаточно велико, причем вместо требуемой в точном решении коллинеарностн в и 3 наблюдалась почти антиколлинеарность этих векторов. В этот момент с.ь эанима- Таблица 21 ла значительную часть гг и никакие ухищрения в рамках описанной выше мегодики не приводили к улучшению решения. (Впрочем, как показали дальнейшие расчеты, для грубых выводов часто бывает до- статочно и полученного таким образом решения.) 0 ! 2 3" 3 5' 5 б 7* 7 В 9 9 !О 11 гг 13" !з 14 15* 15 16 !7 г! 25 29 33 ! 3.83 4.83 3.46 !.772 1.772 0.988 0.864 0.864 0.477 0.4! З 0.4! З 0.307 0.268 0.268 0.24! О.237 0.229 0.212 0.21г олюг О. !97 ОЛ97 0.193 0.19! одвб О.1В4 Онзг 0.18! 10.8 6.0 4.2 З,б 1.8 1.1 !.3 0.67 0,35 0.43 о.гг ол! О.! 5 0.075 0.058 0.072 0,052 О.О45 0,023 О.О1З 0.017 О.ОО9 0.0056 0.0064 0.0029 О.ОО!7 0.0009 0.0005 8.8 !.36 1.69 -0.46 олв олз — 0.002 0.39 0.064 -О,О!В ол! 0.039 -О.О25 0.026 0.004 0.008 0.0! б -0 0014 о.о! ! 0,005 -о.оо! О,ОО4 олюг 0.00! 0,0006 0.0005 0.0002 о.ооог гг 27 гг гг 11 11 8.6 4.3 4.3 3.5 1Л !.7 !.7 0,86 О.вб 0.69 055 0.55 0.28 0.28 О,гз 0.14 ОЛ4 0.14 0.07 0.036 0.029 о.огз ВАРИАЦИОННЫЕ ЗАДАЧИ МЕХАНИКИ 483 Для того чтобы получить более аккуратные и достоверные результаты, пришлось существенно усложнить метод.
Метод линейного программирования. Этот сложный алгоритм мы опишем в самых общих чертах. Ои основан на некоторой аппроксимации приращения недифференцируемого слагаемого Г, более простым, ио тоже недифференцируемым. Используется возможность аппроксимировать круговой конус, как поверхность в трехмерном пространстве (С, В, ~ф+ цз), шестигранным конусом (г„ть Ф(С, ц)), где Ф(г„д) определяется решением следующей задачи линейного программирования: ФД, В) = 054 ш!и (а'+ а" + р'+ р«+ у'+ т") «',...Н« при условиях ('-')М+(2//ЗНК-а") =В, -Е'-К')+( '-Т") =Ч, (а' — а") + (р' — р") + (у' — у") = О, а', а", р', р", у', у" в О. Не будем доказывать этого. Читатель, желающий понять, в чем тут дело, пусть начнет с вопроса о том, почему ) Ц = ппп (а'+ а") при а' — а" =ч, а ~0, а'вО, Используем введенную аппроксимацию конуса.
Введем, кроме Ьт(х, у), новые вспомогательные переменные а'(х, у), ..., /"(х, у), (х, у) ~ С,. В терминах этих переменных задача выбора направле- ния спуска для недифференцируемого функционала ставится следу- ющим образом. Требуется найти Ьт(х, у), а'(х, у), ..., У«(х, у), обеспечивающие ппп ((~ (/З(х, у), Ьт(х, у)) Их Иу+ с + 0.54 ~ ~ /(х, у) (а'(х, у) + а "(х, у) + ... + у" (х, у)) Их а~у) о, при условиях з,(х', у') — ~ ~(Ь,(х — х', у — у'), Ьт(х, у)) с(х Ыу— с — «а' — а")/Л + (2/з/З)(р' — Д"))„, „, = О, зз(х, у ) — ~ ~ (Ьз(х — х, у — у ), Ьт(х, у)) Фх ну+ + «~'-Д") — (у — у"))„ч, =О, «а' — а") + (р' — р«) +(у' — у«))„, „, =О, т'(х', у') Е с .