Стохастический метод внешних аппроксимаций для решения выпуклых задач полубесконечной оптимизации. II. Вычислительный эксперимент
С.К.Завриев, А.В.Федосова* МГУ им. Ломоносова, факультет ВМК
November 4, 1999
1. Постановка задачи и алгоритм для решения выпуклых задач полубесконечной оптимизации.
Рассмотрим задачу полубесконечной оптимизации
Р(У°) : найти ж £ Д^,
Xopt = {х ex0 \ f(x) = min f(x')}
X0 = {xeX°\g(x,y)<OVyeY0},
X° С 'Rk. У0 С 'R1. Функции f(x) и g(x,y) предполагаются непрерывно дифференцируемыми и выпуклыми по х на выпуклом компакте Х° для любого у £ У0, 'R'" - т-мерное пространство.
В данной постановке множество У0 задает бесконечную систему ограничений, определяющих допустимое множество Х° задачи Р(У°).
Для решения задачи Р(У°) применим метод внешних аппроксимаций, который заключается в сведении исходной постановки задачи Р(У°) к решению последовательности аппроксимирующих задач P(Yn), п = 1.2...., следующего вида:
Р(У) : найти ж £ Xopt[Y],
X&IY] = {xG X\Y] | f(x) = min f(x>)},
x'€X\t \
*19alina0cmc.msu.ru
Х[У] = {х £ Х° I д(х, у) < О VУ е У}.
На основании указанного метода в [4] построен идеальный вариант алгоритма SMETH.ACTIV.sip и доказана теорема о сходимости. Рассмотрим реальный вариант этого алгоритма для конкретной выпуклой задачи полубесконечной оптимизации Р(У°).
Вспомогательная процедура SPROC.ACTIV.sip призвана формировать множество А¥п = {у*} для добавления к множеству ¥п:
^те+1 := Уп и АУте,
если вспомогательная процедура не установит, что точка х" уже является решением исходной задачи (команда ОСТАНОВ основного алгоритма).
Основная процедура метода SMETH.ACTIV.sip при помощи данных, полученных SPROC.ACTIV.sip, на каждой итерации решает задачу Р(¥п) с конечным числом ограничений.
Для ОСТАНОВА используется параметр 1тах. Если 1гпах раз не нарушается неравенство на Контрольном шаге Процедуры 8РН0С.-АСТ1У^р. то управление передается в основной алгоритм 8МНТН.-АСТП^^р на Шаг ОСТАНОВ и алгоритм заканчивает свою работу.
Метод SMETH.ACTIV.sip
Параметры. гтах < ос.
Шаг 0. (Начальный шаг)
YÏ := 0-
Шаг 1. Находим х" — решение задачи P(Yn). Шаг 2. Выполняем процедуру SРROC.ACTIV.sip с входными параметрами хп ж Yn. Получаем AYn либо выход на ОСТАНОВ. Шаг 3. Формируем новое множество ограничений
^те+1 := Yn U AlÇj.
Шаг 4.
п := п + 1. Переход к шагу 1.
ОСТАНОВ. Из процедуры SPROC.ACTIV.sip получена точка х" — решение исходной задачи Р(У°).
Процедура SPROC.ACTIV.sip
Входные параметры : х G X, Y £ Á4f(Y°). Выходные параметры : AY£Á4f(Y°)
либо ОСТАНОВ.
Параметры процедуры : 7 > 0, £ > 0. Шаг 0. (Инициализация).
i := 1.
Шаг 1. Если i > imax, тогда возврат в SM ET H.ACT IV.sip
на ОСТАНОВ. Шаг 2. (Пассивный случайный поиск критического ограничения).
Найти точку y¡ G У". используя равномерное
вероятностное распределение на У0.
Шаг 2. ACTIV. (Активный поиск критического ограничения).
Применить алгоритм локального поиска для решения
задачи д(х,у) —>• max (метод проекции градиента),
y€Y°
стартуя из точки y¡ для получения точки у* G У0, такой что :
Vi е YstaAx)i 9(х, Vi) < д(х, У*). Шаг 3. (Контрольный шаг). Если
i ■ maх(д(х, у\), ...,д(х, у*)) < 7,
то i := i + 1 и переход к Шагу 1. Шаг 4.
где г* = arg max д(х, у*).
1 <j<i ■'
Выход.
2. Численное исследование алгоритма для задачи упруго-пластического кручения стержня.
Изложим результаты тестирования стохастического метода внешних аппроксимаций для задачи упруго-пластического кручения стержня [1, 2, 3].
Рассмотрим тонкий цилиндрический стержень из упруго-идеаль-нопластического материала, подчиняющегося критерию Мизеса в прямоугольной декартовой системе координат. Пусть ось жз параллельна образующей цилиндрической поверхности стержня, и в Х\Х2 лежит поперечное сечение, которое задается открытой, ограниченной односвязной областью О в К2. Пусть в исходном состоянии стержня отсутствуют напряжения. Стержень с одним закрепленным концом подвергается воздействию возрастающей крутящей пары сил.
Деформированное состояние стержня характеризуется углом закручивания [л на единицу длины стержня. Обозначим через и(х\. х->) функцию напряжений. С возрастанием приложенных внешних сил (и угла закручивания //) в стержне проявляются пластические свойства материала (появляются зоны пластичности). В таких зонах поведение материала перестает быть обратимым (после снятия напряжений сохраняются остаточные деформации).
Исходная краевая задача эквивалентна следующей вариационной задаче [3]: найти функцию и(х\. х-)), удовлетворяющую во всей области О неравенству | Vu| < 1 (для определенности ограничение сверху для |Vu| положено равным единице), на границе области - условию и\у = 0. Г = 0Q. и доставляющую минимум функционалу
Таким образом, вариационная постановка задачи о кручении заключается в определении поля напряжений в деформированном стержне:
где - пространство Соболева 1-го порядка с нулевым следом
на границе Г [2]. Производные понимаются в обобщенном смысле.
Заметим, что функция и(х\. ж-;) будет удовлетворять неравенству |\7у/| < 1 в области упругости, а равенству |\7</| = 1 в области пластичности.
Сведем вариационную постановку (2) к задаче полу бесконечной оптимизации. Для этого исходное бесконечномерное пространство
(1)
(2)
К = {и е НЦ&) |Vu| < 1 п.в. в i)}
аппроксимируем по схеме Ритца, например, на базе ортогональных рядов Фурье. Соответствующее разложение выбираем с учетом граничных условий. При этом число переменных становится конечным (так как пространство из бесконечномерного стало конечномерным), а число ограничений остается бесконечным, ибо они были операторными и должны выполняться V.7' £ Q.
Таким образом, подставив соответствующее разложение в функционал J, и проинтегрировав его по сечению О, получим выпуклую (это следует из выпуклости функционала J) целевую функцию для задачи 57 Р. которая зависит только от коэффициентов ряда Фурье.
Допустимая область задачи полубесконечной оптимизации получается из подстановки разложения функции и(х\. х->) в ограничения
\Vu(xij жг)|2 — 1 < 0 У(ж1,жг) £ О.
Расчеты упруго-пластического кручения проведем для двух видов сечений стержня:
1) О = (0,1) х (0,1) - квадратная область;
2) О - единичный квадрат с "вырезанной" правой верхней четвертью.
Для первого случая используем следующее разложение:
t t
и(х 1,^2) = Л a>ijsin(7rixi)sin(7rjx2). i=ij=i
После его подстановки в функционал (1) получим целевую функцию:
7Г
2 t t ,, t t
ii << ((-!)'-1Ж-1У-1)
о ,¿=1 7Г ,1=1 ¿=1 г 2
где А - квадратная матрица основных переменных (/¡¡. г = 1,..., ] = 1, ...,£, которые являются коэффициентами разложения функции и(х 1, Х2) в ряд Фурье.
Ограничения задачи выглядят следующим образом:
г г
д(А,х) = (X] Л щ^соз{ттгх1)згп{тт]Х2)'кг)2'-\-<¿=11=1
í г
щ^гп^гх^соз^2x2)1^])2 — 1 < О
¿=1¿=1
Ухх £ (0,1) и ж2 £ (0,1).
Таким образом, вариационная постановка задачи упруго-пластического кручения стержня с квадратным сечением (0,1) х (0,1) была сведена к выпуклой задаче полубесконечной оптимизации
Р(О) : Найти А е Д^,
Л& = {АсХ°1/(А)=тшо/(А%
= {АеХ° I д(А, х) < О Уж е О},
с О С к2.
К задаче -Р(П) применим алгоритм SMETH.ACTIV.sip, который был реализован для // = 1, 2,10, 20.
В методе проекции градиента шаг выбирался по правилу Армихо. Для экономии времени счета и получения наиболее точных результатов, сначала задача считалась для 25 переменных с любой начальной точки с высокой заданной точностью 7 (например, 7 = 0.05). Далее решение этой задачи выбиралось в качестве начальной точки для задачи со 100 переменными (т.е. матрица А имела размерность (10x10)). Для улучшения эффективности также учитывались полученные ранее активные точки у £ О, чтобы в задачах с требованием большей точности сократить время счета.
Задача считалась решенной, если 1гпах раз подряд значение
та х(д(х,у1),...,д(х,у*тах))
не превышало г = 1, .... /„и,.,--
Параметр %тах положили равным 500, поэтому -¿Р— = 0.0001 (для
Ътах
варианта с 25 переменными). В критерии останова в задаче со 100 переменными 7 уменьшалось постепенно от 50 до 0.5, поэтому для полученного решения в любой точке сечения выполняется д(А, х) < 0.001.
По полученному решению были сделаны графические рисунки зон пластичности (там, где д(А,х) = 0) и упругости (где д(А,х) < 0). Было получено значение функции и{\11) = 0.412, а тестовое и{\11) = 0.413 для ц = 10 [1].
Случай дг = 1 и дг = 2 соответствует чисто упругому кручению и задача решается за 1 - 2 итерации. Остальные значения // характеризуются наличием зон пластических деформаций. Зона упругости имеет вид креста исходящего из углов квадрата. Чем больше значение дг, тем тоньше становятся области упругости, и тем больше имеем соответствующих активных ограничений. На рисунке 1 приводятся
распределения зон для // упругих деформаций.
10, черная область соответствует зоне
Рисунок 1.
Для сечения с невыпуклой областью О - квадратом (0.1) х (0.1) с вырезанной правой верхней четвертью - разобьем О на три части как показано на Рисунке 2 и зададим и(х\. следующим образом
и(х)
если х £ Ох, и 2- если х £ ГЬ. щ + и2- если х £ О3, где
Х2
о,
í í
у/1 = Л аг]3т{тх\)зт{2/к3x2).
»=17=1
í í
(¡2 = Л Л а^згп{2/кгх\)згп{/к]х2)'
%=\з=\
03
Рисунок 2.
Заданная таким образом ж2) будет непрерывной на всей области определения.
Полученная целевая функция представляется следующим образом:
.... тг (вт (| - 21)) вт^и + Щ
™ = 2 £ § § )-2,--¡+21 I +
]ф21 1ф2]
2гзгп(Щ^) (зт (|(2^ - I)) вт (|(2^ + ¿))\ 15 ( 22^ 1 +
тт.jl (sin (f(i - 2О) Цf(i + 20)\
""""Х ( i - 2¿ + i + 2¿ J +
4jiein(Ç) (sin^(2j-l)) sin (§1 „ 7
a^aik 15¿ (-27^7-+-2ÎT7-) ) +
1 Л Л Л Л / 2ksin(s-f) (sin (^(j-21)) sin + 2¿))' Ö ajiakl-77- -:-7Г,---^ПТТ- +
2 i=2 к i=i k=1 '=1 V 15 V э - 21 3+Ы
j¿2l 1ф2j
irk2 (sin (f(2j - l)) sin (f (2j + l))\
I, 2j2 JTl J +
а^аы Ш [ J^2l + JT21 J +
irjl (sin (f(2j-0) sin (}(2j + 0)^
T ( 2i -1 + 2i + ¿ JJ +
Xf f f f „ vik (sin (U2i-k)) sin(U2i + k))\
2 £ ¿ S £ (^Х ( XX + XX J +
гф2 к кф2г
4ifein(f) ísin (|(¿ — 2fc)) sin (f(i + 2fc))\
—ш— i—ï=2fc— + —TTTk—) +
7r/2 /sin (f (2¿ - k)) sin (1(2¿ + fc))\
а*0и"Г i-2í~k---2ÏTk-j +
2Шп(Щ (sin (f (¿ - 2fc)) sin (f (¿ + 2i-))\\ а,,а1к 15 ^-—^----j j +
1 « « ' f / Trj2 (sin (|(¿ - 2fc)) sin (|(¿ + 2fc))^
2 S g S Si Г^Т ---Tï-~2k J +
¿#2fc fe#2¿
2jsin(^-) (sin (f (2г - A;)) sin (f (2г + fc))\
°*аы 15 [ X^ J +
шк (sin (f (г - 2fc)) sin (f (г +
t l-2k + i + 2k )+
4iksin(3-Y) /sin (|(2i — к)) sin (|(2¿ +
°*аи-1:,, 1--+-STTfe-j j +
«EEEE \ajtakl—+ atjalk-V2 ; 2 ) {"2 1 >2
z г=1 j=2lk=2il=\ \ 4 ZZDÍ2
1 * * * * / 7Г2 16sin(^)sin(%¿)\ ,;2 ,2.
xEEEE [atjalk-+ ajtakl-2 +
* i=2fc j=l jfe=l 1=2 j \ 4 ZZD/C? /
1 * I I I f sin(^)(k2+ 412) sin(^)(4k2+ l2)\
« E EEE ojíow—^^-¿ + oy-o/fc—^^-¿ | +
¿ i=2ki=2lk=ll=l\ K 1
1 ^ ^ ^ ± ( sin(^f)(4i2+j2) sin(^)(i2 + 4j2)\
^EEE E ajiakl--—:-+ aijaik--—;- +
¿ i=l j=í k=2il=2j \ J 1 J
1 ' ' ' ' / ik (sin (l{2i-k)) sin (f (2¿ + k))\
2 hhhh ri¡ ИТ l. + 2¿Tfc J
гф2к кф2г 1ф2у
'sin (f (i - 20) sin (f (i + 2Í))\ ik (sin (f (i - 2fc))
+ «ú'^fcTT -;—^-+
j-2l j + 2l J 4 2 \ г -
sin(f(i + 2&))\ (sin (%(2j - l)) sin
i + 2k J \ 2 j-l 2 j + l
jl (sin (f(2i - к)) sin (f (2i + k))\ ^ (sin (f (j - 21)) ~2 [ 2 i-к 2 i + k
+
aji^M
V l i-2/
sin (f (i + 2Z))\ jl (sin (f (i - 2fc)) sin (f (г + 2fc))\
r—J + —^ r^ J
'sin (}(2j-Q) sin (f(2j + Q)\\
2 j-i + 2i + ¿ + .¿'¿.....
Ограничения в этом случае представляются как:
g(A,x) = (JZ E aijCos(7rixi)sin(27rjx2)7ri) + i=ij=i
t t 2 (E E aijsin(7rixi)cos(27rjx2)27rj) — 1 < 0, если x G Oi; i=ij=i
t t 2 i=ij=i
t t 2 (E E ajisin(27rixi)cos(7rjx2)7rj) — 1 < 0, если x G O2;
¿=ij=i
t t ¿=ij=i
t t 2
t=lj=l
t t
t=lj=l
t t 2
E E a>ijsin(7rixi)cos(27rjx2)27rj) — 1 < 0, если ж G O3. t=lj=l
Подобная постановка задачи полубесконечной оптимизации соответствует Р(П) с рассматриваемой областью О, и, следовательно, для ее решения применим алгоритм SMETH.ACTIV.sip.
Замечание.
Точки, лежащие на отрезках х\ = 0.5, могут относиться к области Q\ или О3, т.к. функция и(х\. Х2) является непрерывной по всей области задания О. Аналогично для точек, где х-> = 0.5.
В связи с тем, что в данном случае область О не выпуклая, метод проекции градиента реализовался каждый раз на области Qj. i = 1,2,3, куда попадала точка у-,. выброшенная равномерным вероятностным распределением (Шаг 1 SPROC.ACTIV.sip). В этом случае "более нужная" точка у* (Шаг 1 ACTIV из SPROC.ACTIV.-sip) попадает на ту же самую Q;. которой принадлежала и у,.
Стохастический метод внешних аппроксимаций для данного примера исследовался для размерности t = 12, т.е. 144 переменных, начиная с t = 5, постепенно уменьшая требуемую точность выполнения ограничений -Р— с 0.1 до 0.001. Данная задача считалась с теми
^шах
же остальными параметрами для // = -у. 10, 20, 40.
Распределение областей упругости и пластичности для // = 10 показано на Рисунке 3. Геометрические результаты расчетов согласуются с решениями, полученными в [3, 1] другими методами.
Рисунок 3.
Здесь также с увеличением // наблюдается увеличение множества активных ограничений и, соответственно, времени счета и числа итераций.
Выполненное исследование показало, что предложенный алгоритм стохастического метода внешних аппроксимаций, использующий активизированный "двухфазный" поиск ограничений, дает хорошее приближение для решения выпуклых задач SIP и наглядные результаты с возможностью варьирования физических параметров задачи. Алгоритм удобен в реализации, обладает хорошими свойствами сходимости и может быть рекомендован для решения выпуклых задач полубесконечной оптимизации.
References
[1] Гловински Р., Лионе Ж. Л., Тремольер Р. Численное исследование вариационных неравенств. М.: Мир, 1979.
[2] Дюво Г., Лионе Ж. Л. Неравенства в механике и физике. М.: Наука, 1980.
[3] Черноусько Ф. Л., Баничук Н. В. Вариационные задачи механики и управления. М.: Наука, 1973.
[4] Завриев С.К., Федосова А. В. Стохастический алгоритм внешних аппроксимаций для решения выпуклых задач полубесконечной оптимизации. Обоснование сходимости / Исследовано в России, (в печати).