CC BY
13ISSN 2782-2001 Системы анализа и обработки данных том 86, № 2, 2022, с. 105-120
http://journals.nstu.ru/vestnik Analysis and data processing systems Vol. 86, No. 2, 2022, pp. 105-120
КОМПЬЮТЕРНЫЕ НАУКИ И ИНФОРМАТИКА
COMPUTER SCIENCE
AND INFORMATION SCIENCE
УДК 621.3.019.34
DOI: 10.17212/2782-2001-2022-2-105-120
Расчет температурного поля многослойных несущих конструкций
численными методами
Д. В. ОЗЕРКИН
634050, г. Томск, пр. Ленина, 40, Томский государственный университет систем управления и радиоэлектроники
ozerkin.denis@yandex.ru
Современные подходы к проектированию электронных средств включают в себя моделирование теплофизических процессов, происходящих в разрабатываемой аппаратуре. При этом распространенной задачей в теплофизическом моделировании выступает задача определения температурного поля несущей конструкции электронного средства (платы, подложки). Знание количественных показателей температурного поля несущей конструкции напрямую влияет на прогнозирование надежности всего электронного средства. Несущие конструкции электронных средств всё чаще выполняются в виде многослойных конструкций, где могут присутствовать проводящие, изолирующие и полупроводниковые слои. Указанные особенности предъявляют новые требования к разработке математических моделей, алгоритмов, пакетов прикладных программ для вычисления температурных полей несущих конструкций электронных средств.
В статье рассмотрен процесс моделирования температурного поля для несущих конструкций электронных средств с помощью распространенных программных комплексов Ma1:hCAD и БоШШогкз. Всё многообразие несущих конструкций электронных средств поделено на три категории: условно одномерные, условно двумерные и трехмерные. Показано, что температурные поля условно одномерных (стержни) и условно двумерных (однослойные платы и подложки) конструкций целесообразно вычислять методом конечных разностей в программном комплексе Ma1:hCAD. Температурные поля трехмерных несущих конструкций, в том числе, сложной геометрической конфигурации, целесообразно вычислять методом конечных элементов в программном комплексе БоШШогкз. Разработанные алгоритмы расчета температурных полей проверены посредством решения тестовой задачи. Сравнительный анализ показал, что расхождение результатов расчетов относительно тестовой задачи не превышает 0,8 %. Разработанные алгоритмы расчета температурного поля могут быть практически полезны в инженерной деятельности разработчика электронной аппаратуры.
Ключевые слова: уравнение теплопроводности, математическая модель, метод конечных разностей, метод конечных элементов, коэффициент температуропроводности, граничные условия, пространственная сетка, временной слой, теплопроводность
*
Статья получена 29 января 2022 г.
ВВЕДЕНИЕ
Аналитические методы расчета температурного поля позволяют реализовать лишь простейшие одномерные модели. По существу, удается получить аналитические выражения обозримого вида лишь для решений отдельных линейных уравнений [1]:
1) уравнение для одномерного температурного поля стержня;
2) уравнение для плоской стенки, сплошных или полых цилиндров и шаров с температурным полем, изменяющихся только по толщине стенок;
3) уравнение для одномерного потока теплоносителя.
При переходе к системам таких уравнений круг задач, для которых возможно эффективно использовать аналитические методы, еще более сужается.
Применение вычислительной техники и численных методов расчета значительно расширяет классы исследуемых задач теплообмена, позволяя получать приближенные решения многомерных, нелинейных, нестационарных задач, для которых использование точных и приближенных аналитических методов не представляется возможным.
За основу рассматриваемой ниже математической модели распределения температурного поля приняты известные численные методы [2-6] решения дифференциальных уравнений. Однако, в отличие от универсальных методик расчета температурного поля, предлагается разновидность математической модели для решения узкоспециализированной задачи: расчет температурного поля плоской многослойной конструкции - платы или подложки (далее - твердотельная система) с соответствующими допущениями и ограничениями.
Проблема надежности современной радиоэлектронной аппаратуры неразрывно связана с расчетом температурного поля несущих конструкций электронных средств. Общая задача расчета температурного поля несущих конструкций численными методами возможна в различных постановках, а именно: в одномерной, двумерной и трехмерной постановке с точки зрения пространственных координат. Исходя из этого для предстоящего исследования анизотропности теплового поля несущих конструкций электронных средств целесообразно ввести их классификацию:
- условно одномерные (стержни);
- условно двумерные (однослойные платы и подложки);
- трехмерные (многослойные платы и подложки).
1. ПОСТАНОВКА ЗАДАЧИ
Эта задача сводится к анализу обобщенной теплофизической модели твердотельной системы, созданной по объемной многослойной технологии [7, 8]. Численное моделирование температурных полей объектов иной геометрической конфигурации представлено во многих отечественных и зарубежных работах, например в [9, 10]. Графическое представление указанной задачи представлено на рис. 1. Структура параллелепипеда из г слоев имеет размер 1Х и 1у по осям абсцисс и ординат соответственно. По оси аппликат отмечены слои твердотельной системы размерностью И^, Иг,.., И, Иг . Общая толщина твердотельной системы будет Н = ^ И . В каждом слое И твердотельной системы может находиться источник тепловыделения с плотностью тепловы-
деления qf (х, у), причем координаты расположения источника по оси абсцисс и ординат могут быть произвольными, а по оси аппликат - фиксированы.
На шести поверхностях твердотельной системы (верхней, нижней и четырех боковых) заданы условия конвективно-радиационного теплообмена с окружающей средой, которые в общем случае различны: ав, ан, «бъ аБ2, аБ3, аБ4 . Уравнение для теплопроводности для /-го слоя твердотельной системы выглядит следующим образом:
дТ
ср-
дг
г = Л
( д2Tj (x, y, z) д2Tj (x, y, z) d1Ti (x, y, z)Л
дх
ду 2
дz 2
F
-Z q(x, y), (i) f=i
где с - коэффициент теплоемкости материала /-го слоя твердотельной системы; р - плотность материала; X - коэффициент теплопроводности материала.
ly
0
Рис. 1. Графическое представление теплофизической модели твердотельной системы
Fig. 1. Graphical representation of the thermophysical model of a solid state system
Граничные условия III и IV рода для слоев 1, r, i выглядят так:
дТ1
-Л1 дzl
Л дТг
л r-
= аН (Тср - Т1) ;
z=0
дzr
= а В (Тср - Tr ) ;
z=Н
T\z=hj Ti +1 lz =hj+1 '
-л, дТ
дzi
= -Л
дТ
i+1
i+1
z=h
дz-¡
i+1
z=h+1
дТ[
= аБ1(Тср - Ti);
дx,
(2)
(3)
(4)
(5)
(6)
x
дТ-
Ь дХХ~ = аБ2(Тср -Т); (7)
= аБз(Тср -Т); (8)
дуг дТ
Ь ^ = а Б4 (Тср - Т), (9)
дУг
где Тср - температура окружающей среды. Начальное условие имеет вид
Т| т=0 = Тср (х, у, z). (10)
Искомой в задаче (1)—(10) является функция Т (х, у, z, т), заданная в непрерывной области:
0 = {0 < х <1Х }х {0 < у <1у }х{0 < z < Н }х {0 < г <^тах}.
Известны два основных численных метода решения уравнения (1) в частных производных [1, 11]:
— метод конечных разностей (МКР);
— метод конечных элементов (МКЭ).
2. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ ТЕПЛОПРОВОДНОСТИ, ПРИГОДНЫЕ ДЛЯ ПРАКТИЧЕСКОГО РЕШЕНИЯ ЗАДАЧ РАСПРЕДЕЛЕНИЯ ТЕМПЕРАТУРНОГО ПОЛЯ
2.1. МЕТОД КОНЕЧНЫХ РАЗНОСТЕЙ ДЛЯ РЕШЕНИЯ ЗАДАЧ РАСПРЕДЕЛНИЯ ТЕМПЕРАТУРНОГО ПОЛЯ
Наиболее приемлема для расчетов тепловых режимов плоской конструкции локально-одномерная схема [12, 13], отличительной особенностью которой является сочетание сильных сторон явных схем (малые затраты машинного времени на один шаг по времени) и неявных схем (безусловная устойчивость).
Процедуру численного расчета уравнения теплопроводности (1) по локально-одномерной схеме разделим на следующие этапы.
1. Замена исходной области непрерывного изменения переменных пространственно-временной сеткой.
2. Построение разностной схемы.
3. Решение системы разностных уравнений.
4. Многократный повтор решения системы разностных уравнений вдоль координат х, у, z.
При расщеплении уравнения теплопроводности (1) точное решение в конце малого промежутка времени [Tj-1, Tj ] трансформируем в приближенное решение ю(х, у, т ■).
Перед построением локально-одномерной схемы для /-го слоя твердотельной системы примем, что пространственная Xх, Xу, Xг и временная Ах
сетки - эквидистантные, после чего проведем дискретизацию задачи (1)-(9) по пространственным переменным с использованием неявных схем [1, 12, 13]. В результате получим систему конечно-разностных уравнений. Применительно к пространственной сетке Xх дискретизация будет выглядеть как решение и(х, у, х/ ) на данном промежутке времени:
л j _ лj _
п,т,к п,т,к _ а (а / _ 2Л j + А j )-
2 I п+1,т,к п,т,к п_1,т,к I х х
Ах
Я/ (х, у) .
ср
(11)
а0 _ т / _ 1.
**п,т,к 10п,т,к? J
Н,т,к =&]п-п,к, j > 1;
й-/ _ О,./ л1, т, к л0, т, к
х х
аБ1 (Тср _Л0,т,к ) ;
(12)
(13)
^Ы,т,к N _1,т,к х х
_ аБ2 (Тср _ ,т,к
) •
(14)
где п, т, к - индексы пространственной сетки дискретизации по координатам х, у, г; / - индекс временной сетки дискретизации; N - крайний номер индекса пространственной сетки дискретизации по координате х. / - индекс номера слоя твердотельной системы; а _ —— коэффициент температуропровод-
ср
ности материала /-го слоя твердотельной системы.
Дискретизация по пространственным сеткам Xу, Xг выглядит аналогично. Применительно к пространственной сетке хг появляются граничные условия 4-го рода:
(
ю
}
п,т,К юп,т,К _1 X г
Л
(
/+1
Л
ю
}
п,т,1 Юп,т,0
X г
Л
(15)
/+1
ю / 1 _ ю / П п,т,1 п,т,0
_аН Мер
(Тср юп,т,0 );
г _0
ю/ у — Ю/ V 1
_аВ (Тср _<т,К ),
г _Н
X
г
г
где г - крайний индекс слоя твердотельной системы; ю - решение на данном промежутке времени по пространственной координате г; К - крайний номер индекса пространственной сетки дискретизации по координате г.
Искомое решение имеет вид
<т,к = Т(х, У, г, ^]) + 0(Дт2), (17)
где 0 (Дт ) - погрешность численного расчета (погрешность конечно-разностной математической модели).
На рис. 2 представлен алгоритм расчета температурного поля плоской конструкции по локально-одномерной схеме. Решение системы разностных уравнений (11)—(14) сводится к расчету на каждом шаге по времени набора одномерных задач методом правой прогонки.
2.2. МЕТОД КОНЕЧНЫХ ЭЛЕМЕНТОВ ДЛЯ РЕШЕНИЯ ЗАДАЧ РАСПРЕДЕЛНИЯ ТЕМПЕРАТУРНОГО ПОЛЯ
В МКЭ твердотельная плоская конструкция представляет собой непрерывную область О с границами Г. Проводится разбиение пространства ^-го слоя твердотельной системы на конечные элементы - тетраэдры. Каждый элементарный тетраэдр имеет порядковый номер ц и узлы (вершины): к, I.
Объем конечного элемента обозначим как Vе. Полагаем, что текущая температура Т в каждом конечном элементе ц может быть аппроксимирована как функция координат х, у, г. Пространственное распределение температуры внутри любого конечного элемента аппроксимируется суммой произведений координат на коэффициенты, равные приближенным значениям температуры в узловых точках, принадлежащих данному элементу.
Применим вариационный подход к решению уравнения (1). Известно [14], что задача решения уравнения (1) эквивалентна задаче определения функции Т (х, у, г, т), минимизирующей функционал I [Т (х, у, г, т)]. Чтобы функционал I принимал минимальное значение, его подынтегральная функция должна в соответствии с методом Ритца [14] удовлетворять уравнению
I [Т(х, у, г, I)] = 0,
или
Ш
О
атл
Г л2
д Тл (х, у, г) д Тл (х, у, г) д Тл (х, у, г)
^ Л
■ Е ч/(x, у) /=1
дх2
ду 2
Г
Г дТ
+"а(Тср " ТЛ)
IV_
дг 2
Л
с1Г = 0,
(18)
где п - внешняя нормаль к поверхностям твердотельной многослойной системы.
С
х
3
Установка начальных значенийj := О
I
Установка граничных значений по х
I
Инициализация счетчика числа «продольных» стержней m := 1
К
Решение одномерной задачи в направлении координаты х (метод прогонки стержня т )
I
Инкремент числа «продольных» стержней m := m + 1
да
X
m < M
I
нет
Установка граничных значений по y
X
Инициализация счетчика числа «поперечных» стержней n := 1
Решение одномерной задачи в направлении координаты у (метод прогонки стержня п )
Инкремент числа «поперечных» стержней n := n + 1
да
n < N
1
нет
Установка граничных значений по z
I
Инициализация счетчика числа «вертикальных» стержней к := 1
г
Решение одномерной задачи в направлении координаты г (метод прогонки стержня к )
Инкремент числа «вертикальных» стержней к := к + 1
да
I
к < K
нет
Инкремент временного слоя j :=j + 1
X
да
j<J
< ^max)
С
нет
J
Рис. 2. Обобщенный алгоритм расчета температурного поля по локально-одномерной схеме
Fig. 2. A generalized algorithm for calculating the temperature field using a locally one-dimensional scheme
После интегрирования и суммирования выражение (18) будет представлять систему алгебраических уравнений относительно температуры всех ц узлов. Представим уравнение в векторно-матричной форме:
Т Р +1 _ т Р , С--Л • ТР+1 _ Е . (19)
Для получения явной расчетной схемы выразим из (19) температуру Т Р+1:
С
—ТР + Е т р+1 _ А_
_ С _ л
Аt
на новом временном слое.
3. МОДЕЛИРОВАНИЕ РАСПРЕДЕЛЕНИЯ ТЕМПЕРАТУРНОГО ПОЛЯ ЧИСЛЕННЫМИ МЕТОДАМИ В РАЗЛИЧНЫХ ПОСТАНОВКАХ
Программный комплекс МаШСАБ - это компьютерное программное обеспечение, предназначенное в первую очередь для верификации, валидации, документирования и повторного использования инженерных расчетов.
Используя функциональные возможности МаШСАБ по программированию, становится возможным разработать программный блок, реализующий конечно-разностный метод расчета температурного поля на основе математической модели (11)—(17).
Практическое применение разработанного программного блока продемонстрируем на примере решения нестационарной задачи теплопроводности в одномерной и двумерной постановке (см. задачи № 1 и № 2).
Задача № 1. Требуется найти распределение температуры вдоль стержня длиной 1х = 0,05 м в течение времени Хфин = 100 с. Материал стержня - стеклотекстолит. Теплофизические параметры стеклотекстолита с = 840 Дж/(кг • °С); р = 1800 кг/м ; X = 0,3 Вт/(м • °С). Внутри стержня в пределах х1 < х < х2 находится источник тепла с интенсивностью тепловыделения Q = 1 • 105 Вт/м, где х1 = 0,0197 м; х2 = 0,0297 м. Начальная температура стержня Т0 = 20 °С. Температура окружающей среды на левом и правом торцах стержня Тср.лев = = Тср пр = 120 °С. На торцах стержня заданы граничные условия 3-го рода с коэффициентом конвективной теплоотдачи ак лев = ак пр
= 10 Вт/(м2 • °С). Количество узлов пространственной сетки N = 100. Шаг временной сетки Ах = 1 с.
Задача № 2. Требуется найти распределение температуры по двумерной пластине длиной и шириной 1х _ 1у = 0,05 м в течение времени Хфин = 100 с.
Материал пластины - стеклотекстолит. Теплофизические параметры стеклотекстолита с = 840 Дж/(кг • °С); р = 1800 кг/м ; X = 0,3 Вт/(м • °С). Внутри пластины в пределах х1 < х < х2 и у1 < у < у2 находится источник тепла с интен-
5 2
сивностью тепловыделения Q = 1-10 Вт/м , где х^ = 0,02 м; Х2 = 0,03 м; У1 = 0,02 м; У2 = 0,03 м. Начальная температура пластины Т = 20 °С. Температура окружающей среды на торцах пластины Тср = 20 °С. На торцах пластины заданы граничные условия 3-го рода с коэффициентом конвективной теплоотдачи ак = 10 Вт/(м - °С). Количество узлов пространственной сетки N = М = 100. Шаг временной сетки Ах = 1 с.
Результаты численного решения задачи № 1 представлены в виде температурного профиля вдоль координат х в момент времени х = 100 с (рис. 3, а). Результаты численного решения задачи № 2 представлены в виде поверхности - температурного поля - в момент времени х = 100 с (рис. 3, б).
а б
Рис. 3. Результаты решения задач в одномерной и двумерной постановке Fig. 3. Results of solving problems in one-dimensional and two-dimensional formulations
Программный комплекс SolidWorks - это компьютерная программа для твердотельного моделирования и автоматизированного проектирования, которая работает под управлением операционной системы Microsoft Windows.
Практическое применение модуля SolidWorks Simulation Thermal Analysis продемонстрируем на примере решения задачи (см. задачу № 3) теплопроводности в трехмерной постановке для стационарного и нестационарного режима, а также для случая температурной зависимости коэффициента теплоотдачи а и коэффициента теплопроводности X.
Задача № 3. Требуется найти распределение температуры в течение времени Тфин = 100 с для нестационарного режима с температурной зависимостью
коэффициентов теплоотдачи и теплопроводности в подложке с источником тепловыделения (кристаллом). Габаритные размеры подложки: lx = ly = 0,01 м;
lz = 0,001 м. Материал подложки - кремний. Теплофизические параметры кремния: с = 648 Дж/(кг • °С); р = 2330 кг/м ; X = 150 Вт/(м • °С). В пределах Х1 < x < Х2 , У1 < y < У2, Z1 < z < Z2 находится кристалл с мощностью тепловыделения P = 0,1 Вт, где x1 = 0,0045 м; x2 = 0,0055 м; y1 = 0,0045 м; У2 = 0,0055 м; Z1 = 0 м; Z2 = 0,001 м. Тепловая мощность кристалла увеличивается от нуля во время т = 0 и достигнет максимального значения (0,1 Вт)
через 100 с. Габаритные размеры кристалла: 0,001 х 0,001 х 0,001 м. Материал кристалла - кремний. Начальная температура подложки Tq = 20 °С. Температура окружающей среды ТСр = 20 °С. На поверхностях подложки заданы граничные условия 3-го рода с коэффициентом конвективной теплоотдачи ак = = 15 Вт/(м2 • °С).
Используя функциональные возможности дополнительного модуля инженерного расчета SolidWorks Simulation Thermal Analysis, проведем расчет температурного поля подложки. Подложка в терминах SolidWorks представляет собой сборку, состоящую из собственно подложки и одного кристалла. Кристалл размещен симметрично относительно продольной и поперечной оси подложки.
В результате в основном окне построения сборки отобразится температурное поле конструкции в виде термограммы - цветового представления распределения температур (рис. 4, а).
Применив функциональную возможность программного комплекса по зондированию произвольного узла пространственной сетки, можно вывести на экран график в декартовой системе координат - зависимость температуры выбранного пространственного узла от времени (рис. 4, б).
*
а б
Рис. 4. Результаты решения задачи в трехмерной постановке
Fig. 4. Results of solving the problem in a three-dimensional formulation
Для оценки погрешности расчета температурного поля рассмотренными выше численными методами воспользуемся решением тестовой задачи. Тестовую задачу решим тремя способами:
1) аналитическим способом;
2) численным способом (МКР) с помощью программного комплекса MathCAD;
3) численным способом (МКЭ) с помощью программного комплекса SolidWorks.
В качестве тестовой задачи воспользуемся задачей № 26 из [15]. Задача № 26. Начальная температура одинакова. Около одной поверхности задана постоянная температура среды, на другой поверхности тепловой поток равен нулю (рис. 5, а). Коэффициент температуропроводности твердого
тела a = 5,72 • 10_6 м2/с; коэффициент теплопроводности X = 20 Вт/(м • К);
длина тела I = 0,05 м; начальная температура тела Г0 = 200 К; температура окружающей среды ТСр = 400 К.
Требуется определить температуру поверхностей.
Аналитический способ
Аналитический способ заключается в определении температур поверхностей с использованием расчетных графиков (номограмм). Определим температуры поверхности при ц = 0. Здесь ц - относительная безразмерная простран-
х
ственная координата ц = у. На графике, представленном на рис. 5, б [15], произвольно выберем одну из функций безразмерной температуры и соответствующий ей аргумент - критерий Фурье.
Расчетная формула для определения температуры на поверхностях тела [15]:
т = т+е(7ср - 7Ъ), (20)
где е - относительная безразмерная температура.
По формуле (20) определим температуру на поверхности при ц = 0:
Т = 200 + 0,6 (400 - 200) = 320 К.
Определим температуры поверхности при ц = 1. На графике, представленном на рис. 5, б [15], выбираем аналогичные значения функции безразмерной температуры и соответствующего ей аргумента: Б1 = 2 и Го = 0,3. По формуле (20) определим температуру на поверхности при ц = 1:
Т = 200 + 0,18 (400 - 200) = 236 К.
Тср T
Рис. 5. Тепловая схема тестовой задачи Fig. 5. A thermal diagram of the test task
Численный способ (МКР) с помощью программного комплекса MathCAD
Для решения тестовой задачи в MathCAD применим программный блок по методу правой прогонки, рассмотренный выше. Результаты численного решения тестовой задачи представлены в виде температурного профиля вдоль координат х в момент времени т = 131 с (рис. 6, а). Определим в MathCAD
значения температур на поверхностях твердого тела (температуры узлов № 1 и № 100):
Т1 = 320,7 К;
Т100 = 234,5 83 К.
Численный способ (МКЭ) с помощью программного комплекса
SolidWorks
Применим для решения тестовой задачи в программном комплексе SolidWorks рассмотренную ранее методику решения задач для нестационарного теплового режима. В результате решения задачи в программном комплексе в основном окне построения сборки отобразится температурное поле конструкции в виде термограммы в момент времени 131 с. С помощью инструмента Зондирование «щелкнем» в десяти точках вдоль длины твердого тела, причем первая и десятая точки соответствуют правому и левому торцу твердого тела соответственно. В результате на экране появится график температурного профиля - зависимость температуры от пространственной координаты x (рис. 6, б). Определим в программном комплексе SolidWorks Simulation точные значения температур на поверхностях твердого тела (температуры узлов № 3240 и № 985):
Т3240 = 319,4 К;
Т985 = 234,2 К.
Для удобства результаты всех рассмотренных выше способов решения тестовой задачи сведены в таблицу. Можно видеть, что в среде MathCAD понижение шага дискретизации пространственной сетки в 10 раз приводит к уменьшению относительной погрешности нахождения температуры на 0,04.. .0,1 %. Вариация размеров конечного элемента в среде SolidWorks дает стабильные результаты нахождения температуры. Так, уменьшение размеров конечного элемента в 4 раза не приводит к сколько-нибудь заметным изменениям в результатах нахождения температуры поверхностей. Численные способы решения тестовой задачи в средах MathCAD и SolidWorks дают относительную погрешность нахождения температуры в диапазоне от 0,5 до 0,8 %.
а б
Рис. 6. Результаты решения тестовой задачи Fig. 6. Results of solving the test problem
Сравнительный анализ решения тестовой задачи Comparative analysis of the test problem solution
№ п/п Способ решения Шаг пространственной сетки дискретизации, мм Температура левой поверхности, К Относительная погрешность, % Температура правой поверхности, К Относительная погрешность, %
1 Аналитический - 320 - 236 -
2 Численный 0,5 320,7 0,2 234,6 0,6
в среде MathCAD 0,1 320,5 0,16 234,7 0,55
0,05 320,5 0,16 234,8 0,5
3 Численный 0,74 319,4 0,2 234,2 0,8
в среде SolidWorks 0,41 319,4 0,2 234,2 0,8
0,18 319,4 0,2 234,2 0,8
ЗАКЛЮЧЕНИЕ
В результате проведенного исследования разработана методика, практическая пригодная для нахождения значений абсолютных и относительных рабочих температур несущих конструкций электронных средств. Аналитический способ решения уравнения теплопроводности (1) для решения современных задач теплопроводности в двух- и трехмерной постановке в нестационарном режиме является чрезвычайно сложной проблемой, которую практически невозможно представить в обозримом виде.
Эффективным способом решения уравнения теплопроводности (1) в различных постановках и в различных режимах являются численные методы, реализуемые с помощью современных программных комплексов, таких как МаШСАБ и 8оШ1№огк8.
Современные программные средства решения уравнения теплопроводности численными методами обладают высокой степенью адекватности. Сравнительный анализ решения тестовой задачи теплопроводности аналитическим и численным методами позволяет сказать, что погрешность численного счета не превышает 0,8 %.
Разработанные программные блоки для решения уравнения теплопроводности в одно-, двух- и трехмерной постановке в среде МаШСАБ целесообразно применять для случаев конструкций электронных средств в форме параллелепипедов. Разработанный программный блок для решения уравнения теплопроводности в трехмерной постановке в среде SolidWorks применим для случаев конструкций электронных средств сложной геометрической конфигурации.
СПИСОК ЛИТЕРАТУРЫ
1. Селютин В.А. Машинное конструирование электронных устройств. - М.: Советское радио, 1977. - 384 с.
2. ДульневГ.Н., СемяшкинЭ.М. Теплообмен в радиоэлектронных аппаратах. - Л.: Энергия, 1968. - 360 с.
3. Дульнев Г.Н., Парфенов В.Г., Сигалов А.В. Применение ЭВМ для решения задач теплообмена: учебное пособие для теплофизических и теплоэнергетических специальностей вузов. -М.: Высшая школа, 1990. - 207 с.
4. Крейт Ф., Блэк У. Основы теплопередачи: пер с англ. - М.: Мир, 1983. - 512 с.
5. Берковский Б.М., Ноготов Е. Ф. Разностные методы исследования задач теплообмена. -Минск: Наука и техника, 1976. - 144 с.
6. Самарский А.А. Теория разностных схем. - М.: Наука, 1989. - 616 с.
7. Овчинников С.В., Ляшенко А.В. Обобщенная модель для аналитического расчета температурных полей и тепловых сопротивлений элементов твердотельной радиоэлектроники // Ге-теромагнитная микроэлектроника. - 2012. - № 12. - С. 11-18.
8. Кандалов П.И., Мадера А.Г. Моделирование температурных полей в многослойных структурах // Программные продукты и системы. - 2008. - № 4. - С. 46-49.
9. Иванцивский В.В. Численное моделирование температурных полей в материалах при упрочнении с использованием концентрированных объемных источников нагрева // Научный вестник НГТУ. - 2004. - № 2 (17). - С. 161-172.
10. Присекин В.Л., Расторгуев Г.И., Белоусов А.И. Анализ температурных полей и напряжений в элементах конструкций при воздействии периодических тепловых потоков // Научный вестник НГТУ. - 2001. - № 2 (11). - С. 133-142.
11 Базилевич Р.П. Декомпозиционные и топологические методы автоматизированного конструирования электронных устройств. - Львов: Вища школа, 1981. - 168 с.
12. Штейн М.Е., Штейн Б.Е. Методы машинного проектирования цифровой аппаратуры. - М.: Советское радио, 1973. - 296 с.
13. Морозов К.К., Одиноков В.Г. Использование ЭЦВМ при конструировании некоторых узлов радиоэлектронной аппаратуры. - М.: Советское радио, 1972. - 104 с.
14. Сегерлинд Л. Применение метода конечных элементов. - М.: Мир, 1979. - 393 с.
15. Пехович А.И., Жидких В.М. Расчеты теплового режима твердых тел. - Л.: Энергия, 1976. - 352 с.
Озеркин Денис Витальевич, декан радиоконструкторского факультета Томского государственного университета систем управления и радиоэлектроники (ТУСУР), доцент кафедры конструирования и производства радиоаппаратуры (КИПР). Основное направление научных исследований - экспериментальные методы и средства исследования технических систем, автоматизированное проектирование радиоэлектронных систем, проектирование термостабильных радиотехнических устройств. Имеет более 50 научных и научно-методических работ. E-mail: ozerkin.denis@yandex.ru
Ozerkin Denis Vitalievich, Dean of the Radio Engineering Department in the Tomsk State University of Control Systems and Radioelectronics (TUSUR), associate professor at the Department of Design and Production of Radio Equipment (KIPR). The main field of his scientific research is experimental methods and means for studying technical systems, computer-aided design of radio electronic systems, and design of thermostable radio engineering devices. He has published more than 50 research and scientific-methodical works. E-mail: ozerkin.denis@yan-dex.ru.
DOI: 10.17212/2782-2001-2022-2-105-120
Computation of the temperature field of multilayer load-bearing structures by numerical methods *
D.V. OZERKIN
Tomsk State University of Control System and Radioelectronics, 40 Lenin Prospekt, Tomsk, 634050, Russian Federation
ozerkin.denis.@yandex.ru
Abstract
Modern approaches to the electronic means design include modeling of thermophysical processes occurring in the apparatus being developed. At the same time, a common task in thermophysical modeling is the problem of determining the temperature field of the supporting structure of an electronic device (board, substrate). Knowledge of the quantitative indicators of the supporting structure temperature field directly affects the prediction of the entire electronic device reliability. The load-bearing structures of electronic devices are increasingly made in the form of multilayer structures where conductive, insulating and semiconductor layers can be present. These features impose new requirements on the development of mathematical models, algorithms, application software packages for calculating the temperature fields of the supporting structures of electronic means.
The article considers the process of modeling the temperature field for the supporting structures of electronic means using common software systems MathCAD and SolidWorks. The whole variety of supporting structures of electronic means is divided into three categories: conditionally one-dimensional, conditionally two-dimensional and three-dimensional. It is shown that the temperature fields of conditionally one-dimensional (rods) and conditionally two-dimensional (single-layer boards and substrates) structures should be calculated by the finite difference method in the MathCAD software package. Temperature fields of three-dimensional load-bearing structures, including complex geometric configurations, should be calculated using the finite element method in the SolidWorks software package. The developed algorithms for calculating temperature fields are verified by solving a test problem. Comparative analysis has shown that the discrepancy between the calculation results relative to the test problem does not exceed 0.8%. The developed algorithms for calculating the temperature field can be practically useful in the engineering activities of the developer of electronic equipment.
Keywords: heat equation, mathematical model, finite difference method, finite element method, thermal diffusivity, boundary conditions, spatial grid, time layer, thermal conductivity
REFERENCES
1. Selyutin V.A. Mashinnoe konstruirovanie elektronnykh ustroistv [Machine design of electronic devices]. Moscow, Sovetskoe radio Publ., 1977. 384 p.
2. Dul'nev G.N., Semyashkin E.M. Teploobmen v radioelektronnykh apparatakh [Heat transfer in radio electronic devices]. Leningrad, Energiya Publ., 1968. 360 p.
3. Dul'nev G.N., Parfenov V.G., Sigalov A.V. Primenenie EVMdlya resheniya zadach tep-loo-bmena [The use of computers for solving problems of heat transfer]. Moscow, Vysshaya shkola Publ., 1990. 207 p.
4. Kreith F., Black W. Basic heat transfer. New York etc., 1980 (Russ. ed.: Kreit F., Blek U. Osnovy teploperedachi. Moscow, Mir Publ., 1983. 512 p.).
5. Berkovskii B.M., Nogotov E.F. Raznostnye metody issledovaniya zadach teploobmena [Difference methods for studying heat transfer problems]. Minsk, Nauka i tekhnika Publ., 1976. 144 p.
*
Received 29 January 2022.
6. Samarskii A.A. Teoriya raznostnykh skhem [Theory of difference schemes]. Moscow, Nauka Publ., 1989. 616 p.
7. Ovchinnikov S.V., Lyashenko A.V. Obobshchennaya model' dlya analiticheskogo rascheta tem-peraturnykh polei i teplovykh soprotivlenii elementov tverdotel'noi radioelektroniki [Synthesized model for analytical calculation of temperature field and thermal resistance elements of solid state electronics]. Geteromagnitnaya mikroelektronika = Heteromagnetic Microelectronics, 2012, no. 12, pp. 11-18.
8. Kandalov P.I., Madera A.G. Modelirovanie temperaturnykh polei v mnogosloinykh struktu-rakh [Modeling of temperature fields in multilayer structures]. Programmnye produkty i sistemy = Software and Systems, 2008, no. 4, pp. 46-49.
9. Ivantsivsky V.V. Chislennoe modelirovanie temperaturnykh polei v materialakh pri uprochne-nii s ispol'zovaniem kontsentrirovannykh ob"emnykh istochnikov nagreva [Numerical modeling of temperature fields in materials during hardening using concentrated volumetric heating sources]. Nauchnyi vestnik Novosibirskogo gosudarstvennogo tekhnicheskogo universiteta = Science bulletin of the Novosibirsk state technical university, 2004, no. 2 (17), pp. 161-172.
10. Prisekin V.L., Rastorguev G.I., Belousov A.I. Analiz temperaturnykh polei i napryazhenii v elementakh konstruktsii pri vozdeistvii periodicheskikh teplovykh potokov [Analysis of temperature fields and stresses in structural elements under the influence of periodic heat flows]. Nauchnyi vestnik Novosibirskogo gosudarstvennogo tekhnicheskogo universiteta = Science bulletin of the Novosibirsk state technical university, 2001, no. 2 (11), pp. 133-142.
11. Bazilevich R.P. Dekompozitsionnye i topologicheskie metody avtomatizirovannogo kon-struirovaniya elektronnykh ustroistv [Decomposition and topological methods of automated design of electronic devices]. Lviv, Vishcha shkola Publ., 1981. 168 p.
12. Shtein M.E., Shtein B.E. Metody mashinnogoproektirovaniya tsifrovoi apparatury [Methods of machine design of digital equipment]. Moscow, Sovetskoe radio Publ., 1973. 296 p.
13. Morozov K.K., Odinokov V.G. Ispol'zovanie ETsVMpri konstruirovanii nekotorykh uzlov radioelektronnoi apparatury [The use of digital computers in the design of some components of radio-electronic equipment]. Moscow, Sovetskoe radio Publ., 1972. 104 p.
14. Segerlind L.J. Applied finite element analysis. New York etc., 1976 (Russ. ed.: Segerlind L. Primenenie metoda konechnykh elementov. Moscow, Mir Publ., 1979. 393 p.).
15. Pekhovich A.I., Zhidkikh V.M. Raschety teplovogo rezhima tverdykh tel [Calculations of the thermal regime of solids]. Leningrad, Energiya Publ., 1976. 352 p.
Для цитирования:
Озеркин Д.В. Расчет температурного поля многослойных несущих конструкций численными методами // Системы анализа и обработки данных. - 2022. - № 2 (86). - С. 105-120. -DOI: 10.17212/2782-2001-2022-2-105-120.
For citation:
Ozerkin D.V. Raschet temperaturnogo polya mnogosloinykh nesushchikh konstruktsii chis-len-nymi metodami [Computation of the temperature field of multilayer load-bearing structures by numerical methods]. Sistemy analiza i obrabotki dannykh = Analysis and Data Processing Systems, 2022, no. 2 (86), pp. 105-120. DOI: 10.17212/2782-2001-2022-2-105-120.
ISSN2782-2001, http://journals.nstu.ru/vestnik Analysis and data processing systems Vol. 86, No 2, 2022, pp. 105-120
