УДК 51-74
МАТЕМАТИЧЕСКИЕ МЕТОДЫ МОДЕЛИРОВАНИЯ ТЕПЛОВЫХ ПОЛЕЙ В ТРЕХМЕРНОЙ СБОРКЕ ИНТЕГРАЛЬНЫХ СХЕМ
В.Ф.Барабанов,М.Н. Аралов,Н.И.Гребенникова
В данной статье рассматриваются математические методы, применимые для расчета тепловых параметров и моделирования тепловых полей в трехмерной сборке интегральных схем
Ключевые слова: трехмерная сборка, тепловой анализ, функция Грина, разрывный метод Галёркина
Введение. Технология трехмерной сборки интегральных схем или 3Б-1С (рис. 1) является одним из наиболее перспективных методов, позволяющих повышать плотность упаковки, уменьшать задержку сигнала и обеспечивать конструкторскую гибкость. Однако с увеличением плотности упаковки и удельной мощности возникают тепловые проблемы, которые становятся критическим фактором в определении производительности, надежности и стоимости микросхемы. Известно, что чрезмерно высокая температура может значительно ухудшить надёжность межсоединений и устройства, и даже вызвать функциональные отказы посредством электротермической связи[1]. Следовательно, тепловой анализ необходим на системном этапе разработки устройст-ва[2]. Кроме того, тепловой анализ важен для размещения датчиков температуры и расчёта электромиграции.
Рис.1. Трехмерная сборка интегральных схем
В данной статье рассматриваются способы расчета температурных характеристик 3Б-1С на основе функции Грина и разрывного метода Галёркина.
Применение функции Грина.Управляющее уравнение в явлении теплопередачи для 3Б-1С выглядит следующим образом[3]:
теплопроводность в направлении У, К2 -теплопроводность в направлении 2, р- плотно-стькг/м3), с- удельная теплоёмкость материала (Дж/(кг-К)).
Для постоянной теплопроводности уравнение принимает следующий вид:
'д2Т д2Т д2Т ,дх2 ду2 dz2
дТ
= pcTt
или
д2Т д2Т д2Т дх2 ду2 дz2
1 дТ
(2)
Где рс- объёмная теплоёмкость, а=К/рс- коэффициент тепловой диффузии.
Обычно тепловой анализ выполняется на основе моделирования теплового распределения сложной структуры 3Б-1С (Рис. 2) и постепенного вычисления переходных процессов [3]. Через равные промежутки времени строятся карты распределения температуры для поверхности микросхемы и результаты сохраняются. Конечная модель может быть довольно большой, а для анализа переходных процессов требуется тысячи таких итераций, что приводит к сложностям разработки схемы.
Рис. 2. Тепловая модель в 3D-IC
¿(^К^Э + ^Э^ (1)
где T = T(x,y,z;t) -температурное распределение, Kx-теплопроводность в направлении Х (Вт/(м-К)), Ky -
Барабанов Владимир Федорович - ВГТУ, д-р техн. наук, профессор, тел. (473)257-28-50, e-mail: [email protected] Аралов Михаил Николаевич - ВГТУ, аспирант, e-mail: [email protected]
Гребенникова Наталия Ивановна - ВГТУ, канд. техн. наук, e-mail:[email protected]
Температура в определенной точке микросхемы может быть вычислена через интеграл свёрт-ки[3]при использовании функции Грина[4].
Интеграл свертки:
С1 дР(т)
Т(х,у,2,1) = I Т3(х,у,2,1 — т)—-—#т. (3) к дт
Здесь Т - температура, как функция расположения и времени 1и т; Р - распределение энергии, которое меняется в зависимости от параметра вре-
мени т; Т8 - значение функции Т при включении питания.
Зная Т8 , для данногоР температура Т в любой области (х,у,ъ) в любое время 1 может быть вычислена прямым способом. На рис. 3 показано сравнение результатов при использовании функции Грина и сценария среды ЛК8УБ[3].Алгоритм расчета тепловых параметров в среде ЛКБУБподробно изложен в [5].
34
32
о
о 30
я'
г<
р*. н 23
с:
А 26
В
й
он 24
22
га
¿Л
[/ 1 ------АТШГЯ
СтГннП
Э 50 1СО
Время, г
Рис. 3. Сравнение результатов вычислений
Метод на основе функции Грина даетопреде-ленную погрешность, но расчет ускоряется в десятки раз [3].
Разрывный метод Галёркина.Для оценки распределения температуры микросхема разделяется на более мелкие элементы[6]. Физические параметры, такие как теплопроводность и коэффициент теплопередачи, зависят от определенных свойств материала и примененных методов охлаждения. Граничные условия определяются операционной средой. Средство моделирования должно использовать геометрию расположения, граничные условия и физические тепловые параметры как начальные значения, чтобы сформулировать систему дифференциальных уравнений в частных производных, которые приводятся к системе уравнений в полных дифференциалах при помощи разрывного метода Галерки-на. Эти уравнения численно интегрируются последовательным способом, используя метод Рунге-Кутты [6].
Схематическое представление микросхемы и его тепловой модели показано на рис. 4.
Микросхема разделена на ячейки равномерно в х, у, и ъ направлениях, здесь 5х, 5у и 5ъ- размеры стороны каждой ячейки. Уравнение зависимости диффузии тепла через тепловую проводимость:
т)т] -о = о,
(4)
где Q - тепловой источник, Т - температура во время 1 Су-теплоёмкость объёма V, дгаёТ = [5Т/5х,5Т/5у, 5Т/5ъ], и матрица g является матрицей проводимости материала.
Для материалов, у которых есть три ортогональных направления различных тепловых прово-димостей, g=diag(g1); 1=х, у, gy, и gz являются тепловыми коэффициентами проводимостей.
Для решения этих уравнений, при применении численных методов, используется конечная дискретизация, т.е. модель микросхемы анализируется в многочисленные трехмерные элементы, где смежные элементы взаимодействуют через диффузию тепла. Каждый элемент достаточно маленький, чтобы выражать его температуру как дифференциальное уравнение, рассматривать его материальные характеристики, рассеяние мощности и температуры его соседних элементов. Каждой ячейке присвоена удельная теплоемкость. Если двойная сетка будет сформирована, присоединением к центрам смежных ячеек, то каждый край двойной сетки пересечет точно одну поверхность основной сетки. Теплопроводность присвоена краю двойной сетки. Если эти две ячейки по обе стороны от поверхности принадлежат тому же материалу, присвоенная теплопроводность - теплопроводность материала.
X
й X
Ъу
ы
V:
У.
Аж
Ьу
Ъх
Рис. 4. Дискретизация микросхемы
Интеграция частями в случае многомерного интеграла обобщена в теореме о дивергенции:
| + | (gгadз)*(gгad Т)
- I + | зк(Т - Т6)#7 = 0,
#4
(5)
где Та - температура окружающей среды, И - поверхностный коэффициент теплопередачи,определенный как Ь=1/(Лей-Я), где Лей- - эффективная область, нормальная к направлению потока тепла, и Я - эквива-
лентное тепловое сопротивление. Граничное условие Дирихле формы Т=0 (абсолютная температура, равная температуре окружающей среды). Граничное условие в 7=шт(7), как предполагается, имеет смешанный тип д25Т/57-ЬТ=0, где gz - теплопроводность в z направлении. Физически это соответствует тепловым потерям, являющимся пропорциональными различию между абсолютной температурой и температурой окружающей среды.
Если компонент градиента температуры вдоль z направления незначителен, то первоначально трехмерная модель упрощается до двух активных координат. Так как температура не меняется в зависимости от z, интегралы в (5) могут быть упрощены при предварительном интегрировании в направлении толщины, ёУ = AzdS и ёБ = AzdC. Объемные интегралы тогда оцениваются по площади поперечного сечения. Поверхностные интегралы вычислены как интегралы по контуру поперечного сечения Сс
^ЗСудТ-АгйБ + I (gгad з)g(gгad Т)тАг#Б — БС д Ьс
— I + I зк(Т — Та)Аг#Б = 0 (6)
Здесь з(х) = 0, для х Е Сс. Система уравнений полного дифференциала после дискретизации в пространстве (решения методом конечных элементов в форме Галеркина):
@СА-дд- + @С-АТ- = 1Н1+1=1+1К1
-=1 -=1
Где
С-А = I М(А)СуМ(-)Аг#Б ^<
С
= I
(7)
N(j)QAzdS
Lgj = - @ f (gradN(j))hg(grad Nj)TAzdS
С^, Gji, матрицы ёмкости и проводимости соответственно, Ьд - внутреннее тепловыделение.
Граничное условие в конечном виде дано как
i=Nf+l
Уравнение (7) может быть модернизировано и решено методом Рунге-Кутты[6]. Такой способ расчета тепловых параметров ресурсоемок, но дает хорошую точность.
Заключение.В данной статье описаны подходы к расчету тепловых параметров интегральных схем на основе функции Грина и разрывного метода Галёркина. Каждый подход может быть реализован программно и применим при создании специальных средств проектирования электронных компонентов. Материал статьи может помочь разработчику выбрать наиболее подходящий способ расчета тепловых параметров для решения конкретной задачи.
Литература
1. Брагин, Д.М. Математическое моделирование процесса размещения разногабаритных элементов на монтажном пространстве с учетом тепловых полей[Текст]/ Д.М.Брагин, В.Ф.Барабанов, А.М. Нужный //Вестник Воронежского государственного технического университета. - 2007.- Т. 3. - № 5. - С.76 -80.
2. Барабанов В.Ф., Аралов М.Н. - Структура программной модели цифрового устройства - Информационные технологии моделирования и управления (ИТМУ): Научно-технический журнал, № 6 (84), с. 583-589.
3. S. Pan and N.Chang,3D IC Dynamic Thermal Ana-lysi With Hierarchical And Configurable Chip Thermal Mode Slides L IEEE 3DIC 2013.
4. A. D. Polyanin and V. F. Zaitsev, Handbook of Exact Solutions for Ordinary Differential Equations (2nd edition), Chapman & Hall/CRC Press, Boca Raton, 2003. ISBN 1-58488-297-2
5. http://cae.ustu.ru/cont/soft/ansys.htm
6. A. Zjajo, N.P. van der Meijs and T.G.R.M. van Leuken; Thermal Analysis of 3D Integrated Circuits Based on DiscontinousGalerkin Finite Element Method; InProceedings of ISQED 2012.
7. Барабанов В.Ф.,Подвальный С.Л. Интерактивные средства моделирования сложных технологическихпро-цессов. Воронеж: ВГТУ, 2000.
8. Брагин Д.М. Интегрированный программный комплекс моделирования и проектирования электронных средств с учётом тепловых полей [Текст] / Д. М. Брагин, В. Ф. Барабанов, А. М. Нужный // Системы управления и информационные технологии. - 2007. - №3(29). - С. 63-66.
N
\L
1 — AI -X1 L
N(j)cvN(i)AzdS
dT9 ~dt
Воронежский государственный технический университет
MATHEMATICAL METHODS OF SIMULATION OF THERMAL FIELDS IN 3D-IC V.F.Barabanov,M.N.Aralov, N.I. Grebennikova
In this article mathematical methods, applicable for calculation of thermal parameters and simulation of thermal fields in three-dimensional integrated circuits packaging are considered
с
G
Keywords: three-dimensional assembly, thermal analysis, Green's function, discontinuous Galyorkin'smethod