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

АЛГОРИТМЫ ПРОГРАММ ДЛЯ РЕШЕНИЯ ПРЯМЫХ И ОБРАТНЫХ ЗАДАЧ..

7

ТЕПЛОФИЗИКА И ТЕОРЕТИЧЕСКАЯ ТЕПЛОТЕХНИКА

УДК 536.62

АЛГОРИТМЫ ПРОГРАММ ДЛЯ РЕШЕНИЯ ПРЯМЫХ И ОБРАТНЫХ ЗАДАЧ ТЕПЛОПРОВОДНОСТИ ПРИ ИСПОЛЬЗОВАНИИ ДИФФЕРЕНЦИАЛЬНО-РАЗНОСТНЫХ МОДЕЛЕЙ

К.В. Кириллов, Н.В. Пилипенко

Исследовано применение различных модификаций цифрового фильтра Калмана (ФК) для решения граничных и коэффициентных обратных задач теплопроводности. Приведено описание как математических моделей теплопереноса и измерений, так и алгоритмов вычислительных подпрограмм. Представлены результаты тестирования разработанных программ.

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

Введение

Одной из наиболее проблемных задач теплометрии при исследовании промышленных объектов и технологических процессов является определение нестационарных условий теплообмена с помощью приемников теплового потока (ПТП) по измеренным в них температурам или их разностям в отдельных точках. Такие задачи относятся к нестационарным граничным обратным задачам теплопроводности (ОЗТ). Если теплофизические характеристики (ТФХ) ПТП известны лишь приблизительно, то необходимо решать комбинированную ОЗТ: граничную ОЗТ - по восстановлению входящих тепловых потоков и коэффициентную ОЗТ - по идентификации соответствующих ТФХ.

Решение прямой задачи теплопроводности

В качестве математической модели для описания одномерного теплопереноса в ПТП различных типов применяются дифференциально-разностные модели (ДРМ), подробно описанные в работах [1-3], которые в векторно-матричной форме для линейных стационарных систем обыкновенных дифференциальных уравнений (СОДУ) имеют вид:

^ Т (т) = Р • Т (т) + О • и (т) ,

где Т (т) и и (т) - векторы состояния и управления; Р и О - матрицы обратных связей и управления. Общее решение СОДУ (1) имеет следующий вид:

Т (т) = Ф (т, т0 )• Т (т0) +1 Ф(т, 9) • О(в) • Щеу 9,

где Ф(, то ) = ехр(((-то)) - переходная матрица состояния (матрица Коши) системы; то - начальный момент времени. Для программной реализации решения (2) вводится дискретное время т^ = кДт, а также дискретные векторы Тк = Т(тк) и и к = и(тк). Тогда дискретная переходная матрица Ф = Ф к+1 к = Ф(тк+1, тк) может быть вычислена с требуемой точностью путем суммирования необходимого числа членов следующего бесконечного ряда:

Ф = I + РДт + - Р 2 Дт2 +... +—Рт Дтт +..., 2! т!

где I - единичная матрица. Решением прямой задачи теплопроводности (ПЗТ) в этом случае является последовательное применение для каждого момента времени следующей известной формуле расчета Тк+1 по значениям Ф и Тк:

Тк+1 =Ф-Тк +1( I+ Ф) •О• и-Дт.

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

Ук = н • Тк + £ к,

где Ук и £к - векторы измерений и случайных погрешностей; Н - матрица измерений.

Решение обратной задачи теплопроводности

В работах [1, 2] показана целесообразность использования метода параметрической идентификации для решения ОЗТ, так как последний удовлетворяет общепринятым требованиям устойчивости и сходимости вычислительных процедур, точности конечных результатов, универсальности, простоты программной реализации и др. Сущность же метода сводится к предварительной параметризации задачи и последующему нахождению оптимальной несмещенной оценки либо вектора состояния, либо вектора искомых параметров системы, дающей минимум нормы вектора невязки между измеренными в опыте температурами и прогнозами измерений температуры, рассчитанными по модели. Для получения оценок используется рекуррентная вычислительная процедура цифрового ФК. Рассмотрим подробнее два наиболее распространенных ФК: линейный ФК по расширенному вектору состояния системы и нелинейный ФК по вектору искомых параметров.

Под параметризацией ОЗТ понимается априорная кусочно-линейная аппроксимация подлежащего восстановлению теплового потока на всем интервале измерений, где в качестве системы базисных функций применяются 5-сплайны 1-го порядка. Тогда на 2-ом участке аппроксимации значение теплового потока находится по следующей формуле:

= 4а2 ■ ^гА + 4ЬЕ ■

где 4а2 и 4ь2 - значения теплового потока на левой и правой границах участка соответственно; Бр2-1 и

Бр^1 - 5-сплайны. Линейный ФК по расширенному вектору состояния системы (ФК-1) основан на введении расширенного вектора состояния Ял:

R *

т

lzk

IXzi t2zk "' tnzk qaz qbz ]

где Qz = [ (¡а2 ] - вектор искомых параметров, а также на соответствующем расширении ДРМ за счет очевидных уравнений ¿¡а2 = 0, ¿¡ьг = 0 и простейшей коррекции правой части модели измерений.

Алгоритм ФК-1 для одного участка сплайн-аппроксимации описывается следующими уравнения-

ми:

= ®k+1,k -R+ + -2( + фк+u)-Gr • Uk • AT ; Pk+1 =фk+i,k •Pk+•ф^+l,k; K+1 = Pk+1 • HR -(k+i HR + ^);

R++1 = R *+1 + Kk+1 -(Yk+1 - hr R -+1);

P+ = P - K HP-k+1 k+1 k+111 r k+1 ;

где P - ковариационная матрица ошибок оценок; K - весовая матрица; N - ковариационная матрица случайных погрешностей измерений; индексы «-» и «+» обозначают априорные и апостериорные значения, соответственно. Алгоритм ФК-1 обеспечивает нахождение несмещенной оценки Rk, т.е.

^ Rk^ = Rk), дающей минимум дискретной квадратичной функции невязки:

ф(k ) = Z(Yk - HrRk ) -N-1 -(Yk - HrRk ) .

k =1

ФК-1 был реализован в виде программного комплекса «Heat Identification», который непосредственно восстанавливает как температуры, так и входящий тепловой поток, следовательно, его целесообразно использовать в тех случаях, когда начальное распределение температур по толщине ПТП известно лишь приблизительно.

Нелинейный ФК по вектору искомых параметров (ФК-2) основан на введении вектора Qz = [Qqz QAz ] = [qa z qbz Xz ]T , для которого выполняется условие Q = const. Тогда модель ПТП имеет следующий вид:

Q = о, (1)

а модель измерений

Yk = Yk(Qo) + Sk , (2)

где Yk (Q0) - модельный вектор измерений; Q 0 - истинное значение вектора искомых параметров.

АЛГОРИТМЫ ПРОГРАММ ДЛЯ РЕШЕНИЯ ПРЯМЫХ И ОБРАТНЫХ ЗАДАЧ...

К модели (1), (2) может быть применен алгоритм дискретного нелинейного ФК, позволяющий получать рекуррентные оценки (5к+1 вектора искомых параметров О и ковариационную матрицу Рк+1 их ошибок по найденным на предыдущем к-ом шаге (, Рк и известному вектору измерений Ук+1. Алгоритм имеет следующий вид:

K+1 = Pk • H+ • (( k+1 PkH T+1 + ж

<3 k+1 = <3 k + Kk

I Yk+1 Yk+1

) • (Q k ));

P = P - K H P

k+1 1 k ^k+111 k+11 k-

где Нк+1 - матрица функций чувствительности; Ук+1 ((к) - модельный вектор измерения, рассчитываемый по модели теплопереноса в ПТП для момента времени к+1 с использованием предыдущей оценки (к вектора Ок.

Матрица функций чувствительности Нк+1 имеет следующий вид:

>1,к (О) дуик (О) к (О)'

Hk+1 =

3Yk (Q)

5Q

Q = <3 k

Здесь

ty» (Q)

94a

Q = Q k

- функции чувствительности j-го

измерения к искомому параметру qa, q^ и X в k+1 момент времени.

Алгоритм ФК-2 обеспечивает нахождение несмещенной оценки Qk, т.е. E^ Qk^ = E (Qk), дающей минимум дискретной квадратичной функции невязки:

Ф (Qk ) = £ (н - Yk (Qk ) • ж1. (Yk - Yk (Qk )) .

k =1

ФК-2 был реализован в виде программы «Heat Conduction», который непосредственно восстанавливает как тепловой поток, так и теплопроводность, следовательно, его целесообразно использовать в тех случаях, когда теплопроводность материала ПТП известна лишь приблизительно.

Результаты имитационного моделирования

Ниже представлены результаты математического моделирования для градиентного однородного ПТП типа вспомогательной стенки толщиной h = 0,005 м и со следующими ТФХ: X = 15 Вт/(м-К); c = 485 Дж/(кг-К); р = 8000 кг/м3. Входящий в ПТП тепловой поток изменялся по закону q1(x) = [10000sin(0,1x)+10000] Вт/м2, на тыльной стороне q2(T) = 0 Вт/м2. Задавались температуры поверхности ¿1 и второго блока /2 при уровне погрешностей в измерениях О =0,1 °С; длине участка сплайн-аппроксимации Дz = 10.Дт (Дт = 0,01 с); начальном распределении T0 =[30 — 30]T °С.

Результаты восстановления теплового потока и температурного поля по толщине тепломера с помощью ФК-1 представлены на рис. 1. Начальные оценки принимались вдвое меньше эталонных:

R0 =[15 — 15 5000 5000]T, а начальное значение ковариационной матрицы

P0 = diag(100,..., 100,1012,1012).

Результаты восстановления теплового потока и уточнения теплопроводности материала ПТП с помощью ФК-2 представлены на рис. 2. Начальные оценки принимались, как и в предыдущем случае, в

двое меньше эталонных: Q0 =[5000 5000 7,5]T, а начальное значение ковариационной матрицы:

: diag (к

P0 = diag (1012,1012,100)

q, Вт/м2

20000

10000

t, °C 30

25

20

15

/1

-2 .3

4

0,2

0,4 0,6 а

0,8 т,с

0,2 0,4

0,6 0,8 б

Рис. 1. Эталонный (1) и восстановленный (2) тепловые потоки (а); заданная на поверхности первого блока (1') и восстановленные на блоках 1-5 температуры ПТП (б)

q, Вт/м2

20000

10000

0

X, Вт/(м • К) 14 12 10 8 6

0,2 0,4 0,6 0,8 т,с а

0,4 0,6 0,8 т,с б

Рис. 2. Эталонные (1) и восстановленные (2) значения теплового потока (а) и теплопроводности ПТП (б)

Заключение

В статье приведено описание математических моделей, как процесса теплопереноса, так и измерений в различных типах сенсоров нестационарного теплового потока; рассмотрены алгоритмы программ для решения прямых и обратных задач теплопроводности. Для получения оценок значений теплового потока разработаны программы двух модификаций ФК, которые позволили оценить поток в реальном времени.

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

0

т,с

0

Литература

1. Пилипенко Н.В. Методы параметрической идентификации в нестационарной теплометрии (ч. 1) // Изв. вузов. Приборостроение. - 2003. - № 8. - С. 50-54.

2. Пилипенко Н.В. Методы параметрической идентификации в нестационарной теплометрии (ч. 2) // Изв. вузов. Приборостроение. - 2003. - № 10. - С. 67-71.

3. Pilipenko N. Parametrical Identification of Differential-difference Heat Transfer Models in Non-stationary Thermal Measurements // Heat Transfer Research. - 2008. - V. 39. - № 4. - Р. 311 -315.

Кириллов Кирилл Валерьевич - Санкт-Петербургский государственный университет информационных

технологий, механики и оптики, аспирант, kirill.kirillov@gmail.com Пилипенко Николай Васильевич - Санкт-Петербургский государственный университет информационных

технологий, механики и оптики, доктор технических наук, профессор, pilipenko38@mail.ru