CC BY
14
УДК 681.2; 536.083
О ДОСТОВЕРНОСТИ НЕКОТОРЫХ МАТЕМАТИЧЕСКИХ МОДЕЛЕЙ НА БАЗЕ ВТОРОЙ ИНТЕГРАЛЬНОЙ ФОРМЫ УРАВНЕНИЯ ТЕПЛОПРОВОДНОСТИ
Рассмотрен общий вид второй интегральной формы уравнения теплопроводности применительно к наиболее распространенным моделям объектов теплофизических измерений: прямой цилиндр, параллелепипед и полуограниченное тело, из которого получены интегральные уравнения для указанных моделей. Приводятся результаты проверки их достоверности на различных непрерывных и дискретных тепловых моделях.
Ключевые слова: интегральная форма, уравнение теплопроводности, теплопроводность, объемная теплоемкость, конвективный теплообмен, кондуктивная, аккумуляционная составляющие.
Термины «интегральная форма уравнения теплопроводности (ИФУТ)», и «уравнение теплопроводности в интегральной форме» были применены академиками А.Н.Тихоновым и А. А. Самарским в работе [1] для названия соответствующих уравнений, описывающих баланс тепла за промежуток времени на некотором отрезке координаты или для объема. Из интегральной формы можно получить дифференциальную, как это показано в [1], которая представляет общеизвестное дифференциальное уравнение теплопроводности.
Интегральная форма представления уравнений, описывающих различные физические процессы, нашла применение при разработке интегро-интерполяционного метода построения консервативных разностных схем [2] для численного решения определенных задач математической физики.
Кроме того, ИФУТ используется для решения некоторых измерительных задач, например, на ее основе получена функции измерения теплоемкости [3]:
где 0(т), t2(t) - температура на противоположных гранях образца площадью S; m - масса образца; ту- время, начиная с которого достигается приемлемая точность определения интеграла по приближенной формуле; а так же функция измерения плотности теплового потока, при условии теплоизоляции противоположной поверхности, которая применяется при разработке энтальпийного тепломера [3]: q1 (т) = CLt'T(x), где t'T(x) - скорость
изменения средней на участке [0, L] температуры.
36
Ю.И. Азима, С.И. Сидельников
У
2S J [qi( т) - q2( т)]^т
У
Более перспективным, для решения задач измерения теплофизических свойств веществ, материалов и изделий, является вторая интегральная форма уравнения теплопроводности (ИФУТ-2) [4]. Она представляет собой баланс средних на интервале выбранной координаты количеств тепла: прошедшего через границы объекта, прошедшего вследствие теплопроводности и аккумулированном за счет теплоемкости. Для наиболее применяемых моделей объектов теплофизических измерений: прямой цилиндр, параллелепипед и полуограниченное тело, ИФУТ-2 для составляющей вектора теплового потока по выбранной координате % на интервале [b, B] может быть представлена в следующем виде:
B dXi ГГ^с , оЛ. -R^(S)Л. _П ч? dXi _
S (х i)
JJQ(X, T)dS + S(хi _ e, т)Q(sЧ(хi _ в,T)J
в s (хi )(г) в
_ х J [ (S % т)-1(s )(? т )]dc + С J Sd(Xx_ J s (х i )t(s )(х г-, т) dx i,
(1)
где х - вектор координат в прямоугольной, цилиндрической или сферической системе: х1= х, х2= у, Хз= z - для прямоугольных координат, х1= г, х2= z-для цилиндрических, х1= г- для сферических; b, В - границы интервала на выбранной оси координат О%-; S(% = b), S(%) - площади поверхностей, перпендикулярных составляющей вектора теплового потока по координате % и проходящих соответственно через переменную точку % и левую границу интервала % = b; Г- поверхность, перпендикулярная поверхностям площадью S(% = b), S(%); JJQ(г) (Xt, т)dГ - количество тепла, проходящего че-
(Г)
рез поверхность Г; t(s)(%, т) - средняя температура на площади S(%); 1, С -теплопроводность и объемная теплоемкость.
Или в сжатом виде:
1Q<%)_ Q
^, т)+ Qx(Lxi }(С, т),
X
(2)
где Q^' >(Х, т); Q,^1
(С, т)- среднее на интервале [b, В] количество тепла,
соответственно, проходящего через слой толщиной Lх _ B - в в направление оси координат О%-, и аккумулированного в нем; ^ Qx(r) (т) - среднее
на интервале [b, В] количество тепла, проходящего через граничную поверхность слоя толщиной LX , перпендикулярного составляющей теплово-
х i
го потока по выбранной координате.
Данная форма представления уравнения теплопроводности отличается по своей структуре от уравнения баланса тепла присутствием кондук-тивной составляющей и наличием одного, интегрального по времени, теплового потока Q(s)(хi _ в, т), который проходит через одну из границ, на-
пример b, и совпадает по направлению с выбранной составляющей теплопередачи в объекте исследования по оси координат О%-.
Такая структура уравнения теплопроводности (1) открывает перспективу идентификации коэффициентов 1 и С, определяющих теплофизические свойства, а также теплового потока или количества тепла на выбранной границе. Это, в свою очередь, позволяет проектировать новые методы измерения теплофизических величин.
На базе интегрального уравнения (1) была получена интегральная форма уравнения теплопроводности для различных моделей объектов теплофизических измерений. Данные уравнения являются оригинальными математическими выражениями, и корректность их получения может вызвать определенные вопросы, обусловленных их существенным отличием от привычных и понятных решений краевых задач теплопроводности. Поэтому имеет смысл показать достоверность полученных интегральных уравнений на различных тепловых моделях, что будет наглядным, убедительным и понятным, чем доказывать корректность тех или иных преобразований. После анализа данных подтверждающих фактов останется признать правильность выполненных преобразований.
Суть исследований заключалась в том, что для выбранной модели вычисляются левая и правая части уравнения баланс средних количеств тепла (2), и путем их сравнения определяется величина ошибки, по которой можно судить о достоверности исследуемого уравнения. Для выполнения этой процедуры использовалось, предварительно полученное для данной тепловой модели, аналитическое или численное решение краевой задачи теплопроводности. Вследствие того, что данное решение удовлетворяет дифференциальному уравнению теплопроводности для непрерывной модели или системе алгебраических уравнений - для дискретной, то оно должно удовлетворять и ИФУТ-2, которая является одной из форм представления уравнения теплопроводности.
Для определения величины ошибки использовалось выражение:
5( т)
Е £Т(т)-of->(*, т)+дЬ >(с, т)
е дг)(т)
• 100
(3)
Необходимо отметить, что источниками ошибки 5 (t) тождественности (2), кроме недостоверности интегральной формы, вследствие некорректности выполненных преобразований, являются: ограниченность функционального ряда, описывающего решение краевой задачи, величина шага дискретности модели по координате и времени, приближенное определение интегралов по координате и времени, значения параметров, входящих в критерий Fo. Будем считать, что погрешность порядка 5=(10'4^1)% является признаком того, что исследуемая ИФУТ корректно получена и является достоверной, а полученная ошибка тождественности обусловлена вышеперечисленными причинами.
Для упрощения исследований были выбраны одномерные тепловые модели, для которых известны аналитические решения краевой задачи теплопроводности. Для двухмерных объектов были разработаны дискретные модели в среде Mathcad, представляющие систему консервативных разностных уравнений, полученных интегро-интерполяционным методом [2]. Для их решения применялась встроенная функция Radau. Определение температуры в заданных точках дискретной модели, объекта как непрерывной функции времени, осуществлялось с помощью встроенной функции csplaine, выполняющей кубическую сплайн интерполяцию по дискретным временным точкам. Шаг дискретизации по времени составлял 0,1с.
Сначала приведем результаты исследований одномерных моделей в прямоугольной, цилиндрической и сферической системах координат. В прямоугольных координатах ИФУТ-2 для интервала [b, В] объекта в виде одномерной пластины имеет вид:
q(pт) 0 Q(b, т) 0
—J [t (в, т) -1 (В, т)]<2т +— J J t (х, T)dxdx
В x
X
L
, т) -1(В, T)]dT-J J t(x, T)dxdx
в в
в x
0,
0,
BB
т
0
т
0
(4)
(5)
где Q(x, т) 0 = -X J dz (х = в, В); L = В- в.
0 dx
Проверка достоверности проводилась на тепловой модели в виде пластины со следующими краевыми условиями:
t(x,0) = 0, - X dt(0, т) = q, t(h, т) = 0. d x
Для нее температурное поле описывается следующим выражением
[3]:
t( x, т) = ^~ X
qh (1 xX - h
V
2q v ( — > cos
X n=
x X
P nj V h)
exp
(* n021) (1
V
h
0
x
h )
Л f cos
x
P n J
V h)
dx , (6)
где pn = 0,5n(2n -1), h - толщина пластины, 1, а - теплопроводность и температуропроводность.
Для определения погрешности 5(т) выполнения интегральной формы (4) и (5) подставим решение (6 ) в левую и правую части данных выражений. Полученные погрешности при n =100 для различных: t, 1, а, в,Be [0, h], попадают в нижеприведенные границы адекватности. Иллюстрация данного исследования представлена на рис.1.
39
Рис. 1. Результаты исследований ИФУТ на одномерной модели пластины толщиной L=2,5 мм из материала с ТФС: X = 0,2 Вт/(мК), а =10 м /с; а - относительно количества тепла Q(t) (Дж/м ), поступившего в пластину через левую границу x=0 (уравнение (4)): Qx(t), Qc(t)(Дж/м ) - кондуктивная и аккумуляционная составляющие ИФУТ, б - относительно количества тепла QL(t) (Дж/м ), поступившего в пластину в правую границу x=L (уравнение (5)): QL1(t), QLc(t)(Дж/м ) - кондуктивная и аккумуляционная составляющие ИФУТ, 8(t),Sl(t) (%) - относительная погрешность тождественности, соответственно уравнений (4), (5), т(с) - время
В случае одномерного полого цилиндра с внутренним и внешним радиусами r0 и R и действующим источником тепла плотностью q(r0t) ИФУТ-2 описывается следующим выражением:
У
Q(r0, т)|T = — j [t(r, т) -1( R, т)] d
C R dx
L
т+L j
r
r0 r0
/
j rt(r, T)d \
T
0 ,
T
0
(7)
Q(R, t)T = L* j [t(r,, т) -1(R, т)] dT - С j — j rt(r, T)d r
L 0 L R r R
ro r
C 0 dx
T
0,
(8)
где
L = r0 j r'ldx = r0 ln
R
, l' = R in
R
- приведенная толщина цилиндриче-
r
r
0
0
ской стенки со стороны, соответственно, внутреннего r0 и внешнего радиуса R; Q(r0, t), Q(R, t) [Дж/м ]- количество тепла, проходящего через внутреннюю и внешнюю границу полого цилиндра.
40
Достоверность уравнений (7) (8) проверялась на тепловой модели одномерного полого цилиндра со следующими краевыми условиями:
t(x,0) = 0, -Xdt(r)’T) = q, t(R, т) = 0. о x
Решение данной краевой задачи имеет вид [5]:
t( x, т) = q X
r0ln
fR Л
+ П
n=i
V r J
, e.J „г,)J02(Rvn)[J0(rv„)i1 (r0Vn)-Y0(rVn)Ji(r0Vn)]
+л^ехр------------------—T(----) 2(a v,---------
Vn [Ji (r0 Vn ) J0 (RVn )]
где J0 (x), Ji (x), Y0 (x), Yi (x) - функции Бесселя первого и второго рода соответственно нулевого и первого порядков; vn - положительные корни уравнения: Ji (r0 v)Y0 (Rv) -Y i (r0v)J0 (Rv) = 0 .
Вычислительный эксперимент проводился при следующих исходных данных: n=30, r0=5h, R=r0+40h (h=0,i25, 0,25, 0,5 мм), l и а - соответствуют ТФС низкотеплопроводных материалов; Т=50^200 с. Очевидно, что точность вычисления составляющих уравнений (7) и (8) определяется
аТ
количеством n членов функционального ряда. Для n=30 и Fo = — > 0,i по-
x
грешность выполнения данных уравнений не превышает 0,05%, что также подтверждает их достоверность.
Для точечного источника, действующего в полуограниченной среде, ИФУТ-2 имеет вид:
0(т) = т* J [t(P, т) - t(R т)]йТ + J Ч J r2t(r, т)dr,
1 т 1 р r 0
(9)
где
I* =— fr2dx= — 2п^ 2п
(
i i
л
Vp r
Q(t) - количество тепла, выделяемого в
точечном источнике за время т; р е [0,R].
Точное решение уравнения теплопроводности для точечного источника постоянной тепловой мощности д(Вт), действующего в полуограниченной среде, при нулевых начальных условиях представляет собой следующее выражение [6]:
t(r, т) = q erfc 2nXr
f
r
Результаты вычислений погрешности 5(т) представлены на рис.2. Они также показывают достоверность уравнения (9).
4i
Рис. 2. Результаты исследований ИФУТ для полуограниченного тела с точечным источником: Qx(t), - кондуктивная
и аккумуляционная составляющие ИФУТ, Q(t) (Дж) - количество тепла, выделившегося в точечном источнике, t (с) - время
Рассмотрим двухмерные модели объекта в виде пластины и цилиндра. Прямоугольник имеет высоту Ly и образован прямыми: у = 0, y = Ly, x = 0, x = Lx . ИФУТ-2 для составляющей вектора теплового потока по оси х и по оси у, соответственно, имеют вид:
Lx x Lx x
LxQ
o[y V т)+2
L
J dxJ Qy (x, 0, т)dx + J dxJ Qy (x, Ly, т^)dx
0 0
= X J t L )(0, т) - t L)(Lx, т) dт + C J dx J t^Ly \ x, T)dx,
00
7(by),
(10)
L у Q У )(0, т) + -L-
L
у у
00
L
у у
J dy J Qy(0, у, т )dy + J dy J Qy (Lx, у, т )dy
0 0
0 0
(11)
т Г 1 Ly у
X J [t{Lx) (0, т) - t{Lx) (Ly, т)Jdт + C J dy J tL) (у, т)dy.
где Q{4 >(0, т) = Ll J dT J dy = -xJ &
dt L >(0, т)
ёт - среднее по высоте
у т 0
Ly прямоугольника количество тепла, поступившего в объект через границу х=0 за время т; Q (x, у, т) = -Х JЭ t (x,у’т) dx - количество тепла, проте-
0 Э у
кающего в направление оси Оу через сечения с координатами у= 0 или
_(т ) 1 Ly
у = Ly; Г у (x, т) = — J t(x, у, i)dy- средняя по высоте Ly температура
Ly 0
границ прямоугольника с координатами x = 0, x = Lx.
Lx x
0
0
0
Дискретная модель двумерной пластины состоит из n*m элементарных ячеек, имеющих форму прямоугольника hx х hy, длиной hx и высотой hy где n - количество слоев, уложенных по толщине пластины (Ly = nhy); m - количество ячеек уложенных в каждом слое по длине пластины (Lx = mhx).
Исследования выполнялась на дискретной модели размером 30hx □ 10hy со следующими краевыми условиями:
t( x, у, о) = 0,
q(0, у, т) = -a t(0, у, т), q(Lx, у, т)= -аtL, у, т), q(x, Ly, т)= -а t(x, Ly, т),
q (x, 0, т ) =
0т
~l Vx^" у ■
x e [0, hx ],
а t(x,0, т) x e [hx, Lx 1
где 0 - скорость изменения теплового потока; а- коэффициент теплоотда-
чи.
Для данных краевых условий уравнение баланса (10) на площади
Lx х Ly = 20,5hx х 10hy имеет следующий вид:
т
-(Ly N a
( L
- at[Ly>(0,т)-
LyLx
Lx x Lx x
I 11(x,0, т)dxdx + 111(x, Ly, т)dxdx
0 0
dт +
0т2 hx (Lx - 0,5 hx) X t Г-L ).A . -(Ly )/r J, CLf.f-(Ly), w
+------ту—-------= — J [t y (0, т) - r y!(Lx, т)Jdт + — J dx J Г y’(x, т)dx .
2 LyLx
L
L
x 0 x 0 0
Вычисление двукратных интегралов по координате выполнялось по квадратурной формуле на базе интерполяционного многочлена Лагранжа [7]:
Xmdx x
J -у J xst(x)dx = S Pkt(xkX
x
x0 x0
k=0
где s= 0 - для прямоугольных координат, s= 1- для цилиндрических; pk -весовые коэффициенты, определяемые путем интегрирования многочлена Лагранжа степени K-1, совпадающего с функцией распределения температуры в Kточках:
Pk=
x
s x - x
xm i x K
j — J n
l xs»0 '=2 ^ xk - x
dx
xm dx x K-1( s x - x..
\
J dx J n
x
x0 x x0 i=1 V xk xi J
x m x k -1
J -j J П
x
s x - xi
dx
K
если k = 1,
если k = K,
n
x0 x x0 i=1 V xk xi Ji=k+1
x
s x - xi
V xk xi
dx если k Ф1, k Ф K
где x = x0 +(i - 0,5) h, (i = 1, K) - координаты точек, в которых рассчитывается температура.
0
V hxhx
m
Определение однократных интегралов по координате выполнялось по формуле прямоугольников с шагом соответственно hx или hy. Результаты исследований представлены на рис. 3. При расчетах были приняты следующие исходные данные: 1=0,5 Вт/м-К, а=5 10-7м2/с, hx= hy= 0,2 мм, длина пластины L=30hx, толщина H=\0hy, коэффициент теплоотдачи a=50 Вт/м2 ■ К,
Рис. 3. Результаты исследований ИФУТ: а - по координате х, б - по координате у, на дискретной модели двухмерной пластины: Qx(t), Qc(t), Qv(t) (Дж/м ) - кондуктивная и аккумуляционная составляющие ИФУТ;
Qv(t),Qk(t) (Дж/м ) - составляющие ИФУТ, определяющие конвективный теплообмен пластины с окружающей средой и действие источника тепла; 5 (%)- погрешность тождественности
ИФУТ, t (с)- время
Круглый прямой цилиндр имеет радиус основания R и высоту Lz. ИФУТ-2 для составляющей вектора теплового потока по оси Oz:
2 Lz z
LQy\z = 0,T)+ - jdz jQr (R,z,т)dz = XJ[t{s-)(0,т) -1{S-)(Lz,т)]dr +
0 0 t (12)
Lz z
C j dz j t(Sr) (z, т )dz,
+ i
00
где QzSr)(z = 0, т) = —22X j dr j r ^t(r’z’T) dr - среднее на площади
D2 J J 7) t
t 0
основания с координатой z=0 количество тепла, поступившего в цилиндр;
2
t (Sr )(zi ) = — j rt (r, z )dr - средняя температура основания z=0 или z=Lz;
D2 J
R
R 2
Qr (R,z, Т ) = -Х j t
dt (R, z, т) Э r
—т
количество тепла,
поступившего по
координате r через боковую поверхность цилиндра; Sr - площадь основания цилиндра.
ИФУТ-2 для составляющей вектора теплового потока по оси r:
L
R drr
j — j r[Qz (r,0, t) + Qz (r, Lz, t)] dr
0 r 0
R d
= X j [t(Lz) (0, t) -t ^) (R, t)] dT + C j — j rt(Lz )(r, t) dr,
(13)
r
00
T
где Qz (r,0, t); Qz (r, Lz, t) - количество тепла, поступившего в цилиндр
1 Lz
через основания цилиндра; t(Lz) (r, т) = — j t(r, z, t)—z - средняя по высоте
Lz 0
цилиндра температура, как функция радиуса.
Для исследования достоверности интегральных уравнений (12) и (13) применялась дискретная модель цилиндра, состоящего из n*m элементарных ячеек, имеющих форму полого цилиндра с толщиной стенки h (нулевая ячейка представляет цилиндр радиуса h) и высотой h, где n -количество слоев по высоте цилиндра; m - количество ячеек в каждом слое по радиусу цилиндра. Были использованы следующие граничные условия:
t(r, z,0) = 0, q (r, z = H, T) = -a z t (r, H, t ), q(R, z, T) =-a rt (R, z, T),
q(t) = 0t r e [0, kh\
- a zt(r,0, t) r e[kh, R], n = 1, R / h
где 0 - скорость изменения теплового потока; az, ar - коэффициент теплоотдачи, соответственно по оси Oz и Or, h - высота и толщина стенки цилиндрической трубки (ячейки) с внутренним радиусом r=ih о = 1, n -1), составляющих дискретную модель цилиндра; m, n - количество ячеек соответственно по высоте и радиусу цилиндра; k - количество ячеек, на которые воздействует тепловой поток q(t).
Результаты вычисления ошибки 5(т) для уравнений (12), (13) представлены на рис.4. При расчетах были приняты следующие данные: 1=0,2 Вт/м-К, а=1,2 10-7м2/с, h=0,2 мм, радиус цилиндра R=20h - для уравнения (12) и R=30h - для уравнения (13), высота цилиндра H=10h, радиус нагревателя r=3h, коэффициент теплоотдачи, соответственно в направлении координаты r: ar= 200 Вт/ м ■ К, и z: az= 20 Вт/м ■ К.
45
q(r, z = 0, т)
1.5 “i 1 1 1 1
1.25 t 1 ъ 1
1 t -
5z(t) ■ ■ ■ ■
У*'* 0.75 ог 1 1 " '■ \ "
0.5
ъ
0.25 %*****_
0 1 \ 4гг**1+ *** Т***ч
О---------1-------1--------1-------1------^
О 20 40 60 80 100
т
Рис. 4. Относительная погрешность dz, 8r (%) тождественности ИФУТ для цилиндра, соответственно по координате z и радиусу r,
как функция времени т(с)
Таким образом, проведенные исследования подтверждают достоверность всех приведенных уравнений, представляющих собой вторую интегральную форму уравнения теплопроводности. Практическая ценность данных уравнений заключается в том, что они создают теоретическую базу для проектирования новых нестационарных методов измерения ТФС и тепловых потоков, обеспечивающих приборную реализацию с улучшенными метрологическими и эксплуатационными характеристиками по сравнению с существующими.
Кроме того, данный метод получения второй интегральной формы уравнения теплопроводности может быть применен к уравнению диффузии, дифференциальная форма которого имеет аналогичную структуру с уравнением теплопроводности [1]. Это также создает перспективы разработки новых методов измерения коэффициентов диффузии и пористости.
Список литературы
1. Тихонов А.Н., Самарский А. А. Уравнения математической физики. Москва: Наука, Изд. МГУ, 1977. 799 с.
2. Самарский А.А.Теория разностных схем. М.: Наука, 1977. 656 с
3. Теплофизические измерения и приборы / под ред. Е.С. Платуно-ва. Л.: Машиностроение, 1986. 256 с.
4. Азима Ю.И. Представление уравнения теплопроводности в интегральной форме для задач измерения теплофизических величин // Известия Тульского государственного университета. Технические науки. Тула: Изд-во Тул.ГУ. 2014. Вып. 4. С. 123-134.
5. Карслоу Г. Теплопроводность твердых тел / Г. Карслоу, Д. Егер. М.: Наука, 1964. 487 с.
6. Филиппов П. И. Приложение теории теплопроводности к теплофизическим измерениям.; под ред. В. П. Ушкалова. Новосибирск, Наука, 1973. 61 с.
7. Никольский С.М. Квадратурные формулы. М: Наука, 1979. 256 с.
Азима Юрий Иванович, канд. техн. наук, доц., iuia@yandex.ru, Россия, Новомосковск, Филиал Российского химико-технологического университета им. Д.И. Менделеева,
Сидельников Сергей Иванович, канд. техн. наук, доц., sidsergl 1@mail.ru, Россия, Новомосковск, Филиал Российского химико-технологического университета им. Д. И. Менделеева
ABO UT THE RELIABILITY OF SOME MA THEMA TICAL MODELS BASED ON THE SECOND INTEGRAL OF THE FORM THE HEAT EQUATION
Y.I. Azima
The General form of the second integral form for the heat equation with respect to the most common models of objects thermophysical measurements: straight cylinder, parallelepiped and semi-bounded body from which the received integrated equations for these models. The results of validation for various discrete and continuous thermal models
Key words: integral form, the equation of heat conduction, thermal conductivity, volumetric heat capacity, convective heat transfer, conduction, storage components
Azima Yuriy Ivanovich, candidate of technical sciences, docent, iuia@yandex. ru Russia, Novomoskovsk, Novomoskovsk’s Institute (subdivision) of the Mendeleyev Russian Chemical-Technological University
Sidelnikov Sergey Ivanovich, candidate of technical sciences, docent, sidsergl 1 @mail. ru, Russia, Novomoskovsk, Novomoskovsk’s Institute (subdivision) of the Mendeleyev Russian Chemical-Technological University
