УДК 519.6
ПРИБЛИЖЕННОЕ РЕШЕНИЕ ЗАДАЧИ ОБ ОПТИМАЛЬНОМ ВЫБОРЕ ИСТОЧНИКОВ ТЕПЛА
А.Г. Брусенцев, О.В. Осипов
Белгородский государственный технологический университет им. В.Г. Шухова, ул. Костюкова, 46, Белгород, 308012, Россия, e-mail: [email protected], [email protected]
Аннотация. Предлагается и обосновывается метод численного решения задачи об оптимальном выборе плотности источников тепла. Приводится описание основанных алгоритмов и результаты численных экспериментов.
Ключевые слова: плотность источников тепла, обратная задача теплопроводности, функция Грина, конечномерная аппроксимация, симплекс-метод.
Задача об оптимальном расположении источников тепловых полей всегда являлась актуальной задачей при проектировании в строительстве, металлургии и других областях техники. Она допускает ряд постановок, которые не эквивалентны из-за различий в критериях оптимизации. Здесь мы рассматриваем задачу нахождения плотности источников тепла минимальной мощности, которая обеспечивает заданный температурный режим в некоторой области в условиях ее стационарного теплового баланса с окружающей средой. Возможные постановки и пути решения такой задачи обсуждались в [1]. С математической точки зрения эта задача относится к задачам оптимального управления для эллиптических краевых задач. Существование решений и общие свойства подобных задач для квадратичных целевых функционалов, а также приближенные методы их решения изучались рядом авторов (см. [2, 3], и приведенную там библиографию). Нашу задачу можно отнести также к обратным задачам теплопроводности (ОЗТ), методы приближенного решения которых рассмотрены в [4]. Однако и здесь речь идет, в основном, о квадратичном целевом функционале. В данной работе целевой функционал является линейным. Трудности в решении задачи порождены входящими в формулировку неравенствами, описывающими заданный температурный режим. Кроме того, отсутствие коэрцетивности у линейного функционала порождает трудности в установлении существования решения задачи. Мы видоизменяем формулировку задачи и строим конечномерную аппроксимацию этой задачи в виде последовательности задач линейного программирования, а также показываем, что эта аппроксимация обладает особым свойством регулярности по функционалу. Последнее свойство позволяет считать решение конечномерной задачи с достаточно большим номером, приближенным решением исходной задачи. В работе приводится также описание основанных на построенной аппроксимации вычислительных алгоритмов и результаты численных экспериментов.
1. Формальная постановка задачи и ее конечномерная аппроксимация
В ограниченной связной области DcRm требуется определить функцию f (x) ^ 0, доставляющую минимум линейному функционалу
J{f} = J f (x) dVm ^ min , (1)
D
при следующих условиях
X^u + f = 0,
, (2)
(du/dn + ßu) |aD = 0 ,
M(x) — T0 ^ u(x) ^ m(x) — T0 . (3)
Здесь X,T0 — известные константы, ß = ß(x) — известная функция, определенная
на границе dD, m(x), M(x) — задаваемые в области D минимальный и максимальный профили температуры, которые считаются непрерывными функциями. Вводя функцию Грина G(x, £) для оператора Лапласа с краевым условием (2), условия (2), (3) можно заменить следующими
m(x) — To ^ —1 iG(x,i)f(i)dVm ^ M (x) — To, f(x) ^ 0 , f(x) E S (D). (4)
X i&D
В оптимизационной задаче (1), (4) S(D) является классом функций, среди которых разыскивается плотность источников f (x). Ниже в качестве S(D) берется вещественное гильбертово пространство L2(D) квадратично интегрируемых в D функций.
Построим конечномерную аппроксимацию задачи (1), (4) в виде задачи линейного программирования. Разобьем область D на n частей <D = У Dj >. Определим
I j=1 J
подпространство Sn(D) С S(D) кусочно-постоянных функций вида f (x) = fj, x E Dj (j = 1, 2,...,n). Введем в Sn(D) базис, состоящий из функций ej (x) = 1, x E Dj, и
ез(Х) = 0, х/Ву. Тогда f (х) =^2¡зе^(х). Рассмотрим на 5(В) оператор
3 = 1
=-(1/х) [ С(х,ОЦОЛ'т .
Введем также обозначения аз = (Сег, ез), (т(х) — Т0, ег(х)) = аг, (М(х) — Т0, ег(х)) = Ьг, где (• , •) — скалярное произведение в Ь2(В). Построим конечномерную аппроксимацию задачи (1), (4)
Jn{f} = J2(mes Dj )fj ^ min,
j=i
n
ai ^ ^ aij fj ^ bi, fi ^ 0 (г = 1, 2,...,n). j=i
В одномерном случае (т =1) такая аппроксимация строится легко, т.к. нетрудно найти функцию Грина 0(Х) в явном виде. Область В здесь является интервалом (а,Ь), а краевая задача (2) примет форму
{ (х) = 0 / (Аи/АХ - в1и) 1 х=а =0;
йХ2 ’ \ (йи/йх + в2и) 1 х=ъ= 0 ,
где вг,в2 > 0 — некоторые константы. Функция Грина запишется в виде
Иаъ(х - 1ъ)(£ + 1а), при а ^ ^ х;
[
С(х,£)
Яаь(х + ¡а)(£ - ¡ь), При X ^ ^ а,
где ЯаЬ = вів2/(вів2 (ь - а)+ ві + в 2) , ¡а = (1 - віа)/ві, ¡ь =(1+ в2Ь)/132.
При т > 1 для построения задачи (5) численными методами находится функция иі = Ові, которая является решением краевой задачи (2) с f = ві. Далее задача (5) решается симплекс методом.
Отметим несколько естественных модификаций задачи (5). В случае достаточно большой функции М( X) - Т0, условия сверху в (4) на решение задачи не влияют. После конечномерной аппроксимации получаем задачу без ограничений сверху (Ьі = ж). Такую задачу будем называть односторонней задачей вида (5) и обозначать ^+(п). Двустороннюю задачу (5) обозначим Z0(n).
При постановке задачи до дискретизации, кроме условия (4), возможно появление и других условий, которые могут привести к некоторым новым модификациям задачи (5). Одна из естественных модификаций состоит в требовании невозможности расположения источников тепла в некоторой части области В, т. е. возникает дополнительное требование: f (X) = 0 при X Є В0 С В Оно приводит к тому, что в задаче (5) уменьшается количество переменных. Исчезают некоторые слагаемые в целевой функции и ограничениях. Количество ограничений остается равным п.
Вторая модификация связана с присутствием некоторых фиксированных источников до оптимизации. Другими словами, функция f (X) в (1) состоит из двух слагаемых, одно из которых известная функция, а второе подлежит определению. После дискретизации для этого второго слагаемого получим задачу вида (5) с другими значениями чисел аі,Ьі. В целевой функции появится свободный член. И в этом случае вполне удобен симплекс метод. Разумеется, общая схема решения задачи останется неизменной и в случае совмещения последних двух модификаций, для которых возможен как двусторонний, так и односторонний варианты.
Замечание 1. Плотность источников f (X) в задаче (1), (4) можно варьировать без изменения значения целевого функционала 3{f}. Действительно, используя формулу Гаусса-Остроградского и условия (2), получим
3^} = J f (X)dVm = -х ! ^пйУт = -х J = X / ви^8 ■
Б О дО дО
Таким образом, значение целевого функционала полностью определяется значениями на границе dD функции п(х), которая удовлетворяет краевым условиям (2), также заданным лишь на границе области. При этом функцию п(х) можно выбирать внутри области D произвольно, лишь бы удовлетворялось неравенство Ап ^ 0. То есть f ( х) = -хАп(х) можно изменять, не меняя значения J{f}. Эти соображения показывают, что если задача (1), (4) имеет точное решение, то оно далеко не единственно. Кроме того, плотность источников f ( х), удовлетворяющую (4), при достаточно широких условиях можно изменить без изменения значения функционала J{f} так, чтобы ограничения сверху в (4) стали строгими неравенствами.
2. Уточнение постановки задачи.
Регулярность конечномерной аппроксимации по функционалу
Ответ на вопрос о существовании точного решения задачи (1), (4) существенно зависит от класса функций S(D) и здесь мы этот вопрос не рассматриваем. В настоящей работе мы рассматриваем задачу отыскания плотности источников суммарной мощности сколь угодно близкой к нижней грани функционала при условиях (4). Обозначим y = inf J{f}. Пусть также (Jn)min — минимальное значение целевой функции f е(4)
задачи (5) с n переменными.
Определение. Конечномерную аппроксимацию (5), (т.е. последовательность задач Z0(n) или Z+ (n)) назовем регулярной по функционалу в классе S(D), если справедливо неравенство
Иш (Jn)min ^ y ■ n—
При наличии регулярности решение конечномерной задачи при достаточно большом n можно считать приближенным решением нашей уточненной задачи. Действительно, система ограничений в (5) означает, что неравенства в (4) удовлетворяются в среднем по Dj. При этом значение (Jn)min при достаточно больших n приближенно не превосходит y■
Условия регулярности для последовательности задач Z+ (n) дает
Теорема. Пусть последовательность разбиений области D удовлетворяет условиям:
1. lim max(diam Dj) = 0;
n—<x> j
2. ßn/un ^ C, где ßn = max(mes Dj), vn = min(mes Dj),
jj
а C — независящая от n константа. Если при всех х ED и некотором k0 > 0 выполнено неравенство т(х)—Т0 ^ k0, то последовательность односторонних конечномерных задач Z+(n) регулярна по функционалу в классе L2(D).
Рассмотрим сначала задачу Zs(n), которая отличается от задачи Z0(n) тем, что a и bi в ней заменены на ai — 5 и bi + 5 при некотором 5 > 0. Обозначим через In(5) минимальное значение целевой функции задачи Zs (n).
Лемма 1. Если задача 2о(п) имеет решение, то при любом 8 > 0 справедливы неравенства
2п
(Зп)шт ^ !п{8) ^ (З'п)шт 8 ^ ^ Уг ,
г=1
где у г(г = 1, 2п) — координаты какой-нибудь точки максимума задачи, двойственной по отношению к Z0(n).
□ Поскольку точка минимума задачи Z0(n) удовлетворяет системе ограничений задачи Z¡(п), то справедливо неравенство (3п)т\п ^ 1п(8). Целевая функция задачи, двойственной по отношению к задаче Zs (п), имеет вид
п п 2п
^(у) = ^ агУг ЬгУг+п — 8^Уг-
г=1 г=1 г=1
Согласно первой теореме двойственности, если обозначить через 3Т (у) целевую функцию задачи, двойственной к Z0(n), а через у0 точку максимума задачи двойственной к Z¿ (п), получим
2п 2п 2п
1п(8) = 3 (У°) — 8^уг0 ^ —3п (уг) — 8^2/Уг = (3п)т\п — 8^^уУг,
г=1 г=1 г=1
что завершает доказательство леммы 1. В
Лемма 2. Пусть (3п)т;п и 1п(8) — минимальные значения целевых функций задач Z+ (п) и Z+ (п) соответственно. Тогда справедливы неравенства
\Jn)mm ^ In(ß) ^ (Jn)min ( 1
1 mm ai
i
(1 - .
у min ai J
□ Задачу Z+ (n) можно считать частным случаем задачи Z0(n) при достаточно больших bi. При этом, для точки минимума задачи Z0(n) ограничения сверху в (5) являются строгими неравенствами, а это, в силу условий дополнительной нежёстко-сти, означает, что в точке максимума двойственной задачи yi = 0 (i = n, 2n). Поэтому
n 2n
(Jn)min = Jn (У) = Yl aiyi ^ (min aj)^2i yi. Отсюда и из леммы 1 вытекает справедли-
n2n
JT
П
i=l i i=l
вость леммы 2. В
Лемма 3. Для произвольного числа а > 0 и любой функции ф(х) £ Ь2(П) найдется такое число е > 0, что при выполнении для разбиения области В неравенства
П
шах(йгатВу) <е можно выбрать кусочно-постоянную функцию /( х) = ^ /*г* ( х), для 3 3=1
которой выполнены условия 3{/} = 3{ф}, ||/ — ф|| < а, где || • || — норма в Ь2(П).
□ Хорошо известно, что для любой функции ф(х) £ Ь2(Б) можно найти такую функцию ф(х) £ О0(С), что ||ф — фЦ < а/2. Ввиду равномерной непрерывности ф(х)
существует такое е > 0, что при условии max(diam Dj) < е можно найти кусочно-
j
постоянную функцию f1(x), для которой ||/i — ф\\ < а/2 или f — ф|| < а. Зафиксируем разбиение области D и рассмотрим задачу на минимум для функции n переменных
П
Ф(/) = \\f — ф||2 = J2(fi ■ (mes Di) — 2fi ■ (ф, ei)) + ||ф||2. В стационарной точке справед-
i=1
ливы равенства fi ■ (mes Di) — (ф,еi) = 0, i = l,n. Поэтому выбирая fi = (ф,е^) /mes Di, получим кусочно-постоянную функцию f (x), для которой J{f} = J{ф}, ||f — ф|| < а.
□ Доказательство теоремы. Пусть а — произвольное положительное число. Выберем функцию 0 ^ ф(х) £ Ь2(П), удовлетворяющую левому неравенству в (4), такую, что 3{ф} — 7 < а. Согласно лемме 3 вследствие условия 1) теоремы, для всякого разбиения области В, связанного с задачей Z+(n) при достаточно большом п, найдется кусочнопостоянная функция /( х), для которой выполнены условия 3{/} = 3{ф}, I/ — фЦ < а. Эта функция может не удовлетворять левому неравенству в (4). Однако она удовлетворяет системе ограничений задачи Z+(n) при достаточно большом п и некотором 8 > 0. Действительно Оф = О/ + О(ф — /) ^ т(х) — Т0. Умножим скалярно в Ь2(0) обе части неравенства на ег(х). Получим
п
(О/, гг) = ^2 агз/ ^ аг — (О(ф — /), гг) ^ аг — ||О||(ше8 А)||ф — /1| ^ аг — ЦОЦ^а.
3=1
Тем самым показано, что /(х) удовлетворяет системе ограничений задачи Z+(n) с 8 = ||О||^па. Отсюда, если учесть лемму 2 и условия теоремы, получим
3 {/} ^ 1п(8) ^ (3п)т1п(1 — ^Уп (^пко) 1а) ^ (3п)т1п(1 — |О|Ск0 1а) •
Таким образом, для любого сколь угодно малого а > 0 найдется такой номер Ыа, что при п ^ Ыа выполнено неравенство (3п)т\п ^ (^ + а)(1 — ЦОЦСкО^)01, которое и доказывает теорему. В
Эта теорема точно также доказывается для других, упомянутых выше, односторонних модификаций рассматриваемой задачи.
Замечание 2. Регулярность по функционалу последовательности двусторонних задач Z0(n) легко устанавливается по приведенной схеме, если для задач Z0(n) и Z(п) удается доказать справедливость леммы 2. Ее доказательство сводится к установлению существования точки минимума задачи Z0 (п), для которой правые неравенства в (5) являются строгими. В пользу существования такой точки минимума говорит замечание 1.
3. Описание алгоритмов и результаты численных экспериментов
Для приближенного решения задачи (1), (4) с т = 1, 2, 3 создан программный комплекс ИеаЮоге. Он использует сведение задачи к задаче (5) при равномерном разбиении
области по каждой координатной оси. На блок-схеме на рис. 1 приведен общий алгоритм решения -мерной задачи с использованием численного метода для вычисления коэффициентов ау. Прямая задача (2), сведённая к разностному уравнению Аи = ej, составленному интегро-интерполяционным методом [5] (методом баланса), решается методом Гаусса на сетке, частота которой по каждой осп в и раз больше количества разбиений по соответствующей оси. В одномерном случае предусмотрена возможность получения значений коэффициентов аналитическим путем. Сравнение результатов при различных способах (аналитическом и приближенном) определения а^ показывает приемлемость приближенного нахождения этих коэффициентов. При т > 1 для достижения приемлемой точности определения необходимо значительное машинное время. Для случая т = 2 проведено большое количество экспериментов на различных областях. Удалось произвести также просчет нескольких вариантов при т = 3.
Heatcore (mix) M(x),f3(x),x,T0,D,v)
)
j := l,n
>
fix) =
Решение уравнения (2) с 'l, fix) e Dj,
[0, fix) i D} численно на сетке. При этом каждая подобласть Dj заменяется сеткой с количеством точек (v + 1)™. Результат - значение температуры Uj(x).
При численном решении данного уравнения результат интерпретируется как приближённое значение интеграла Uj(Xi) = -± f G(Xi,£)f(£)dVm te-Dj в узлах сетки
1,п
>
;Hj = J uj(x)dVm. Di
Решение задачи (-5) с использованием симплекс-метода..
Результат - оптима льная плотность источников fmin
^Возврат^™^
Гис. 1. Блок-схема алгоритма решения m-мерной задачи.
Вычисления показывают, что с ростом п последовательность (3/г)т;п довольно быстро стабилизируется. При этом оптимальная плотность источников имеет несколько от-
дельных зон локализации возрастающей остроты. На рис. 2 представлен результат численного решения задачи (5) для отрезка (?п= 1). Для изображённых здесь результатов решения задачи при аналитическом вычислении а^, уже при 7 разбиениях минимум полностью стабилизируется и равен /т;п = 65 К°- м/с.
Отметим, что при численном нахождении этих коэффициентов /тіп также оказывается равным 65, независимо от частоты сетки. Для иллюстрации эффективности разработанного метода на рис. 2 приведено сравнение результатов при выбранном произвольно расположении источников (кривые серого цвета) с найденным методически (чёрного цвета). Оба распределения источников соответствуют одному температурному коридору. Оптимизация распределения источников в данном примере экономит (76, 25 — 65)/65 • 100% ~ 17, 3% тепловой энергии.
а 71 °К
X = 5 м2/с ДО) = 4 м1 т(х) = 5 °К п = 40 Jmш = 65 °К-м/с
Т= 0 /?т = 9 м1 М(х) = 6.5 °К ./ = 76.25 °К-м/с
0,1 0,2 0,3 0,4 0,5 0,6 0,7 0,8 0,9 л, м
Рис. 2. Оптимальное и не оптимальное распределение плотностей источников / и соответствующих температур Т на отрезке.
На рис. 3 представлен результат решения задачи (5) для квадрата (т. = 2).
о
Рис. 3. Распределение оптимальной плотности источников / и соответствующей температуры Т в квадрате.
В таблицу сведены результаты численных экспериментов для области, представленной на Рис. 3. Для решения задачи область разбивалась на п частей, и расчеты производились на различных сетках с коэффициентами температуропроводности, равными 5 м2/с и 100 м2/с. Как показывает эксперимент, с ростом коэффициента температуропроводности значение функционала стабилизируется быстрее, и увеличение частоты сетки не оказывает большого влияния на точность.
Таблица 1
Зависимость значения функционала ^п)тт от п
п V = 2 г/ = 4 V = 6
X = 5 м2/с Х = Ю0 м2/с Х = 5 м2/с Х = ЮО м2/с Х = 5 м2/с Х = ЮО м2/с
3x3 30,56674 30,02876 30,59955 30,03027 30,60560 30,03055
5 х 5 30,18947 30,00955 30,20090 30,01010 30,20302 30,01020
7x7 30,09682 30,00487 30,10256 30,00515 30,10362 30,00520
10 х 10 30,04754 30,00239 30,05032 30,00252 30,05083 30,00255
15 х 15 30,02116 30,00106 30,0224 30,00112 — —
20 х 20 30,01191 30,00060 30,01260 30,00063 — —
30 х 30 30,00530 30,00027 — — — —
40 х 40 30,00298 — — — — —
На следующем рисунке представлена область (?п = 3) с оптимальными источниками, закрашенными тем темнее и непрозрачнее, чем выше их плотность.
Рис. 4. Распределение оптимальной плотности источников / в кубе.
Литература
1. Брусенцев А.Г., Брусенцева В.С. Задача об оптимальном выборе источников тепла // Сб. трудов XXIII международной конференции «Математические методы в технике и технологиях». Т.2. - 2010. - С.43-46.
2. Лионс Ж.-Л. Оптимальное управление системами, описываемыми уравнениями с частными производными / М.: Мир,1972. - 412 с.
3. Федоренко Р.П. Приближенное решение задач оптимального управления / М.: Наука, 1978. - 497 с.
4. Алифанов О.М. Обратные задачи теплообмена / М.: Машиностроение, 1988. - 280 с.
5. Самарский А.А. Теория разностных схем / М.: Наука, 1977. - 656 с.
APPROXIMATE SOLUTION OF THE OPTIMAL CHOICE PROBLEM
OF HEAT SOURCES A.G. Brusentsev, O.V. Osipov
Belgorod State Technological University named after V.G. Shuchov,
Kostyukova st., 46, Belgorod 308012, Russia, e-mail: [email protected], [email protected]
Abstract. Method of numerical solution of optimal choice problem connected the heat sources density is proposed and justified. Description of algorithms based on the method and numerical results are given.
Key words: heat sources density, inverse problem of heat conduction, Green’s function, finitedimensional approximation, simplex method.