CC BY
43
НАУЧНОЕ ИЗДАНИЕ МГТУ ИМ. Н. Э. БАУМАНА
НАУКА и ОБРАЗОВАНИЕ
Эл № ФС77 - 48211. Государственная регистрация №0421200025. ISSN 1994-0408
электронный научно-технический журнал
Сравнительный анализ оценок коэффициента теплопроводности композита с шаровыми включениями # 07, июль 2013 Б01:10.7463/0713.0569319
Зарубин В. С., Кувыркин Г. Н., Савельева И. Ю.
УДК 541.124
Россия, МГТУ им. Н.Э. Баумана fn2@bmstu.ru
Введение
Существуют различные подходы к построению математических моделей процесса переноса тепловой энергии в композите с шаровыми включениями, предназначенных для получения оценок эффективного коэффициента теплопроводности неоднородного твердого тела, которым является такой композит [1, 2, 3, 4, 5, 6, 7, 8]. Проведем сравнительный количественный анализ этих оценок применительно к композиту, состоящему из изотропных по отношению к свойству теплопроводности матрицы и шаровых включений с коэффициентами теплопроводности Ат и Ао соответственно. В силу электротепловой аналогии [9] такие оценки могут быть использованы для расчета эффективных значений электропроводности, диэлектрической и магнитной проницаемости композита рассматриваемой структуры. Именно расчету электропроводности неоднородной среды с шаровыми включениями была посвящена первая известная в этом направлении работа, опубликованная К. Максвеллом [10] еще в 1873 году.
1. Основные соотношения
Выделим в композите некоторый объем V*, осреднение в пределах которого свойств матрицы и включений дает представление о характеристиках композита в целом. Такой объем принято называть представительным. Точка М Е V* может быть расположена как в матрице, занимающей объем У^* С V*,так и в одном из включений, занимающих объем V* = V* \ ут. Отношение У0*/V* определяет объемную концентрацию Су включений.
Связь между векторами g(M) и q(M) соответственно градиента температурного поля и плотности теплового потока устанавливают вытекающие из гипотезы Био — Фурье [3] соотношения
q(M ) = -Л(М^(М), g(M) = ^(М )/Л(М),
где Л(М) = Л0, если М е V)*, и Л(М) = Лт, если М е ^• После осреднения этих соотношений по объему V* запишем
(1)
(я) = V-*/ Ч(М) ¿V (М) = - V-*! Л(М)ё(М) ¿V (М) = (Л§);
V * V *
(ё) = V-*/ 8(М) ¿V(М) = -V-! ЛМ) ¿V(М) = (я/Л).
V * V *
Введем локальные отклонения (флуктуации) от средних значений
я'(М ) = я(М) -(я) и ё'(М ) = ё(М) -(ё), М е V*.
В линейном приближении примем, что я'(М) = 0(М)(я) и ё'(М) = б(М)(ё), М е V*, где 0(М) и б(М) —тензоры второго ранга, зависящие от положения точки М е V*. Учитывая, что средние по объему V* флуктуации равны нулю, равенства (1) можно представить в виде
(я) = ((Л) + (Л0))(ё) = ЛКё), (ё) = ((1/Л) + (П/Л))(я) = (я)/Л2, (2)
где Л1 = ((Л) + (Лб)) и Л2 = 1/((1/Л) + (б/Л)) — оценки значения Л* эффективного коэффициента теплопроводности композита, связывающего осредненные по его представительному объему V* векторы (я) и (ё) плотности теплового потока и градиента температуры соответственно.
Из равенств (2) следует, что построенные оценки значения Л* в общем случае не совпадают со средними значениями (Л) и 1 /(1/Л). Если предположить отсутствие флуктуаций либо вектора градиента температуры, либо вектора плотности теплового потока, то для таких частных случаев получим соответствующие оценки среднего значения Л коэффициента теплопроводности композита в виде
Л = (Л) = + Лт(1 - ^Л2 = (1/Л) = ^/Ло + (1 - С)/Лт. (3)
2. Двусторонние оценки
Определяемые равенствами (3) оценки ограничивают возможные значения эффективного коэффициента теплопроводности композита сверху и снизу, т.е. являются двусторонними оценками значения Л*. Сначала, использовав эти равенства, установим, что Л1 ^ Л2, или
ч2 , (л г< \2 , п (Л П \( Л0 Л
С2 + (1 - ^)2 + ^(1 - ^и + тЧ ^ 1.
V Лт Л0 /
Действительно, неравенство справедливо, так как Л0/Лт + Лт/Л0 ^ 2.
Оценки Л1 и Л2 можно получить, если воспользоваться двойственной вариационной формулировкой стационарной задачи теплопроводности в неоднородном твердом теле [11, 12], включающей минимизируемый функционал
[ё] = 1/ Л(М)(ё(М))2 ¿V(М) (4)
V *
и максимизируемый функционал
IМ = -1 / ^ (М) - / /т (Р )я(Р) ■ п(Р) ¿Б (Р), (5)
V * ят
где п(Р) — единичный вектор внешней нормали в точках Р Е Бт участков Бт С Б поверхности 5, ограничивающей объем V* неоднородного тела, на которых задано распределение температуры /т(Р). Функционал (4) допустимо рассматривать на таких распределениях Т(М) (М € V*) температуры в этом теле, которым соответствуют кусочно непрерывные распределения вектора градиента g(M) = — Л(М) УТ(М) (V — дифференциальный оператор Гамильтона) и которые удовлетворяют заданному граничному условию Т (Р) = /т (Р), Р € Бт. При этом участки Б \ Бт поверхности тела должны быть идеально теплоизолированы. Для функционала (5) допустимыми являются непрерывные распределения вектора я(М), М € V, удовлетворяющие дополнительным условиям V ■ я = 0 в V* и я ■ п = 0 на Б \ Бт [11].
Из экстремальных свойств функционалов (4) и (5) следует цепочка неравенств
3И ^ 3^ I[я], (6)
где g*(M), М € V, — истинное распределение градиента температурного поля в неоднородном теле. На этом распределении функционал (4) достигает своего наименьшего значения, равного [12]
3= — 2 / /т(Р)Л(Р^(Р) ■ п(Р) ¿Б(Р). (7)
Ят
Пусть часть Бт поверхности Б неоднородного тела объемом V* состоит из двух не граничащих между собой изотермических участков БТ и БТ с заданными значениями температур Т' и Т" соответственно. Эти участки отделены друг от друга идеально теплоизолированными участками поверхности Б (рис. 1). В этом случае вместо равенства (7) получим
23 ^*] = Т'д' + Т''д'', (8)
где ^' и Q'' — суммарные тепловые потоки, поступающие в область V через участки БТ и Бт'' ее поверхности соответственно, причем
Q' = — У Л(Р'^*(Р') ■ п(Р') ¿Б(Р'), Q'' = — I Л(Р'')g*(P'') ■ п(Р'') ¿Б(Р'').
В силу закона сохранения тепловой энергии при установившемся процессе теплопроводности в неоднородном теле, занимающем область V, имеем Q' + Q'' = 0. Поэтому, учитывая равенство (8), получаем Q' = —Q'' = 23[Т*]/(Т' — Т''). Отсюда находим термическое сопротивление рассматриваемого тела
_ = Т' — Т'' = (Т' — Т'')2 (9)
Q' 23 [g*] (9)
между изотермическими участками БТ и БТ поверхности Б. Цепочка неравенств (6) и равенство (9) позволяют построить двусторонние оценки значения термического сопротивления
тела в виде
(Г' - г'')2 , R , (Г' - Г'')2 ,]0.
где Ji[g] и Ii[q] — значения функционалов (4) и (5) на каких-либо допустимых распределениях векторов соответственно градиента температурного поля и плотности теплового потока.
Если область V*, соответствующую представительному объему композита, выбрать в достаточно простой форме, то после замены в этом объеме рассматриваемого композита однородным материалом с искомым значением Л* эффективного коэффициента теплопроводности можно однозначно выразить через Л* величину RT в соотношении (10) и тем самым получить двусторонние оценки значения Л*. Выберем представительный объем композита с объемной концентрацией CV включений в виде прямого цилиндра высотой H с идеально теплоизолированной боковой поверхностью и площадью S0 каждого из оснований, т.е. в данном случае V* = HS0. Примем температуру одного из оснований за нуль отсчета, а температуру другого положим равной T0.
Для функционала (4) простейшим допустимым распределением вектора g в таком объеме является постоянное значение g = T0/H единственной составляющей этого вектора, перпендикулярной основаниям цилиндра. Тогда из формулы (4) следует
Ji[g] = (Гц/2)(So/H)(ЛоСу + Лт(1 - CV)). (]])
Для максимизируемого функционала (5) в качестве допустимого распределения вектора q можно принять значение q единственной составляющей этого вектора, также перпендикулярной основаниям цилиндра. В этом случае формула (5) примет вид
q2 fCy 1 - C
11 [9] = -ТЯ5Ч ^ + ^Г^ ) - ^Я). (12)
2 Л0 Лт
Значение 9 найдем из необходимого условия (9) /¿9 = 0 максимума функционала, что дает 9 = -(Т0/#)/(Су/Л0 + (1 - )/Лт). После подстановки этого значения в формулу (12)
получим
11 [9]
(т02/2)(5о/Я)
(13)
Су/Ао + (1 - Су)/Лт
При замене композита в представительном объеме однородным материалом с коэффициентом теплопроводности Л* для термического сопротивления цилиндра находим Кт = = Я/(Л* 50) и с учетом неравенств (10) и формул (3), (11) и (13) получаем двусторонние оценки
А1 — АоСу + Лт (1 - Су) ^ Л* ^
1
Л2
(14)
Су/ло + (1 - Су)/Ат
На рис. 2 представлены построенные по формулам (3) графики зависимостей А1 = А1 /Лт (штрихпунктирные линии) и Л2 = Л2 /Лт (штриховые линии) от Су при различных значениях Л = А0/Ат. При малом отличии значения Л от единицы разность ДЛ = Л1 — Л2 сравнительно мала во всем промежутке изменения Су. По мере отклонения значения Л от единицы несмотря на совпадение оценок при Су = 0 и Су = 1 эта разность для промежуточных значений Су становится значительной. При фиксированном значении Л максимальному значению (ДА)тах = (1 - ^)2 соответствует объемная концентрация включений
1
(Су )тах —
1 + 1/^Л'
^Ф 2
4
%
\
/ %2, %
/
Х=5
Х=2
к=0,5 1=0,2
0,2
0,4
0,6
0,8
Рис. 2. Графики зависимостей безразмерных двусторонних оценок эффективного коэффициента теплопроводности композита от объемной концентрации включений при различных значениях Л
Если в качестве первого приближения к искомому значению Л* = Л* /Лт принять
А + Л2
Л
2
(15)
то возможная относительная погрешность при (Су )тах не превысит значения п(Л)
ДЛт
2Л*
причем п* (Л) = п*(1/Л). На рис. 3 в полулогарифмических координатах штриховой линией изображен график зависимости п(Л) при Л ^ 1.
Рис. 3. Графики зависимостей п(А) и п*(А)
Двусторонние оценки можно сблизить, если использовать вариационный принцип Ха-
шина — Штрикмана [9]: при Л < 1
А! = А! -
су(1 — су )(ат — ао)2 ао(1 - су) + ат(2 + су)
^ А* ^ А! —
СУ (1 — СУ )(Ат — А0)2
ао(3 — су) + АтСу
а2,
(16)
а при Л > 1
а:
А1 -
су(1 — су)(а-ат)2 ^ а* А0 (3 — су) + АтСу
> а! —
су(1 — Су )(Ат — А0)2 ао(1 — су ) + ат (2 + су)
л;
(17)
Эти оценки совпадают с результатами, полученными в рамках обобщенного сингулярного приближения теории случайных функций [6] и учитывающими в формулах (2) влияние флук-туаций в распределениях векторов градиента температурного поля и плотности теплового, т.е. а^ = А^, v = 1, 2. На рис. 2 по формулам (15) и (16) с учетом равенств (3) построены графики зависимостей А! = а!/ат (сплошные кривые) и А2 = а2/ат (пунктирные линии) от су при различных значениях Л.
Разность Да* = А! — а2 при промежуточных значениях су также возрастает по мере отклонения значения Л от единицы. При фиксированном значении Л наибольшего значения (ДА*)тах эта разность достигает при объемной концентрации
(су )т
2Л + Л2 — у/ Л(2 + Л)(1 + 2Л)
Л2 1
Для среднего арифметического значения А* = (а! + А2)/2 возможная относительная погрешность при (су)тах не превысит значения п*(л) = (Да*)тах/(2А*), причем п*(л) = п*(1/Л). График зависимости п*(л) при Л ^ 1 представлен на рис. 3 в полулогарифмических координатах пунктирной линией.
Необходимо отметить, что вычисленные по формуле (15) значения а* выходят за пределы двусторонних оценок а! и А2 при Л < 1 по мере приближения су к нулю, а при Л > 1 по мере приближения су к единице.
3. Модели теплового взаимодействия включения и матрицы
Представленные выше двусторонние оценки определяют границы возможных значений Л* эффективного коэффициента теплопроводности композита без учета особенностей распределения температуры в матрице в непосредственной близости к включению. Учет влияния таких особенностей возможен в рамках математической модели теплового взаимодействия включения и матрицы как в предположении идеального теплового контакта между ними [8], так и при наличии термического сопротивления на контактной поверхности [13]. При химическом взаимодействии материалов включения и матрицы между ними возможно образование промежуточного слоя с коэффициентом теплопроводности, отличиющимся от коэффициентов теплопроводности как включения, так и матрицы. В работах [14, 15] построены математические модели, учитывающие наличие промежуточного слоя с постоянным или непрерывно изменяющимся по его толщине коэффициентом теплопроводности.
Математическую модель переноса тепловой энергии в композите построим в предположении, что шаровые включения не контактируют между собой, т.е. отделены друг от друга слоем материала матрицы. Для учета взаимодействия между включениями и матрицей при идеальном тепловом контакте примем в качестве представительного элемента структуры композита составную шаровую частицу, состоящую из шарового включения с коэффициентом теплопроводности Л0 и радиусом R0 и окружающего его шарового слоя из материала матрицы с коэффициентом теплопроводности Лт и наружным радиусом Rm [8]. Эта частица в тепловом отношении взаимодействует с неограниченным массивом однородного материала, коэффициент теплопроводности Л которого подлежит определению как эффективная характеристика композита. Таким образом, модель композита содержит три фазы: включение, слой матрицы и неограниченный массив однородного материала.
Центр включения поместим в начале сферической системы координат. Примем, что на большом расстоянии r от начала координат задан вектор градиента температурного поля в однородном материале, имеющий модуль g и направленный по оси сферической системы координат, от которой происходит отсчет угловой координаты т.е. при r ^ то распределение температуры в этом материале описывает функция T^(r, = gr cos Эта функция удовлетворяет уравнению Лапласа, которое в сферических координатах имеет вид
В данном случае благодаря коллинеарности заданного вектора градиента температурного поля оси отсчета угловой координаты § распределение температуры симметрично относительно этой оси и не зависит от угловой координаты т.е. d2T/d^2 = 0.
По мере приближения к составной шаровой частице температурное поле в однородном материале претерпевает возмущение, описываемое также удовлетворяющим уравнению Лапласа дополнительным слагаемым [3] AT(r, §) = (B/r2) cos §, где B — подлежащий определению постоянный коэффициент. Таким образом, температурное поле в однородном ма-
1 d2T
r2 sin2 § d^2
0.
териале, удовлетворяющее заданному условию при r ^ то и уравнению Лапласа, описывает функция
T(r, tf) = T^(r, tf) + AT(r, tf) = (gr + B/r2) costf. (18)
Аналогичные зависимости описывают распределения температуры в шаровом включении
To(r, tf) = (Aor + Bo/r2) cos tf (19)
и в слое материала матрицы
Tm(r, tf) = (Amr + Bm/r2) COS tf, (20)
причем B0 = 0 в силу ограниченности температуры включения при r = 0.
Таким образом, в равенства (18)—(20) входят 4 неизвестных коэффициента B, A0, Am и Bm, которые необходимо найти из граничных условий на сферических поверхностях с радиусами Д0 и Rm. При r = R0 из условий непрерывности распределения температуры и радиальной составляющей вектора плотности теплового потока следует
дт
T0(R0,tf) = Tm (R0,tf) и A0 0
dT
m
Am й r=Ro dr
. (21)
r=Ro
дг
Отсюда с использованием равенств (19) и (20) при Во = 0 находим
Ао = Ат + Вт /до и ААо = (Ат - 2Вт/Д3), (22)
где А = Ао/Ат. Из аналогичных условий непрерывности распределения температуры и радиальной составляющей вектора плотности теплового потока на сферической аоверхности радиусом Дт с учетом формул (18) и (20) получим
Вт В , Вт ~ ( 2В \
Ат + дт = 9 + и Ат - 2Дт = Х{9 - , (23)
т т т т
где А = А/Ат.
Последовательным исключением неизвестных из равенств (22) и (23) находим
В 3 г 2 + А+(1 - А)(Во/Дт)3 (24)
—— = ¿9А^-=-^-=--9. (24)
дт (2А + 1)(2 + А) + (А - 1)(1 - А)(До/Дт)3
Замена составной шаровой частицы равновеликим шаром радиусом Дт с искомым коэффициентом теплопроводности А приведет к исчезновению возмущения температурного поля в окружающем ее однородном материале. Тогда в равенстве (18) следует положить ДТ (г, $) = 0, что равносильно условию В = 0, которое с учетом формулы (23) и равенства (До/Дт)3 = Су позволяет записать
А =2 + А- 2(1 -А)Су, (25)
2 +А + (1 - А) Су
Эта формула совпадает с известной формулой Максвелла, приведенной в [3] и полученной на основе более простой двухфазной модели, состоящей из включения в виде сплошного
шара и окружающего его материала матрицы. К виду (25) можно привести формулу, представленную в [1], и зависимость, полученную в [16] для эффективного значения магнитной проницаемости и приведенную в [7], но без упоминания о формуле Максвелла. Результаты расчета по формуле (25) совпадают со значениями А1 при Л ^ 1 и со значениями Л* при Л ^ 1, так как в этих случаях равенство для Л**, вытекающее из соотношения (16), и равенство для Л*, следующее из соотношения (17), можно преобразовать к виду (25).
В [2, 17, 18] со ссылкой на работу [19] приведена функциональная зависимость, которую в использованных выше обозначениях можно представить в виде
/ 1 ч 1/3 Г- Л
1 - = Ы Т-Т • (26)
где Л* = Л*/Лт, Л* — эффективный коэффициент теплопроводности композита, оцениваемый этой зависимостью. Равенство (26) удается преобразовать к неполному кубическому уравнению относительно величины (А*)1/3, имеющее при Су € [0; ] единственный действительный корень, который при Л ^ 0 и Л = 1 является положительным.
Шаровую форму включений можно рассматривать как некоторую статистически усредненную по отношению к включениям произвольной формы, но имеющим близкие размеры во всех трех направлениях. В [5] построены двусторонние оценки для эффективного коэффициента Лп теплопроводности композита с кубическими включениями на основе математической модели, в которой использован прием деления кубического представительного объема композита на фрагменты адиабатическими и изотермическими плоскостями. Полученные при этом оценки можно представить в виде
-<2/3
1 - (1 - Л)Су/3(1 - С1/3) ^ Лт ^ Л + (1 - Л)С1/3
Используем еще один подход к построению оценки эффективного коэффициента теплопроводности композита, имея в виду установленную выше линейность распределения температуры в шаровом включении с модулем А0 градиента. Если поместить такое включение с коэффициентом теплопроводности А0 и радиусом Я0 в начало сферической системы координат и окружить неограниченным массивом однородного материала с искомым коэффициентом теплопроводности Л, то при заданном на большом расстоянии г модуле д градиента температурного поля в этом материале вместо равенств (18) и (19) можно записать
Т*(г, 0) = (дг + В*/г2) сое0 и Т0*(г, 0) = А0г сое0. С учетом этих равенств из граничных условий вида (21) получим
а* = д + Вт/Я3 и АоА* = А(д - 2В*/Я03).
3 о
Отсюда находим модуль А* = -——— градиента температурного поля в шаровом включе-
2 + Ло / Л
нии.
= 1 - (1 - Л)С2/3 . Лп > Л + (1 - Л)(1 - Су/3)Су = ЛП (27) Л1 = :—-—„2/^—::тт::7 ^ т" ^-=—-—-= Л2. (27)
Представим матрицу композита в виде совокупности шаровых частиц различного радиуса, что позволит заполнить любой объем без пустот. Тогда в отдельно взятой частице матрицы при условиях теплового взаимодействия с однородным материалом, аналогичных
взаимодействию его с шаровым включением, возникнет распределение температуры с мол* 3д
дулем Ат = 2 + Лто/Л градиента.
Рассматривая величины АО — д и А*т — д как флуктуации модуля градиента в однородном материале, соответствующем композиту, и проводя осреднение этих флуктуаций по объему композита, получаем равенство Су(А0 — д) + (1 — Су)(А^ — д) = 0. Отсюда следует квадратное уравнение относительно искомого значения Л° эффективного коэффициента теплопроводности композита с такой моделью структуры. Решение этого уравнения можно представить в виде
~ Л° 2 — Л — 3(1 — Л)Су + л/(2 — Л — 3(1 — Л)Су )2 + 8Л
л° = - =-V.-. (28)
% т:
Если по объему композита осреднять флуктуации ЛоА0 — Лд и ЛтАт — Лд плотности теплового потока соответственно во включениях и частицах матрицы, то эта процедура также приводит к формуле (28).
При объемной концентрации включений, близкой к единице, использованная при получении формулы (25) модель структуры композита, предполагающая наличие шарового слоя матрицы, окружающего каждое включение, становится неадекватной в силу того, что некоторые включения будут контактировать между собой. В этом случае можно предположить, что частицы включений и матрицы меняются своими ролями, т.е. модель структуры композита представима в виде шаровых частиц из материала матрицы, окруженных шаровым слоем частиц из материала включений. Тогда после преобразований, аналогичных использованным при получении формулы (25), для эффективного коэффициента теплопроводности Л* композита с близким к единице значением Су получим формулу
Л = ^ = Л2Л + 1 — 2(Л — 1)(1 — Су)
Лт 2Л + 1 + (Л — 1)(1 — Су) . ( )
Результаты расчета по этой формуле совпадают со значениями Л2 при Л ^ 1 и со значениями Л1 при Л ^ 1, поскольку в этих случаях равенство для Л2, которое следует из соотношения (16), и равенство для Л1, следующее из соотношения (17), идентичны равенству (29).
4. Результаты расчетов
Для сравнения результатов расчетов по всем полученным формулам выберем сначала два фиксированных значения: Л = 0, 2 и обратное ему значение Л = 5. Для первого значения результаты представим в виде графиков зависимостей от Су отношения Л = Л/Лт с различными индексами, а для второго значения — в виде графиков зависимостей от Су отношения Л = Л/Л0 также с разными индексами. Это позволит представить все результаты расчетов графиками с ординатами, не превышающими единицу.
На рис. 4 все графики для значения Л = 0, 2 совокупности зависимостей от Су отношения Л* /Ат, обозначенного символом Л с различными индексами, имеют общие точки с абсциссами Су = 0 и Су = 1 и ординатами, равными соответственно Л = 1 и Л = 0, 2, а все графики для значения Л = 5 совокупности зависимостей отношения А*/Ао, обозначенного символом Л с разными индексами, также имеют общие точки с абсциссами Су = 0 и Су = 1 и ординатами, но равными соответственно Л = 0, 2 иЛ =1.
Рис. 4. Графики зависимостей безразмерных эффективных коэффициентов теплопроводности композита от объемной концентрации включений при значениях Л = 0, 2 и Л = 5
Из этого рисунка видно, что всюду в интервале (0; 1) изменения Су справедлива цепочка строгих неравенств А1 > Л° > А1 > Л° > Л* > Л2 и Л1 > А1 > Л° > Л* > Л° > > Л2. Таким образом, оценки, определяемые формулой (27), являются двусторонними по отношению к значению Л, которое следует из равенства (25). Вместе с тем величины А^ и Л°, определяемые формулами (26) и (28) соответственно, выходят за границы этих оценок, но остаются в пределах, устанавливаемых соотношениями (16) и (17). Графики для этих величин при малых значениях Су близки к графику, соответствующему формуле (25), но по мере роста значения Су они приближаются к графику, соответствующему формуле (29). Аналогичная закономерность справедлива и для графиков зависимостей и Л° от Су. Дело в том, что по мере приближения значения к единице материалы матрицы и включений меняются ролями, что отражено в модели структуры композита, использованной при выводе формулы (29). Эту особенность учитывает и модель структуры, использованная при получении соотношения (28).
Поскольку разность ДА* увеличивается по мере отклонения значения Л от единицы, проведем сравнение полученных зависимостей от объемной концентрации Су шаровых включений при достаточно малом значении Л = 0, 01 и при достаточно большом значении
Л = 100. Как и на рнс. 4, на рис. 5 все графики для значения Л = 0, 01 совокупности зависимостей Л = А*/Ат от Су имеют общие точки с абсциссами Су = 0 и Су = 1 и ординатами, равными соответственно А =1 и А = 0, 01, а все графики для значения Л = 100 совокупности зависимостей А = А* /Ао от Су также имеют общие точки с абсциссами Су = 0 и Су = 1 и ординатами, но равными соответственно Л = 0, 01 и Л =1. Обозначения графиков на этих рисунках совпадают, а закономерности во взаимном расположении графиков аналогичны, но при этом разности ординат для соответствующих пар двусторонних оценок на рис. 5 в силу большего различия между значениями коэффициентов теплопроводности матрицы и включений оказались больше, чем на рис. 4.
0 0,2 0,4 0,6 0,8 Су
Рис. 5. Графики зависимостей безразмерных эффективных коэффициентов теплопроводности композита от объемной концентрации включений при значениях Л = 0, 01 и Л = 100
При изменении Су в полуинтервале [0, 8; 1) масштаб рис. 5 не позволяет выявить справедливость строгих неравенств < А/Ат < и < А/А0 < Поэтому в табл. 1 представлены результаты расчетов по формулам (25) и (26). Эти результаты подтверждают выполнение этих неравенств.
Таблица 1
Су 0,80 0,85 0,90 0,95 0,99 0,999
Лп Л2 Х/Хт Л? А/Ао Л? 0,1475 0,1520 0,1564 0,1078 0,1142 0,1213 0,1122 0,1147 0,1171 0,1453 0,1514 0,1580 0,0775 0,0786 0,0797 0,2117 0,2172 0,2229 0,0434 0,0437 0,0440 0,3612 0,3649 0,3687 0,016644 0,016655 0,016666 0,746305 0,746917 0,747530 0,01066341 0,01066352 0,01066363 0,96740509 0,96741530 0,96742549
Из сравнения рис. 4 и 5 следует, что по мере отклонении параметра Л от единицы различие между результатами расчетов по рассматриваемым формулам наиболее существенно возрастает при промежуточных значениях Су объемной концентрации включений. При
этом увеличивается и разность двусторонних оценок значений эффективного коэффициента теплопроводности композита. Ясно, наиболее приемлемую расчетную зависимость можно выбрать лишь с учетом экспериментальной проверки, хотя до настоящего времени опубликованных результатов тщательно проведенных экспериментов еще недостаточно для того, чтобы сделать определенные выводы.
5. Сравнение с экспериментом
Рассмотрим путь сравнительного анализа расчетных зависимостей с использованием экспериментальных результатов, представленных в работах [17, 18]. В этих работах предполагалось наличие неидеального теплового контакта между включением и матрицей, характеризуемого коэффициентом а контактного теплообмена и учитываемого введением эффективного коэффициента теплопроводности включения
А* = = , (30)
А. + 1 Л/в + 1 ' '
аК0
где в = аЯ0/Ат. Действительно, если в формулу (25) вместо Л подставить отношение Л* = А*/Ат, то после преобразований можно получить формулу [13]
А = 2Л + в(2 + Л) - 2(в(1 - Л) + Л)Су
= 2Л + в(2 + Л) + (в(1 - Л) + Л)Су ,
которая непосредственно следует из рассмотрения математической модели теплового взаимодействия шарового включения и матрицы при наличии между ними неидеального теплового контакта.
В работе [18] приведены экспериментальные данные по эффективному коэффициенту А теплопроводности композита с алюминиевой матрицей и шаровыми включениями из карбида кремния при Су от 0,53 до 0,60. Включения принадлежали семи фракциям и имели средние значения диаметров 2Я0 от 8,9 мкм до 167 мкм (табл. 2).
Таблица 2
Су 0,58 0,58 0,60 0,59 0,58 0,55 0,53
2^0, мкм 167 86,4 56,8 37,1 23,4 16,9 8,9
Л, Вт/(мК) 221 209 203 204 194 193 154
При обработке этих экспериментальных данных для алюминиевой матрицы было принято значение Ат = 237 Вт/(м-К). Поскольку формулу (30) можно привести к равенству 1/А* = 1/А0 + (1/а)(1/Я0), результаты этой обработки целесообразно представить в виде линейной зависимости 1/А* от 1 /Я0. В табл. 3 сведены результаты вычислений по формулам, номера которых отвечают нижнему индексу у символа А* (индекс 16 означает использование формулы для полусуммы (А1 + А2)/2 соответствующих двусторонних оценок).
1/R0, 1/мкм 0,0120 0,0231 0,0352 0,0539 0,0855 0,1183 0,2247
1/Л25, мм-К/Вт 4,7652 5,2614 5,4920 5,4700 6,0251 6,2168 10,5000
1/Л26, мм-К/Вт 4,7649 5,2599 5,4893 5,4675 6,0178 6,2070 10,3000
1/Л28, мм-К/Вт 4,7647 5,2584 5,4866 5,4648 6,0103 6,1969 10,1000
1/Л29, мм-К/Вт 4,7644 5,2560 5,4827 5,4610 5,9988 6,1794 9,7341
1/Л15, мм-К/Вт 4,7606 5,2459 5,4703 5,4485 5,9821 6,1624 10,1000
1/Л16, мм-К/Вт 4,7648 5,2587 5,4874 5,4655 6,0119 6,1980 10,1000
Из табл. 3 видно, что ни одна из использованных формул не приводит к предполагаемой линейной зависимости 1/Л* от 1 /R0 при постоянном значении коэффициента 1 /а, который по условиям теплового взаимодействия матрицы и включений не должен зависеть от их радиуса. Из равенства 1/Л* = 1/Л0 + (1/a)(1/R0) следует, что 1/Л* ^ 1/Л0 = const при 1/R0 ^ 0, т.е. графики всех предполагаемых линейных зависимостей, соответствующих использованным формулам, должны при 1/R0 = 0 иметь общую точку с ординатой 1/Л0. Для обработки экспериментальных результатов было принято Л0 = 253 Вт/(м-К) [17], т.е. 1/Л0 ~ 3, 9526 мм-К/Вт. По данным, представленным в табл. 3, при помощи линейной регрессии с использованием метода наименьших квадратов получены значения коэффициента 1 /а и модуля M вектора невязки, приведенные в табл. 4.
Таблица 4
Формула (25) (26) (28) (29) (15) (15)
1/а, мм2-К/Вт M, мм К/Вт 0,02909 1,56438 0,02867 1,53766 0,02824 1,52110 0,02751 1,51852 0,02805 1,51970 0,02827 1,52208
Следует отметить, что использование формул для верхней А1 и нижней А2 оценок, определяемых из соотношения (14), приводит к существенно отличающимся от представленных в табл. 4 значениям 1/а и М. Поскольку отношение А0/Ат ~ 1, 0675 достаточно близко к единице, обработка экспериментальных данных по всем использованным формулам привела практически к одинаковым результатам. Формально предпочтение следует отдать формуле (29), которой соответствует наименьшее значение модуля М вектора невязки при применении линейной регрессии для обработки экспериментальных результатов. Наибольшее значение М получено при использовании формулы (25).
Представляет интерес обработка экспериментальных данных для композитов, у которых отношение Л = А0/Ат в большей мере отличается от единицы как в меньшую, так и в большую сторону. Такие данные приведены в работе [18] для композитов с шаровыми включениями из нитрида алюминия (АШ), имеющего коэффициент теплопроводности А0 = = 155 Вт/(м-К), модифицирующими металлические матрицы: медную (ЛСи ~ 0, 3865), алюминиевую (ЛА1 ~ 0,6540), оловянную (Лэп ~ 2, 3134) и свинцовую (Лрь ~ 4, 4286). В табл. 5 приведены результаты измерения эффективного коэффициента теплопроводности А таких композитов, модифицированных включениями, принадлежащими пяти фракциям с различными средними значениями диаметров 2Я0.
2Яо Су А, Вт/(мК)
мкм АШ/Си АШ/А1 АШ/Бп АШ/РЬ
92,1 0,59 223 ± 7 178 ± 5 98 ± 3 74 ± 2
63,7 0,62 210 ± 6 173 ± 5 97 ± 3 75 ± 2
41,0 0,62 204 ± 6 170 ± 5 92 ± 3 71 ± 2
26,6 0,57 206 ± 6 168 ± 5 82 ± 2 60 ± 2
19,5 0,60 187 ± 6 159 ± 5 73 ± 2 53 ± 2
Обработка этих данных по описанной выше методике привела к значениям модуля М вектора невязки, представленным в табл. 6 и имеющим размерность мм-К/Вт. Для каждого типа матрицы композита при помощи каждой из формул вычислены три значения М, причем среднее из них вычислено по номинальному значению А из табл. 5, а верхнее и нижнее — по значениям А соответственно с верхним и нижним допуском.
Таблица 6
Формула (25) (26) (28) (29) (15) (16)
АШ/Си 0,3821 0,3904 0,7458 0,3405 0,2674 0,6365 0,3642 0,1580 0,5014 0,4801 0,1777 0,3002 0,6068 0,3121 0,2911 0,3454 0,1989 0,5547
АШ/А1 0,2721 0,1515 0,4989 0,2730 0,1423 0,4872 0,2742 0,1327 0,4748 0,2746 0,1179 0,4574 0,3214 0,0872 0,3837 0,2714 0,1348 0,4802
АШ/Бп 0,8077 0,6728 0,7766 0,7561 0,6694 0,8090 0,7168 0,6733 0,8410 0,6841 0,6841 0,8772 0,9511 0,7440 0,7387 0,7349 0,6704 0,8256
АШ/РЬ 1,3245 1,2414 1,2173 0,9446 1,0034 1,1151 0,9405 1,0862 1,2568 1,1019 1,2760 1,4622 1,2575 1,2386 1,2627 0,9057 1,0097 1,1554
Из табл. 6 видно, что не удается выбрать какую-либо одну универсальную формулу, которая по критерию минимума модуля М была бы наилучшей для всех рассмотренных вариантов композита. При Л < 1 (композиты АШ/Си и АШ/А1) наиболее приемлемыми можно считать формулы (28) и (29). Значение М = 0, 0872 мм-К/Вт, вычисленное для композита АШ/А1 с использованием соотношения (14), определяющего полусумму (А1 + А2)/2, является наименьшим из всех представленных в табл. 6, но большая чувствительность результатов при сравнительно малом отклонении от номинального значения А не позволяет считать такой подход надежным. Это замечание, хотя и в меньшей степени, относится ко всем формулам, использованным применительно к указанным композитам. Наоборот, при Л > 1 (композиты АШ/Бп и АШ/РЬ) отклонения от номинального значения А сравнительно слабо влияют на значения модуля М, вычисленного с применением всех рассмотренных формул. Для этих композитов предпочтенее можно отдать формуле (26) и соотношениям (16), определяющим полусумму (А* + А*)/2.
Проведенное сопоставление с экспериментом следует считать предварительным, поскольку оно опирается на весьма ограниченный объем экспериментальных данных в срав-
нительно узком интервале изменения параметра Л и объемной концентрации Су включений. Накопление результатов тщательно проведенных экспериментов с композитами в более широком интервале изменения Л и Су включений позволит сделать более определенные выводы о достоверности той или иной расчетной формулы для оценки эффективного коэффициента А теплопроводности композитов с шаровыми включениями.
Заключение
Проведен сравнительный количественный анализ различных математических моделей теплового взаимодействия шаровых включений и матрицы композита, используемых для построения оценок эффективного коэффициента теплопроводности такого композита. С применением двойственной формулировки вариационной задачи стационарной теплопроводности в неоднородном твердом теле получены верхняя и нижняя границы возможных значений этого коэффициента. Для сближения этих границ использован вариационный принцип Хашина — Штрикмана. Установлена связь между оценками, полученными с применением этого принципа и обобщенного сингулярного приближения теории случайных функций.
Выявлено слабое влияние на эффективный коэффициент теплопроводности композита изменения формы включений на кубическую. При помощи модели теплового взаимодействия однородного материала как с шаровыми включениями, так и с частицами матрицы, выведена расчетная формула, обеспечивающая оценку искомого коэффициента теплопроводности этого материала во всем диапазоне возможного изменения объемной концентрации включений.
С привлечением ограниченного объема экспериментальных результатов проведено сравнение рассмотренных расчетных формул по критерию минимума модуля вектора невязки, возникающей при обработке этих результатов с применением линейной регрессии и метода наименьших квадратов.
Работа выполнена по гранту НШ-255.2012.8 программы Президента РФ поддержки ведущих научных школ.
Список литературы
1. Оделевский В.И. Расчет обобщенной проводимости гетерогенных систем // ЖТФ. 1951. Т. 21, вып. 6. С. 667-685.
2. Чудновский А.Ф. Теплофизические характеристики дисперсных материалов. М.: Физ-матгиз, 1962. 456 с.
3. Карслоу Г., Егер Д. Теплопроводность твердых тел: пер. с англ. М.: Наука, 1964. 488 с.
4. Миснар А. Теплопроводность твердых тел, жидкостей, газов и их композиций: пер. с франц. М.: Мир, 1968. 464 с.
5. Дульнев Г.Н., Заричняк Ю.П. Теплопроводность смесей и композиционных материалов. Л.: Энергия, 1974. 264 с.
6. Шермергор Т.Д. Теория упругости микронеоднородных сред. М.: Наука, 1977. 400 с.
7. Кристенсен Р. Введение в механику композитов: пер. с англ. М.: Мир, 1982. 336 с.
8. Зарубин В.С., Кувыркин Г.Н., Савельева И.Ю. Эффективный коэффициент теплопроводности композита с шаровыми включениями // Тепловые процессы в технике. 2012. № 10. С. 470-474.
9. Зарубин В.С. Математическое моделирование в технике. М.: Изд-во МГТУ им. Н.Э. Баумана, 2010. 496 с.
10. Maxwell C. Treatise on electricity and magnetism. Oxford, 1873.
11. Зарубин В.С. Инженерные методы решения задач теплопроводности. М.: Энергоатом-издат, 1983. 328 с.
12. Зарубин В.С., Кувыркин Г.Н. Математические модели механики и электродинамики сплошной среды. М.: Изд-во МГТУ им. Н.Э. Баумана, 2008. 512 с.
13. ЗарубинВ.С., КувыркинГ.Н., Савельева И.Ю. Эффективный коэффициент теплопроводности композита при неидеальном контакте шаровых включений и матрицы // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. 2012. Спец. выпуск «Прикладная математика и механика». С. 84-94.
14. Савельева И.Ю. Эффективный коэффициент теплопроводности композита при наличии промежуточного слоя между шаровыми включениями и матрицей // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. 2012. Спец. выпуск «Математическое моделирование в технике». С.180-186.
15. Зарубин В.С., Кувыркин Г.Н., Савельева И.Ю. Эффективный коэффициент теплопроводности композита при непрерывном изменении теплопроводности между шаровыми включениями и матрицей // Вестник МГТУ им. Н.Э. Баумана. Сер. Естественные науки. 2012. Спец. выпуск «Прикладная математика и механика». С. 95-102.
16. Hashin Z., Shtrikman S. A variational approach to the theory of the effective magnetic permeability of multiphase materials // J. Appl. Phys. 1962. Vol. 33, no. 10. P. 3125-3131. D0I:10.1063/1.1728579.
17. Molina J.M., Narciso J., Weber L., Mortensen A., Louisl E. Thermal conductivity of Al-SiC composites with monomodal and bimodal particle size distribution // Materials Science and Engineering A. 2008. Vol. 480, no. 1-2. P. 483-488.
18. Kida M., Weber L., Monachon C., Mortensen A. Thermal conductivity and interfacial conductance of Al/N particle reinforced metal matrix composites // J. Appl. Phys. 2011. Vol. 109. Art. no. 064907. DOI: 10.1063/1.3553870.
19. Bruggemann D.A.G. Berechnung verschiedener physikalischer Konstanten von heterogenen Substanzen. I. Dielektrizitatskonstanten und Leitfahigkeiten der Mischkorper aus isotropen Substanzen // Ann. der Phys. 1935. Vol. 416, iss. 7. P. 636-664. DOI: 10.1002/andp. 19354160705.
SCIENTIFIC PERIODICAL OF THE BAUMAN MSTU
SCIENCE and EDUCATION
EL № FS77 - 48211. №0421200025. ISSN 1994-0408
electronic scientific and technical journal
Comparative analysis of estimations of heat conduction
of a composite with ball inclusions
# 07, July 2013
DOI: 10.7463/0713.0569319
Zarubin V. S., Kuvyrkin G. N., Savelyeva I. Yu.
Bauman Moscow State Technical University 105005, Moscow, Russian Federation
fn2@bmstu.ru
The authors consider different mathematical models of thermal interaction of ball inclusions and the matrix of a composite, used for estimating the effective heat conduction coefficient of this kind of composite. Application of a dual formulation of the variation problem of stationary thermal conductivity in a nonhomogeneous solid body allowed to get double-sided estimations of possible values of the coefficient. It was determined that transition from ball inclusions to cubic ones influenced slightly the effective heat conduction coefficient of a composite. The authors derived an evaluation formula that allowed to obtain a reliable estimation of this coefficient within the whole range of possible variations of the volume concentration of inclusions. Due to the electrothermal analogy it was possible to apply the obtained results in order to estimate characteristics of electrical conductivity and of dielectric capacitance of composites modified with ball inclusions (considering nanostructured elements).
References
1. Odelevskii V.I. Raschet obobshchennoi provodimosti geterogennykh sistem [Calculation of generalized conductivity of heterogeneous systems]. Zhurnal tekhnicheskoi fiziki [Technical physics], 1951, vol. 21, no. 6, pp. 667-685.
2. Chudnovskii A.F. Teplofizicheskie kharakteristiki dispersnykh materialov [Thermophysical characteristics of disperse materials]. Moscow, Fizmatgiz, 1962. 456 p.
3. CarslawH.S., Jaeger J.C. Conduction of Heat in Solids. 2nd ed. Oxford University Press, 1959. (Russ. ed.: Karslou G., Eger D. Teploprovodnost' tverdykhtel. Moscow, Nauka, 1964. 488 p.).
4. Misnar A. Teploprovodnost' tverdykh tel, zhidkostei, gazov i ikh kompozitsii [The thermal conductivity of solids, liquids, gases and their compositions]. Transl. from French. Moscow, Mir, 1968. 464 p.
5. Dul'nevG.N.,ZarichniakIu.P. Teploprovodnost' smesei i kompozitsionnykh materialov [Thermal conductivity of mixtures and composite materials]. Leningrad, Energiia, 1974. 264 p.
6. Shermergor T.D. Teoriia uprugosti mikroneodnorodnykh sred [Theory of elasticity of microin-homogeneous media]. Moscow, Nauka, 1977. 400 p.
7. Christensen R.M. Mechanics of composite materials. New York, Wiley-Interscience publication, 1979. 348 p. (Russ. ed.: Kristensen R.M. Vvedenie v mekhaniku kompozitov. Moscow, Mir, 1982. 336 p.)
8. Zarubin V.S., Kuvyrkin G.N., Savel'eva I.Iu. Effektivnyi koeffitsient teploprovodnosti kom-pozita s sharovymi vkliucheniiami [The effective thermal conductivity of composites with spherical inclusions]. Teplovye protsessy v tekhnike, 2012, no. 10, pp. 470-474.
9. Zarubin V.S. Matematicheskoe modelirovanie v tekhnike [Mathematical modeling in engineering]. Moscow, Bauman MSTU Publ., 2010. 496 p.
10. Maxwell C. Treatise on electricity and magnetism. Oxford, 1873.
11. Zarubin V.S. Inzhenernye metody resheniia zadach teploprovodnosti [Engineering methods for solving problems of thermal conductivity]. Moscow, Energoatomizdat, 1983. 328 p.
12. Zarubin V.S., Kuvyrkin G.N. Matematicheskie modeli mekhaniki i elektrodinamiki sploshnoi sredy [Mathematical models of mechanics and electrodynamics of continuous media]. Moscow, Bauman MSTU Publ., 2008. 512 p.
13. Zarubin V.S., Kuvyrkin G.N., Savel'eva I.Iu. Effektivnyi koeffitsient teploprovodnosti kom-pozita pri neideal'nom kontakte sharovykh vkliuchenii i matritsy [The effective coefficient of thermal conductivity of the composite at nonideal contact of spherical inclusions and matrix]. Vestnik MGTU im. N.E. Baumana. Ser. Estestvennye nauki [Bulletin of the Bauman MSTU. Ser. Natural science], 2012, Spets. vypusk <<Prikladnaia matematika i mekhanika>> [Spec. iss. "Applied mathematics and mechanics"], pp. 84-94.
14. Savel'eva I.Iu. Effektivnyi koeffitsient teploprovodnosti kompozita pri nalichii promezhutochnogo sloia mezhdu sharovymi vkliucheniiami i matritsei [The effective coefficient of thermal conductivity of the composite in the presence of an intermediate layer between the spherical inclusions and matrix]. Vestnik MGTU im. N.E. Baumana. Ser. Estestvennye nauki [Bulletin of the Bauman MSTU. Ser. Natural science], 2012, Spets. vypusk «Matematicheskoe modelirovanie v tekhnike» [Spec. iss. "Mathematical modeling in engineering"], pp. 180-186.
15. Zarubin V.S., Kuvyrkin G.N., Savel'eva I.Iu. Effektivnyi koeffitsient teploprovodnosti kom-pozita pri nepreryvnom izmenenii teploprovodnosti mezhdu sharovymi vkliucheniiami i ma-tritsei [The effective coefficient of thermal conductivity of the composite at continuous change of thermal conductivity between the spherical inclusions and matrix]. Vestnik MGTU im. N.E. Baumana. Ser. Estestvennye nauki [Bulletin of the Bauman MSTU. Ser. Natural science],
2012, Spets. vypusk «Prikladnaia matematika i mekhanika>> [Spec. iss. "Applied mathematics and mechanics"], pp. 95-102.
16. Hashin Z., Shtrikman S. A variational approach to the theory of the effective magnetic permeability of multiphase materials. J. Appl. Phys., 1962, vol. 33, no. 10, pp. 3125-3131. DOI: 10.1063/1.1728579.
17. Molina J.M., Narciso J., Weber L., Mortensen A., Louisl E. Thermal conductivity of Al-SiC composites with monomodal and bimodal particle size distribution. Materials Science and Engineering A, 2008, vol. 480, no. 1-2, pp. 483-488.
18. Kida M., Weber L., Monachon C., Mortensen A. Thermal conductivity and interfacial conductance of Al/N particle reinforced metal matrix composites. J. Appl. Phys., 2011, vol. 109, art. no. 064907. DOI: 10.1063/1.3553870.
19. Bruggemann D.A.G. Berechnung verschiedener physikalischer Konstanten von heterogenen Substanzen. I. Dielektrizitatskonstanten und Leitfähigkeiten der Mischkörper aus isotropen Substanzen. Ann. der Phys., 1935, vol. 416, iss. 7. pp. 636-664. DOI: 10.1002/andp. 19354160705.
