Численное нахождение обменной матрицы при решении задачи оптимизации распределения источников тепла Текст научной статьи по специальности «Математика»

Брусенцев А.Г., д-р физ.-мат. наук, проф., Осипов О.В., канд. физ.-мат. наук, ст. препод.

Белгородский государственный технологический университет им. В.Г. Шухова

ЧИСЛЕННОЕ НАХОЖДЕНИЕ ОБМЕННОЙ МАТРИЦЫ ПРИ РЕШЕНИИ ЗАДАЧИ ОПТИМИЗАЦИИ РАСПРЕДЕЛЕНИЯ ИСТОЧНИКОВ ТЕПЛА

ov.osipov@gmail.com

Описывается метод нахождения оптимальной плотности источников тепла в движущейся среде внутри области, находящейся в состоянии стационарного теплового баланса с окружающей средой. Получаемые источники имеют минимальную суммарную интенсивность и обеспечивают приемлемое распределение температуры внутри области. Строятся конечномерные аппроксимации задачи, обладающие особым свойством регулярности по функционалу. Это свойство показывает принципиальную возможность численного нахождения квазирешения исходной задачи. Обсуждаются алгоритмы его приближённого нахождения, а также результаты вычислительных экспериментов, проведённых с помощью специально созданного программного комплекса HeatCore.

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

1. Введение. Одним из видов объектов, широко распространённых в различных областях человеческой деятельности, являются системы источников тепла в области пространства, находящейся в состоянии стационарного теплового баланса с окружающей средой. При математическом моделировании таких систем часто возникает связанная с ресурсосберегающими технологиями инженерная задача оптимального распределения источников тепловых полей. В типичной постановке эта задача состоит в таком выборе распределения источников, при котором создаваемое ими температурное поле наименьшим образом отличается от заданного. В качестве критерия оптимизации здесь обычно выступает квадратичный функционал, обладающий свойством коэрцитивности. С математической точки зрения эта задача относится к задачам оптимального управления для эллиптических краевых задач. Существование решений и общие свойства подобных задач для квадратичных целевых функционалов, а также приближённые методы их решения изучались рядом авторов (см. [1,2] и приведённую там библиографию). Эту задачу можно также отнести к обратным задачам теплопроводности (ОЗТ), методы приближённого решения которых рассмотрены в [3]. Однако и здесь речь идет, в основном, о коэрцетивном квадратичном целевом функционале.

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

тепла, которое обеспечивает заданный температурный режим при минимальной суммарной мощности источников. При численном решении этой задачи возникает ряд трудностей, и до настоящего времени в полном объёме она практически не рассматривалась. Целевой функционал здесь является линейным и в силу наличия поточечных неравенств, определяющих температурный режим, возникают значительные трудности в установлении существования точного решения. Вообще говоря, точного решения этой задачи может не существовать. В работе [4] введено так называемое квазирешение, которое всегда существует и с прикладной точки зрения вполне приемлемо. В работах [4-6] предлагается способ конечномерной аппроксимации задачи, на основании которого разрабатываются основные алгоритмы, и создаётся методика приближённого нахождения квазирешения. Эти аппроксимации образуют последовательность задач линейного программирования и, как показано в [4,6], эта последовательность обладает особым свойством регулярности по функционалу. Это свойство обеспечивает принципиальную возможность приближённого нахождения квазирешения. Однако построение конечномерных аппроксимаций основано на очень трудоёмкой процедуре определения так называемой обменной матрицы. В настоящей работе изложен один из наименее трудоёмких способов нахождения этой матрицы, и приведены результаты численных экспериментов, показывающих эффективность разработанного метода.

2. Формулировка оптимизационной задачи

Введём функцию и(х) = Т(X) — Т0 , где Т (X) — установившаяся температура в точке х

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

хАп-V-(ум) + / = 0, X е Б, (1)

где х - коэффициент температуропроводности среды, V (X) - поле скоростей среды, которое предполагается подчинённым условию V = 0 и потенциальным, а / (X) - плотность источников тепла в области В. Для формулировки краевых условий разобьём границу области В на три части: дБ = Г+ ^Г0 , где Г+ - часть границы, которая является входом среды в область В, Г- - часть границы, являющейся выходом (стоком) среды, а Гд - часть непроницаемой границы. Справедливы следующие соотношения:

(Я, V)!Г= 0, (Я, V)!Г_> 0, (Я,Г+< 0,

где Я - единичный вектор внешней нормали к границе дБ. Последнее из этих трёх условий означает, что в область В поступает вещество внешней среды с температурой Т0. В дальнейшем мы предполагаем, что в нашем процессе присутствует поток тепла через непроницаемую для среды границу, равный

а(X) - п(X) (X еГэ). При этом коэффициент теплопередачи а(X) > 0 является функцией, заданной на Г0. Будем считать выполненными следующие краевые условия:

ды

X--+ аы

дп

= 0,

ды

дп

= 0,

(2)

Г

ХдЫ _ (П, v)u)

дп

= 0.

Физически эти условия выражают тепловой баланс области D с окружающей средой.

Будем считать, что поле скоростей области D является потенциальным и его потенциал (р( X) известен. Нас интересует оптимизационная задача нахождения функции f (x) > 0, принадлежащей пространству L2(D) и доставляющей минимум линейному функционалу

I{f} =j f (X)dVm ^ min, (3)

D

при условии, что выполнены равенства (1) и (2), а также неравенства

M(X) - To > u(X) > m(X) - To, при X e D. (4)

Здесь m(X), M(X) - минимальный и максимальный профили температур, которые считаются непрерывными функциями, задаваемыми

в некоторой подобласти D ^ D, которая в дальнейшем называется областью контроля температуры. Отметим, что случай dD = Г0

(Г+, Г_=0) не исключается. В этом случае естественно p(X) = const. Не исключается и

случай, когда D = D. Существование точного решения задачи (3) при условиях (1), (2), (4) в пространстве L2(D) гарантировать нельзя.

Обозначим через у нижнюю границу значений функционала I{ f } (3), когда f(X) пробегает множество неотрицательных функций из L2(D), удовлетворяющих условиям (1), (2), (4).

Определение 1. Квазирешением оптимизационной задачи (3) при условиях (1), (2), (4) с данным допуском s > 0 назовём такую функцию fo(X) > 0 из L2(D), удовлетворяющую ограничениям задачи, для которой выполняется неравенство I{f0} ^ У + s.

Для приближённого нахождения квазирешения дадим более удобную эквивалентную формулировку оптимизационной задачи. Введём новую неизвестную функцию w(X) с помощью

равенства u = wep/(2x\ Подставляя это выражение в (1), получим следующее уравнение:

- x&w + (| Vp |2 /(4*))w = f • e"p/(2^). (5)

Из условий (2) следует, что функция w должна удовлетворять краевому условию

rcw

дп

- + ow

= 0,

(6)

dD

где функция а(X) > 0 задана на границе . 0 ^ г _ области следующим обра-

dD = Г+^ГП^Г_

зом:

x ) =

а( X)/X, X еГъ

1 dp ^

--, X еГ_;

2х дп

1 dp „

---, X еГ,.

2х дп '

Оператор Ьм = -хАм + (| Vp | /(4х))м , действующий в пространстве Ь2(В) на достаточно гладкие функции X), подчинённые краевым условиям (6), является самосопряжённым и положительно определённым, а поэтому имеет ограниченный обратный оператор, определен-

Г

0

Г

+

ный на Ь2(П) (см., например, [7, с. 129]). Отсюда следует, что краевая задача (5), (6) имеет единственное решение, которое можно выразить через функцию Грина оператора Ь. Последнее обстоятельство играет существенную роль в дальнейшем.

Если учесть описанное выше преобразование краевой задачи, то получим следующую формулировку оптимизационной задачи

g (x) = Дx)e

-V( x)/(2*)

J{g} = Je^(x)/(2^g(x)dV ^ min,

D

62(x) > Gg(x) > в1(x), при x e D; g(x) e L2(D), g(x) > 0, при x e D,

где

(7)

6\( x) = (m( x) - To)e

x)/(2*)

02( x) = (M (x) - 70 )e

x)/(2%)

а О — оператор, обратный по отношению к оператору Ь. Для задачи (7) тоже можно ставить задачу о нахождении квазирешения go (х) , при этом квазирешение исходной задачи

/о(X) = ^(Х)/(2ж)go(X) .

3. Конечномерная аппроксимация оптимизационной задачи и её регулярность по функционалу.

Для приближённого нахождения квазирешения построим конечномерную аппроксимацию задачи (7) в виде задачи линейного программирования. Разобьём область В на п частей

D = У Dj . Определ:

им

подпространство

j=1

кусочно-постоянных функций вида g(х) = gj, х е О^ (| = 1,2,...,п) . Введём в Sn(D) базис, состоящий из функций е^ (х) = 1, х е О|, и е^ (х) = 0, х £ О^ . Тогда

g(x) = Ё gjej (x). Разбиение области D

j=1

мы

считаем и разбиением области О ^ О, т.е. при некотором натуральном р справедливо равен-

~ Р

ство О = ^ О1 . Рассмотрим на S(D) оператор

г=1

Оg и обозначим

а1]=(тв8В1У\Ое],е1), а^твяВ)1 (01,е,), Ьг=(тв8В;)-1 (02,в1), где (•,•) — скалярное произведение в Ь2(В). Заменяя класс функций S(D) подпространством Sn(D), умножая скалярно ограничения на базисные функции е{ (х) , получаем конечномерную аппроксимацию задачи (7)

Jn (f) = Ё Cjgj ^ min:

j =1

(8)

<Ёaijgj < 6,-, gi > 0 (i = 1,2,..., p).

j=1

Здесь cj = Je^(x)/(2x)dV . Самым

трудным

D,

при построении задачи (8) является нахождение матрицы aiJ=(mesD1)^} (Ge^e,), поскольку оператор G явно не задан. Эту матрицу мы и называем обменной матрицей. Её определение равносильно нахождению функций Wj = Ge;, которые являются решениями уравнений -XAw+(|V9|2/(4x)) w=ej при краевых условиях (6). Вопросу о нахождении обменной матрицы в основном и посвящена настоящая работа.

Задачу (8) обозначим через Z0(n). При неограниченном дроблении области D получаем последовательность задач линейного программирования, которую мы и называем конечномерной аппроксимацией задачи (7). Обозначим через (Jn)min минимальное значение целевой функции задачи Z0(n).

Определение 2. Конечномерную аппроксимацию назовем регулярной по функционалу, если справедливо неравенство

lim( Jn )min <У , n^ro

где у определено перед определением 1.

При наличии регулярности по функционалу решение задачи Z0(n) при достаточно малом

max (diam D,) можно считать приближённым

i

квазирешением исходной задачи. Действительно, после дискретизации система ограничений означает, что неравенства в исходной задаче (7) практически выполнены, т. к. они удовлетворяются в среднем по D,. При этом значения (J„)min при достаточно больших n не превосходят у + в.

В работе [6] доказано, что при m < 3 и усло-

виях

1) m(x) - 70 >¿0 > 0; 2) lim max (diam Dj ) = 0

n^-ro j

(9)

конечномерная аппроксимация является регулярной по функционалу.

Назовем разбиение области D, состоящее из подобластей Di, дроблением разбиения, состоящего из подобластей , если это разбиение

n

одновременно является и разбиением каждой тл*

области Д* .

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

* 1 £ £(т^Д.)аг]. (10)

ам =

тв$Д,

'к ДуеД, ДеД*

В дальнейшем мы описываем способ определения обменной матрицы (а,у) для «тонкого» разбиения области В, а затем с помощью формулы (10) переходим к матрице (а*), отвечающей «укрупненному» разбиению, относительно которого решается задача (8).

4. Определение обменной матрицы и поля скоростей среды. Метод последовательного усреднения

Нахождение обменной матрицы. Предположим, что область В, имеющая форму параллелепипеда в К", разбита на ячейки, которые мы для простоты считаем одинаковыми множествами: отрезками при т=1, квадратами при т=2 и

|2

кубами при т=3. В дальнейшем считается, что ячейки возникают в результате разрезания области точками (т=1), двумя семействами параллельных сторонам В прямых (т=2) или тремя аналогичными семействами плоскостей (т=3). Ребро ячейки обозначим через к. Эти ограничения на форму ячеек не обязательны, но сильно упрощают итоговые формулы.

Искомые матричные элементы

ау=(те5Вг)-1^е;,ег). Их определение равносильно нахождению функций м1=Ое, которые являются решениями уравнений -ХАw+(|Vф|2l(4x)) ■w=e;■ при краевых условиях (6). Зафиксируем у и выпишем приближенные соотношения между элементами у-ого столбца матрицы ау. Каждый из этих элементов отвечает своей ячейке В,. Ячейки разделим на внутренние и граничные. Первые не примыкают к границе, а относительно вторых мы предполагаем, что, по крайней мере, одна из «граней» находится на границе области. Рассмотрим вначале внутреннюю ячейку, в которой находится источник с плотностью еу (X) . К данной ячейке примыкают

2т соседних ячеек, имеющих с ней общие грани. Очевидно, что

2т

дм

0 =| (-хАм-еу + (|^|2/4хК)^т =-х£ | — ^-кт +{(VI/4Х)^ . (11)

Д,

к=1 д*Ду

Д

Здесь дкВу - граница с £-ой соседней ячейкой, а w является решением уравнения -ХАw+(|Vф|2l(4x)) ■w=e;■ при краевых условиях

г дм )

(6). Обозначим через

дп

дм

интегральное

среднее по £-ой грани от —, матричный эле-

дп

мент для внутренней ячейки с источником через аив, а матричные элементы, отвечающие сосед-

ним ячейкам - через а®. Получаем следующую формулу:

С Рл \

— | = (а(к) -аив)/к + о(к) , при к ^ 0.

-дп )к

Подставляя последнее выражение в (11) и

I I 2

обозначая значение | VI | в центре выбранной

2

ячейки через V , получим равенство

2т

- Хкт-1 £ (а(к) - аив)/ к - кт + у)аивкт / 4х + о(кт ) = 0.

к=1

Из последнего равенства получим

1

2т

аив = 2т2 I л 2 £

2т + V к / 4х к=1 Обозначим матричный элемент, для внутренней ячейки, свободной от источника через асв. Аналогично предыдущему получим

а

(к)

+

к2

9 9

2тх + у к / 4х

+ о(к ).

1

2т

асв

£ а(к) + о(к2).

2т + у2к2 / 4х2 к=1

Предположим, что ячейка является граничной. Соответствующий матричный элемент обозначим через аг. У такой ячейки в зависимости от т, количество 5 граней, примыкающих к гра-

нице, может быть равно 1, 2, или 3. Для граней, примыкающих к соседним ячейкам, как и раньше, считаем справедливым равенство

= (а(к) - аг)/ к + о(к),

,дп )к

а для 1-ой грани, примыкающей к границе, можно считать выполненным равенство

fdw ^ кдп,

= аг + o(h),

Из (11), как и ранее получим для граничной клетки с источником выражение:

где cl - значение функции c(X), определенной в (4), в центре /-ой грани.

1 2m-s

аиг =

I

а^ > +-

h2

+ o(h2).

.2 k=1

(2m - s) + hl cl + h v / 4х" l=1

а для свободной граничной клетки - выражение:

1

(2m-s)x + hxI&i + h v /4х

l=1

2m-s

асг =■

I a(k) + o(h2) .

(2т — + Н£ а1 + 2 / 4х /=1

Отбрасывая члены порядка о(й2), полученные соотношения между элементами у-го столбца обменной матрицы можно записать в вектор-но-матричной форме а = Та + у , где Т — линейный оператор усреднения, у — вектор, определяемый источником. Если сходится матричный ряд

1п + Т + Т2 + ... + Т^ + ...,

2 k=1

то можно записать

а = (In -T)-1у = у + Ty + T2у +... + T^y + T^+1а

Это равенство показывает, что вектор а можно приближенно заменить вектором

ак =у + Ту + Т2у + ...+ Тыу, (12)

причём допускаемая при этом относительная погрешность (измеряемая в специально выбранной норме) не превышает величины ||Т||^+1, где

Т — операторная норма матрицы Т. Для операторной нормы матрицы Т, при следующем вы-

боре нормы вектора а = max

а,

, справедливо

неравенство

<-

2m

2m + v02h2 /4х2

где Vq = min V^l. Таким образом, обратимость

D

матрицы In-T и сходимость соответствующего матричного ряда гарантированы при условии v0 > 0. Разрешимость системы а = Ta + у и сходимость описанного выше алгоритма можно установить и без этого условия. Результаты вычислений показывают быстрое убывание с ростом N погрешности в приближенной формуле (12).

Метод, сведения задачи к рекуррентной формуле (12), назовем методом последовательного усреднения.

Нахождение потенциала поля скоростей среды также можно произвести методом последовательного усреднения. потенциал ф( х) является решением следующей краевой задачи Лф = 0 ;

дф

дп

= 0, ^ дп

= s(X), %

дп

= -s2( X),

где положительные функции st(x), s2( x) считаются известными и удовлетворяющими

условию Js1(X)dS = Js2(X)dS , которое

дГ_

дГ

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

условие

О

Пусть разбиение области D такое же, как и раньше. Среднее значение потенциала поля скоростей среды по ячейке D]■ обозначим

ф = h

JD ф(X)dVm •

Далее, как и выше, можно написать

2m - дф

0 = JD Аф(XWm = IJ

k =1

—ds.

9kDj дп

дф

Пусть -| - интегральное среднее по kI дп

k

дф

ой грани от — . Рассмотрим вначале внутрен-

дп

нюю ячейку Dj. Обозначим среднее потенциала

для этой ячейки

Ф

,(k)

1в-

а средние для соседних

ячеек — через Ф( ). Можно считать выполненным равенство

l

s

s

s

Г

Г

Г

m

д(р I = (ф^) - ф .в) /h + o(h) при

v дп

'k

Таким образом,

1 2m

=

/в 2m j

£ф/ ) +o(h2).

k=1

Для граничной ячейки среднее потенциала обозначим Фуг . Обозначим количество граней,

примыкающих к границе области, через tj, а количества граней, примыкающих к частям Г-, Г+, соответственно через ^у, t2у . Принимая во внимание (6), получим для ячейки, примыкающей к Г_ или Г+,

Ф =-1-j 2m -1

2m-t ,■

Еф1

(k )

j k=1 Здесь s\p, s(/

'1 j '2 j Е s|/,)-Е s2i)

l =1

l =1

h+o(h ).

средние величины плот-

ностей входящего и выходящего потоков среды по /-ой грани ячейки с номером у.

Для ячейки, примыкающей к Г0 , справед-

ливо равенство

Ф„ =

1 2m-t.

1 s: ф(к)+o(h2).

;г 2m -1. k=1 j

Отбрасывая погрешность аппроксимации o(h2), получим векторно-матричное уравнение

(In - T)ф = S,

(13)

где T - линейный оператор усреднения, определяемый равенством

(тф) / =

1 2m-t/

2m -1

Е ф

,(k )

j k =1

в котором предполагается, что для внутренней ячейки tj = 0. Компоненты вектора 5 сопоставляются ячейкам и отличны от нуля только для ячеек, примыкающих к Г- или Г+, причем для этих ячеек

(S )/ =

Л

Е si/)-Е

j. / ) 2/

/=1

/=1

h.

Отметим, что матрица 1п - Т является вырожденной, поскольку для вектора С с одинаковыми координатами (1п - Т)С = 0. Поэтому

при произвольном векторе 5 система (13) не имеет единственного решения. Отметим также,

что Е (S) / = 0 . Если принять во

внимание

/=1

условие , то мы должны разыс-

Д

кивать решение системы (13) в подпространстве векторов Ф, удовлетворяющих условию

п с

£ (Ф) у = 0 , а в этом подпространстве суще-

у=1

ствует единственное решение. Хотя это подпространство не является инвариантным для оператора Т, можно показать, что итерационный процесс

N

(ф)( N ) = S+Е Tks k=1

сходится к решению системы (13).

5. Численная реализация алгоритмов. Результаты вычислительных экспериментов.

Для численного решения рассмотренной задачи создан программный комплекс HeatCore, написанный на языках C#/C++. Имеется возможность решения как двумерных, так и трёхмерных задач с возможностью задания физических характеристик среды, параметров сеток и количества источников тепла. Краевые условия и функции m(x), M(x) задаются в виде скриптов на языке C#. Визуализация результатов расчёта: положение и мощность источников, поле скоростей, температурное поле осуществляется с помощью графического модуля, использующего OpenGL. Задача линейного программирования решается классическим симплекс-методом с использованием искусственного базиса. HeatCore использует равномерные сетки и позволяет решать задачу оптимизации распределения источников тепла для простых областей: квадрат (m=2) и параллелепипед (m=3). Для решения дифференциальных уравнений, кроме приведённого метода последовательного усреднения, так же реализован набор конечно-объёмных схем и методов решения разрежённых СЛАУ.

На рис. 1 приведена блок-схема алгоритма решения га-мерной задачи для движущейся среды с использованием метода последовательного усреднения для нахождения обменной матрицы. На рис. 2 приведены графики стабилизации значений (Jn)min с ростом n, полученные в результате численных экспериментов в трёхмерном случае, соответственно для неподвижной и движущейся сред с различными значениями коэффициента температуропроводности %.

t

t

Рис. 1. Блок-схема алгоритма решения т-мерной задачи

а) неподвижная среда

б) движущаяся среда

Рис. 2. Зависимость На рис. 3 приведён результат численного эксперимента для 3-мерной области В размером 5^4x3 м с принудительным источником холода. В качестве вещества, заполняющего область, был выбран воздух с коэффициентом теплопроводности х = 0,0244 Вт/(м-К), плотностью р = 1,293 кг/м3 и удельной теплоёмкостью с = 1005 Дж/(кг-К) (расчёт произведён в системе СИ).

значения (Лп)тт от п Границы области выполнены из кирпича с коэффициентом теплопроводности 0,5 Вт/(м-К) толщиной 30 см. В правую грань области входит вещество со скоростью $(ху,1) = 10-4 м/с (Г+ имеет площадь 2 м2). Симметрично с левой стороны имеется участок стока Г такой же площади с такой же скоростью вещества 8(х,у,г) = -10-4

м/с. Минимальный и максимальный профили = 260 К. Область контроля занимает всю область

температур заданы функциями т(х-у,г) = 290 ^ за исключением границы О = О — дО . М{х.у.г) = 320 К, температура внешней среды Т0

к = 0,0244 Вт/(м К) т(х,у,2) = 290 К «=10x10x10 а(лг, у,г) = 5/3 Вт/(м;-К) Т„ = 260 К М{х,у,г) = 320 К Г| = 2

Рис. 3. Оптимальное расположение плотности источников тепла внутри параллелепипеда.

Разноцветными объёмами показана искомая оптимальная плотность источников тепла (красные) и сток (синий). Менее прозрачные объёмы соответствуют более мощным источникам. Источники распределяются вблизи места входа вещества Г+, границы, а также около стока, таким образом его нейтрализуя. Суммарная интенсивность источников в системе СИ имеет значение (/п)шт = 432,739 Вт.

Для проверки эффективности был проведён следующий эксперимент. Было взято некоторое случайное распределение источников (рис. 4а) и решена прямая задача теплопроводности для той же самой области (рис. 3), но без стока тепла. Параметры среды были взяты следующие:

Х(х,у,г) = 2,216х10-5 м2/с, а(х,у,г) = 1,7^10 м/с. На дальней плоскости имеется окно с большим значением а(х, у, г) = 3,4х10-5 м/с. В результате было получено, что случайное распределение источников нагревает область контроля в диапазоне температур 292,8-364,8 К. Затем была решена прямая задача с т(ху,г) = 292,8

К и М(хуг) = 364,8 К при (х, у, г) ^ О . В итоге получилось, что случайное распределение (рис. 4б) имеет суммарную интенсивность 0,02687 К-м3/с, а оптимальное - 0,01476 К-м3/с. Сэкономленная мощность при этом составляет около 45 %.

а)неоптимальное б) оптимальное

Рис. 4. Расположение источников в параллелепипеде

6. Заключение. Для достаточно точного вычисления обменной матрицы требуется численно решать семейство краевых задач с использованием тонких сеток. Авторы так и поступали в работах [4]-[6]. Однако с ростом п в конечномерной аппроксимации Z0(n) количество этих краевых задач, равное количеству ячеек В, становиться большим, особенно в трехмерном случае. Это может приводить к невозможности завершить вычисления за разумное время. Поскольку величина п заранее не известна, возникает потребность в менее затратном способе нахождения обменной матрицы, который позволил бы производить надежные численные эксперименты. Метод последовательного усреднения и является таким способом. Численные эксперименты показывают, что при относительно большом значении отношения % / а последовательность (Лп)тт с ростом п сравнительно быстро стабилизируется вблизи у и можно обойтись небольшим значением п. С уменьшением величины % / а скорость стабилизации этой последовательности снижается. Поэтому нахождение квазирешения в широком диапазоне значений параметров модели может потребовать использования аппроксимации Z0(n) с большим значением п. При этом далее приходится решать задачу линейного программирования с п переменными, что тоже может привести к значительным затратам машинного времени.

Отметим, что при одних и тех же параметрах модели результаты численных эксперимен-

тов, произведенных с помощью решения семейства краевых задач и с использованием метода последовательного усреднения практически совпадают.

БИБЛИОГРАФИЧЕСКИЙ СПИСОК

1. Лионс Ж.-Л. Оптимальное управление системами, описываемыми уравнениями с частными производными. М.: Мир, 1972, 412 с.

2. Федоренко Р.П. Приближенное решение задач оптимального управления. М.: Наука, 1978, 497 с.

3. Алифанов О.М. Обратные задачи теплообмена. М.: Машиностроение, 1988, 280 с.

4. Брусенцев А.Г., Осипов О.В. Приближенное решение задачи об оптимальном выборе источников тепла // Научные ведомости БелГУ. Математика. Физика. №5 (124). Выпуск 26. Белгород, 2012. С.60-69.

5. Осипов О.В. Оптимальное расположение источников тепла в неоднородной среде // Вестник Белгородского государственного технологического университета имени В.Г. Шухова. №1. 2013. С.154-158.

6. Брусенцев А.Г., Осипов О.В. Оптимальный выбор источников тепла при наличии конвекции // Научные ведомости БелГУ. Математика. Физика. №26 (169). Выпуск 33. Белгород, 2013. С.64-82.

7. Михлин С.Г. Вариационные методы в математической физике. М.: Наука, 1970.

Brusentsev A.G., Osipov O.V.

NUMERICAL CALCULATION OF THE EXCHANGE MATRIX WHEN SOLVING THE PROBLEM OF OPTIMIZATION OF DISTRIBUTION HEAT SOURCES

Describes a method for finding an optimal density of heat sources in a moving environment within the region in a state of steady-state thermal balance with the environment. The obtained sources have minimum total intensity and provide an acceptable temperature distribution within the region. Construct finite-dimensional approximations of the problem, have the special property of the regularity of the functional. This property allows you to numerically find quasisolution original problem. Algorithms discussed its approximate location, as well as the results of computational experiments conducted with a specially crafted software package HeatCore.

Key words: density heat sources, the inverse heat conduction problem, the optimal control problem for elliptic boundary value problems, finite-dimensional approximation, convection, heat balance, the simplex method.

Брусенцев Александр Григорьевич, доктор физико-математических наук, профессор кафедры программного

обеспечения вычислительной техники и автоматизированных систем.

Белгородский государственный технологический университет им. В.Г. Шухова.

Адрес: Россия, 308012, Белгород, ул. Костюкова, д. 46.

E-mail: brusentsev@mail.ru

Осипов Олег Васильевич, кандидат физико-математических наук, старший преподаватель кафедры программного обеспечения вычислительной техники и автоматизированных систем. Белгородский государственный технологический университет им. В.Г. Шухова. Адрес: Россия, 308012, Белгород, ул. Костюкова, д. 46. E-mail: ov.osipov@gmail.com