Interface constraint expert hybrid model training method for 3D integrated circuit thermal evaluation
By employing a hybrid model training method based on interface constraints and experts, combined with a global backbone network and residual experts, the problems of heat flow mismatch and temperature oscillation in 3D integrated circuit thermal simulation were solved, achieving efficient and accurate temperature field prediction and improving the thermal evaluation capability of multi-material 3D ICs.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- SHENZHEN BIANGXIN TECH CO LTD
- Filing Date
- 2026-03-17
- Publication Date
- 2026-06-12
AI Technical Summary
Existing self-supervised operator learning methods in 3D integrated circuit thermal simulation suffer from several drawbacks. These include the lack of interface heat flow constraints due to the use of a global PDE loss function that assumes uniform thermal conductivity, the inability to capture the temperature gradient behavior on both sides of the material interface due to reliance on random or uniform sampling strategies, and the inability to express piecewise thermophysical properties due to the use of a single global neural network architecture that forces global smoothing of the temperature field. These issues lead to heat flow mismatch, temperature oscillation, and reduced overall temperature prediction accuracy.
An interface-constrained expert hybrid model training method is adopted, which includes a shared global backbone network and multiple residual experts. Paired sampling points at the interface are generated through a bilateral sampling strategy, and the bilateral normal temperature gradient is calculated. Interface heat flow continuity and temperature continuity loss are introduced into the training loss to construct an expert hybrid architecture to adapt to the piecewise thermal behavior of multi-material 3D ICs.
It significantly improves the accuracy and efficiency of thermal assessment in multi-material 3D ICs, enabling fast and physically consistent temperature field prediction. Compared with commercial solvers, it is 63.9 times faster and the average error is reduced to 0.28K.
Smart Images

Figure FT_1 
Figure FT_2 
Figure FT_3
Abstract
Description
Technical Field
[0001] This invention relates to the field of integrated circuit layout design and provides a hybrid model training method for interface constraint experts for thermal evaluation of 3D integrated circuits. Background Technology
[0002] Three-dimensional integrated circuits (3D ICs) significantly improve integration density by vertically stacking multiple chips, but this also leads to a sharp increase in power density, making the system more prone to hot spots and overheating. Thermal integrity (TI) has become a key challenge in 3D IC design [2-5]. Therefore, early layout planning and design space exploration (DSE) for thermal management have become indispensable parts of the modern chip design flow [6-10]. However, the lack of efficient and accurate thermal simulation tools constitutes a major bottleneck in the current design flow [11-13].
[0003] Traditional numerical solvers based on thermal partial differential equations (PDEs), such as the finite element method (FEM) and the finite difference method (FDM), can provide high-precision simulation results, but their computational cost is extremely high and their running time is too long, making it difficult to meet the efficiency requirements of iterative optimization workflows [14-17]. To accelerate the thermal analysis process, operator learning methods have been introduced into this field [18-21]. Among them, supervised operator learning relies on a large number of input-temperature field paired datasets generated by the complete numerical solver; however, generating such labeled data for 3D ICs is costly and difficult to cover the variable power distribution, material stacking combinations and boundary conditions in actual designs.
[0004] Self-supervised operator learning, by training directly on the control thermal PDE and boundary conditions, avoids dependence on real temperature labels [22-26] and shows potential for application in 3D IC thermal simulation. However, existing methods still face fundamental challenges in dealing with the inherent multi-material discontinuities in real 3D IC structures:
[0005] First, there is a sudden change in thermal conductivity at the material interface, while existing self-supervised methods [23,24] usually use a global PDE loss function that assumes uniform thermal conductivity, ignoring the discontinuity of the temperature gradient at the interface in heterogeneous stacked structures, leading to heat flow mismatch and temperature oscillation. Figure 2 To address this issue, we introduce an interface-aware loss function that explicitly enforces the continuity of heat flow between material boundaries, enabling the model to learn thermally consistent solutions in multi-material regions.
[0006] Secondly, random or uniform sampling strategies
[25] cannot effectively capture the bilateral gradient behavior at the material interface. Points on both sides of the interface have significantly different temperature gradients, which carry the key physical information of heat transfer between materials. Existing methods are insufficient in fitting the interface thermal behavior due to the defects of the sampling mechanism. To overcome this, we propose a bilateral sampling strategy, which draws geometrically aligned point pairs along the interface normal, so that the model can learn the bilateral temperature gradients and force the continuity of heat flow during training.
[0007] Finally, a single global operator cannot represent the piecewise PDE behavior inherent in multi-material 3D ICs. Although the temperature field is continuous, it is only piecewise smoothed: each material region follows its own thermal equation with unique thermal conductivity. The single global neural operator architecture [23,24] forces global smoothing of the temperature field, which is difficult to express the physical nature of “piecewise smoothing” in multi-material 3D ICs: although the overall temperature field is continuous, each material region follows its own thermal conduction equation with unique thermal conductivity. The global smoothing constraint will blur the material-specific thermal patterns and limit the prediction accuracy of the model in heterogeneous structures.
[0008] In summary, existing thermal simulation technologies have significant shortcomings in terms of computational efficiency, physical modeling capabilities of multi-material interfaces, and segmented thermal behavior representation. There is an urgent need for a rapid thermal evaluation method that can balance computational efficiency, physical consistency, and multi-material adaptability to support efficient design space exploration in 3D IC design. Summary of the Invention
[0009] The purpose of this invention is to solve the technical problems of existing self-supervised operator learning methods in 3D integrated circuit thermal simulation, such as heat flow mismatch at the material interface, temperature oscillation, inaccurate modeling of interface heat transfer, and significantly reduced overall temperature prediction accuracy. These problems arise from the use of a global PDE loss function that assumes uniform thermal conductivity, resulting in a lack of interface heat flow constraints; reliance on random or uniform sampling strategies that cannot capture the temperature gradient behavior on both sides of the material interface; and the use of a single global neural network architecture that forces global smoothing of the temperature field, which fails to express piecewise thermophysical properties.
[0010] To achieve the above objectives, the present invention employs the following technical means:
[0011] This invention provides a hybrid model training method for interface constraint experts for thermal evaluation of 3D integrated circuits, comprising:
[0012] Step 1: Obtain material domain information and power density distribution of 3D integrated circuits and boundary condition mapping The material domain information includes M material subdomains. The stacked structure consists of each subdomain Having piecewise constant thermal conductivity ;
[0013] Step 2: Construct a neural network model, which adopts an expert hybrid architecture, including a shared global backbone network. and multiple residual experts Each residual expert corresponds to a material subdomain in a three-dimensional integrated circuit that is separated by interfaces with discontinuous thermal conductivity.
[0014] Step 3: Generate paired sampling points at the interface of adjacent material subdomains using a bilateral sampling strategy. The paired sampling points are located within the adjacent material subdomains and aligned along the interface normal direction. Calculate the bilateral normal temperature gradient based on the paired sampling points.
[0015] Step 4: Calculate the training loss by weighting and summing the various losses according to their respective weights to form the total loss function. The training loss includes:
[0016] Volume loss Use the corresponding thermal conductivity Independent calculation of volumetric PDE loss This is used to constrain the temperature field within each material subdomain to satisfy the heat conduction equation;
[0017] Boundary loss Based on the thermal conductivity of the material at the boundary Calculate the corresponding boundary condition loss to constrain the boundaries of the three-dimensional integrated circuit to satisfy the thermal boundary condition;
[0018] Interfacial heat flow continuity loss The difference in normal heat flow is calculated based on the bilateral normal temperature gradient and the corresponding material thermal conductivity to constrain the continuity of interfacial heat flow and solve the temperature oscillation problem caused by the discontinuity of heat flow at the material interface.
[0019] Interface temperature continuity loss The temperature prediction difference between the paired sampling points is used to constrain the continuity of the interface temperature.
[0020] Step 5: Based on the total loss function Optimize the parameters of the neural network model;
[0021] Step 6: Input the material domain partitioning information, power density distribution data and boundary condition data of the three-dimensional integrated circuit to be evaluated into the trained neural network model, and output the complete thermal field prediction results.
[0022] In the above scheme, step 1 includes the following steps:
[0023] Step 1.1: Obtain material domain information, which includes decomposing the 3D integrated circuit stack into... Material subdomains ,in The total number of material subdomains is a positive integer. Indicates the first Each material subdomain Corresponding to different material regions, and with thermal conductivity fields Relatedly, the thermal conductivity field is a piecewise constant within the material subdomain but discontinuous at the material interface;
[0024] Step 1.2: Obtain power density distribution The power density distribution To describe the spatial distribution of volumetric heat generation within the chip, used to capture the power consumption of functional modules;
[0025] Step 1.3: Obtain boundary conditions The boundary conditions Thermal boundary conditions are defined on the domain boundary, including Neumann boundary conditions and Robin boundary conditions.
[0026] In the above scheme, step 2 includes the following steps:
[0027] In the above scheme, step 2, constructing the neural network model, includes the following steps:
[0028] Step 2.1: Construct a shared global backbone network The shared global backbone network uses spatial coordinates As input, the output basis function The basis functions are used to capture the overall temperature trend that changes smoothly throughout the entire 3D integrated circuit stack;
[0029] Step 2.2: For each material subdomain Building a residual expert The residual expert uses spatial coordinates and the thermal conductivity of the material subdomain As input, it is used to model the local thermal behavior correction amount within the corresponding material region, where For the first A constant thermal conductivity value within a material subdomain;
[0030] Step 2.3: Establish a routing mechanism so that for input spatial coordinates ,according to Material subdomain ,Will Routing to the corresponding residual expert ;
[0031] Step 2.4: Construct a temperature prediction module by adding the output of the shared global backbone network to the output of the routed residual expert to obtain the final temperature prediction value. ,satisfy:
[0032]
[0033] in, This represents the temperature field predicted by the neural network model. Represents the set of model parameters. This indicates that the output of the shared global backbone network is being shared. Indicates the first The output of a residual expert For spatial coordinates, For the first Material subdomains The thermal conductivity (constant within the subdomain). The first is separated by a thermal conductivity discontinuity interface. A material subdomain, , This represents the total number of material subdomains.
[0034] In the above scheme, step 3, which generates paired sampling points at the interface of adjacent material subdomains and calculates the bilateral normal temperature gradient components using a bilateral sampling strategy, includes the following steps:
[0035] Step 3.1: For each pair of adjacent material subdomains and Determine its shared interface and the shared interface Perform two-dimensional parameterization, in parameter coordinates The upper layer generates stratified sampling points, where and These are the two-dimensional coordinates of the interface.
[0036] Step 3.2: Set the parameter coordinates Mapping the layered sampling points on the surface to a three-dimensional Cartesian coordinate system yields the interface sampling points. ;
[0037] Step 3.3: For each interface sampling point Two offset points are generated along the interface normal as paired sampling points, and the two offset points are respectively located in the material subdomain. and Within, and satisfying the following relationship:
[0038]
[0039]
[0040] in, Indicates that it is located in the material subdomain Offset point within, Indicates that it is located in the material subdomain Offset point within, The preset offset distance, To from the material subdomain Pointing to material subdomain The interface unit normal vector;
[0041] Step 3.4: The paired sampling points and The data is input into a neural network model to obtain the corresponding predicted temperature field value. and ;
[0042] Step 3.5: Calculate the neural network model at paired sampling points using automatic differentiation. and Temperature gradient vector at and ;
[0043] Step 3.6: Calculate the two-sided normal temperature gradient components, i.e., calculate... With interface unit normal vector dot product ,as well as With interface unit normal vector dot product .
[0044] In the above scheme, step 4, calculating the training loss, includes the following steps:
[0045] Step 4.1: For each material subdomain ( ), calculate volumetric PDE loss The volumetric PDE loss Calculated using the following formula:
[0046]
[0047] in, Represents the gradient operator, For material subdomains thermal conductivity, For the neural network model in spatial coordinates The predicted temperature Let the volume heat source function be... The norm squared of the residual vector;
[0048] Step 4.2: For each boundary condition (total) (each), calculate boundary condition loss. The boundary condition loss Based on the thermal conductivity of the material at the boundary It can be calculated using the following formula:
[0049]
[0050] in, This represents the derivative of the temperature field along the outer normal of the boundary. For the specified heat flux density, The heat transfer coefficient is... For ambient temperature, and It is a boundary type indicator and satisfies ;
[0051] Step 4.3: Calculate the interfacial heat flow continuity loss The interfacial heat flow continuity loss Calculated using the following formula:
[0052]
[0053] in, This represents the total area of all material interfaces. Represents a material subdomain and Shared interface between ( ), and For the paired sampling points generated in step 3, and The two-sided normal temperature gradient represents the temperature gradient components along the interface normal at the paired sampling points. The interface area is a micro-element;
[0054] Step 4.4: Calculate the interface temperature continuity loss The interface temperature continuity loss Calculated using the following formula:
[0055]
[0056] in, and These represent the paired sampling points of the neural network model. and The predicted temperature value;
[0057] Step 4.5: Sum the losses according to their weights to form the total loss function. :
[0058]
[0059] in, , , , These are the weighting coefficients for the corresponding loss terms. The total number of material subdomains, This represents the total number of boundary conditions.
[0060] In the above scheme, step 5 is based on the total loss function. Optimizing the parameters of the neural network model includes the following steps:
[0061] Step 5.1: Employ a warm-up training strategy, optimizing only the shared global backbone network during the warm-up phase. The parameters are fixed for all domain-specific residual experts. The parameters enable the shared global backbone network to learn the overall temperature trend of smooth changes in the three-dimensional integrated circuit stack;
[0062] Step 5.2: After the warm-up phase is completed, remove the fixation on the domain-specific residual expert parameters, and simultaneously optimize the shared global backbone network. The parameters and all domain-specific residual experts The parameters are set to make the neural network model converge to the total loss function. The minimum value.
[0063] In the above scheme, step 6, which inputs the material domain partitioning information, power density distribution data, and boundary condition data of the three-dimensional integrated circuit to be evaluated into the trained neural network model and outputs the complete thermal field prediction result, includes the following steps:
[0064] Step 6.1: Obtain the material domain partitioning information of the three-dimensional integrated circuit to be evaluated. The material domain partitioning information includes dividing the three-dimensional integrated circuit into... Material subdomains ,in The total number of material subdomains is a positive integer. Indicates the first Each material subdomain It has thermal conductivity ;
[0065] Step 6.2: Obtain the power density distribution of the 3D integrated circuit to be evaluated. and boundary conditions ;
[0066] Step 6.3: Generate a uniformly distributed set of sampling points within the spatial domain of the three-dimensional integrated circuit. ,in The total number of sampling points. For the first The three-dimensional spatial coordinates of each sampling point;
[0067] Step 6.4: Power density distribution and boundary conditions The input is fed into the branch network of the trained neural network model to generate the encoded latent representation;
[0068] Step 6.5: For each sampling point The material subdomain to which it belongs is determined based on the material domain partitioning information. , sampling points and local thermal conductivity Input to shared global backbone network Obtaining basis functions and sampling points Thermal conductivity of its material subdomain Input to the corresponding residual expert ;
[0069] Step 6.6: Perform a single forward pass using the trained neural network model to calculate each sampling point. Temperature prediction at the location Among them, for those belonging to the material subdomain sampling points The predicted temperature value is calculated using the following formula:
[0070]
[0071] in, Indicates at the sampling point Predicted temperature at the location, Indicates sharing the global backbone network in The output at that location, Indicates the first A residual expert inputs spatial coordinates. and thermal conductivity Output at time;
[0072] Step 6.7: Compile all sampling points The temperature prediction value is output, providing a complete thermal field prediction result covering the entire three-dimensional integrated circuit stack.
[0073] Because the present invention employs the above-mentioned technical means, it has the following beneficial effects:
[0074] This invention constructs a neural network model with a hybrid expert architecture (step 2), including a shared global backbone network and multiple residual experts corresponding one-to-one with material subdomains. It also establishes a routing mechanism based on material domain affiliation to route input coordinates to the corresponding residual experts. The final temperature prediction is composed of the sum of the backbone network output and the output of the residual expert of the corresponding material subdomain. This solves the technical problem that a single global neural network architecture implicitly forces global smoothing of the temperature field and cannot express the physical essence of "segmented smoothing" in 3D integrated circuits (i.e., the overall temperature is continuous, but each material subdomain follows a heat conduction equation with different thermal conductivity, and there are abrupt changes in temperature gradient at the interface). The backbone network focuses on learning the smooth overall temperature trend across material stacks, while each residual expert focuses on modeling the local thermal behavior correction amount determined by a specific thermal conductivity within the corresponding material subdomain. The two work together to enable the model to characterize the thermophysical properties of different material regions, effectively avoiding the ambiguity of material-specific thermal patterns caused by global smoothing constraints, and significantly improving the accuracy of expressing segmented thermal behavior in heterogeneous material structures.
[0075] This invention employs a bilateral sampling strategy (step 3) to generate geometrically aligned paired sampling points (located in the two subdomains respectively) along the normal direction at the shared interface of adjacent material subdomains. Based on automatic differentiation, it calculates the bilateral normal temperature gradient components, solving the technical problem that traditional random or uniform sampling strategies cannot explicitly capture the temperature gradient differences on both sides of the interface, resulting in the loss of key physical information on interface heat transfer. This strategy deeply couples the interface geometry (normal direction) with the sampling process, enabling the model to directly perceive and learn the temperature gradient discontinuity caused by abrupt changes in thermal conductivity on both sides of the interface during training. This provides an accurate and aligned gradient calculation basis for subsequent interface physical constraints, enhancing the model's physical modeling ability for cross-material heat transfer processes.
[0076] This invention explicitly introduces interfacial heat flux continuity loss and interfacial temperature continuity loss into the training loss (in step 4). and The former calculates the heat flow difference based on the bilateral normal temperature gradients and corresponding thermal conductivity of the paired points obtained from bilateral sampling, while the latter constrains the difference in predicted temperature values at the paired points. This solves the technical problem of existing self-supervised methods relying on global PDE loss and ignoring the physical laws of heat flow conservation and temperature continuity at the material interface, leading to interface heat flow mismatch and non-physical temperature oscillations. The loss design strictly follows the basic principles of thermodynamics, and forces the model to satisfy the dual constraints of heat flow continuity (energy conservation) and temperature continuity (physical continuity) at the interface during the optimization process. It eliminates the numerical distortion caused by the discontinuity of thermal conductivity from the perspective of the loss function, ensuring the physical consistency and stability of the predicted temperature field in the interface region.
[0077] This invention solves the technical problem of traditional methods that use a uniform thermal conductivity to calculate global loss and cannot adapt to the independent thermophysical equations and boundary conditions of each subdomain in multi-material structures by independently calculating the volume PDE loss for each material subdomain (step 4.1) and calculating the boundary condition loss based on the thermal conductivity of the boundary material (step 4.2). This design ensures that the loss function is strictly aligned with the physical model: the volume loss ensures that the temperature field inside each material subdomain satisfies the heat conduction equation under its own thermal conductivity, and the boundary loss ensures that the boundary heat exchange behavior matches the local material properties. Thus, the physical constraints are accurately applied both inside the subdomain and at the system boundary, comprehensively improving the local accuracy and overall physical rationality of the complete thermal field prediction.
[0078] This invention employs a preheating training strategy (step 5.1) to initially fix the residual expert parameters and optimize only the shared backbone network. Once the backbone network converges to the overall temperature trend, all parameters are then jointly optimized. This solves the training conflict and convergence instability issues that may arise between the backbone network and residual experts in expert hybrid models due to differences in optimization objectives. This strategy, through phased optimization, allows the backbone network to first establish a reasonable basis for the global temperature field, providing a stable reference for the subsequent refined local corrections by residual experts. It effectively coordinates the optimization tension between global consistency and local specificity, significantly improving the stability and convergence efficiency of the model training process.
[0079] This invention achieves significant synergistic effects by organically combining an expert hybrid architecture, a bilateral sampling strategy, an interface constraint loss function, and a preheating training strategy: the expert hybrid architecture provides the structural foundation for piecewise thermal behavior modeling, enabling the backbone network and residual experts to perform their respective functions; the bilateral sampling strategy provides geometrically aligned bilateral gradient inputs for the interface constraint loss, enabling... and The calculations are grounded in physical principles; the interface constraint loss guides gradient learning at bilateral sampling points, ensuring that expert outputs satisfy thermodynamic continuity at the interface; and the preheating training strategy further guarantees the stable convergence of this complex architecture. These four elements are interconnected—the architecture design supports the implementation of physical constraints, the sampling strategy enhances the accuracy of constraint calculations, the loss function drives the model to learn physical laws, and the training strategy optimizes the learning process—together solving complex technical challenges that cannot be overcome by single techniques, such as the contradiction between "global smoothness and local discontinuity," "missing interface physical information sampling," and "conflict between interface constraints and internal PDE optimization." This allows the model to accurately characterize the complete thermophysical behavior of multi-material 3D integrated circuits—characterized by "overall continuity, segmented smoothness, and interface abrupt changes"—while maintaining the physical consistency of the complete thermal field. This represents a fundamental improvement in the physical fidelity and modeling adaptability of thermal evaluation methods. Attached Figure Description
[0080] Figure 1 This is a simplified flowchart of the present invention;
[0081] Figure 2 This diagram illustrates the thermophysical phenomena at multi-material interfaces. Single-domain partial differential equations cannot accurately describe the thermal behavior at interfaces. The IC-MoE model establishes dedicated partial differential equations for each material and introduces interface constraints, thereby achieving physical consistency in the cross-layer heat transfer process.
[0082] Figure 3 This is the IC-MoE framework. Left side: Input data of multi-material 3D integrated circuits processed using a bilateral normal alignment sampling strategy; Top: Branch network encoding power distribution maps and boundary conditions; Center: Domain-aware expert hybrid architecture, which organically combines a shared global backbone network with residual experts for each corresponding material subdomain.
[0083] Figure 4 This represents the temperature distribution along the z-direction under the influence of a uniform heat source.
[0084] Figure 5 The results are qualitative assessments. From left to right: input power distribution, actual temperature field, prediction results from this model, and corresponding relative error distribution.
[0085] Figure 6 Training curves for a single-network model (left) and an expert hybrid model (right). Red: Total loss; Blue: Squared error of interface heat flux. Detailed Implementation
[0086] The embodiments of the present invention will be described in detail below. Although the present invention will be described and illustrated in conjunction with some specific embodiments, it should be noted that the present invention is not limited to these embodiments. On the contrary, any modifications or equivalent substitutions made to the present invention should be covered within the scope of the claims of the present invention.
[0087] Furthermore, to better illustrate the present invention, numerous specific details are set forth in the following detailed embodiments. Those skilled in the art will understand that the present invention can be practiced without these specific details.
[0088] This invention proposes IC-MoE, an interface-constrained, self-supervised operator-learning MoE framework for fast and physically consistent thermal analysis of 3D integrated circuits (ICs). Our method integrates material-specific PDE constraints, interface-aware sampling, and a domain-aware MoE architecture, effectively capturing global temperature trends and local multi-material thermal behavior. The main contributions are summarized below:
[0089] Interface-aware heat flow formula: We derived a physically consistent heat flow loss function, explicitly enforced the continuity of heat flow between material interfaces, and achieved accurate modeling of thermal conductivity discontinuities in multi-material stacks.
[0090] Bilateral Sampling: We introduce a sampling strategy to generate aligned point pairs across material boundaries, enabling true learning of bilateral gradients and cross-material heat transfer.
[0091] Domain-Aware MoE: We propose a MoE architecture that combines a global backbone network with domain-specific residual experts to capture the piecewise smooth temperature field unique to 3D ICs.
[0092] Efficient and accurate thermal assessment: IC-MoE achieves label-free, physically consistent temperature prediction, 63.9 times faster than commercial solvers, while maintaining an average error of only 0.28K in multi-layer stacked cases.
[0093] To facilitate a better understanding of the technical concept of this invention by those skilled in the art, the relevant technologies involved in this invention will be further described in detail:
[0094] 1.1 Thermal Modeling of 3D ICs
[0095] Given heat sources, material properties, and boundary conditions, thermal simulation in 3D ICs aims to predict the temperature distribution throughout the package. The core task is to solve for the temperature distribution at each location within the domain. Partial differential equation for heat conduction:
[0096] (1)
[0097] in It's temperature. Indicates a volumetric heat source. , and These are the material density, heat capacity, and thermal conductivity, respectively. In this work, we focus on the steady-state case, where the time derivative disappears, resulting in the Poisson equation:
[0098] (2)
[0099] In practical IC thermal modeling, two types of boundary conditions are typically used:
[0100] Neumann boundary conditions define the normal temperature gradient, corresponding to a specific heat flux entering or leaving the surface.
[0101] Neumann boundary conditions. Define the normal temperature gradient, corresponding to a specific heat flux entering or leaving the surface.
[0102] (3)
[0103] The adiabatic boundary is Special circumstances.
[0104] Robin boundary condition. This boundary condition simulates the balance between internal conduction and convective cooling of the solid to the surrounding environment.
[0105] (4)
[0106] in It is the heat transfer coefficient. It refers to the ambient temperature.
[0107] 1.2 Operator Learning for Partial Differential Equations (PDEs)
[0108] Operator learning aims to approximate nonlinear mappings between function spaces, rather than point-to-point input-output mappings. Given a family of hot PDEs, the goal is to approximate the solution operators:
[0109] (5)
[0110] in It is the input function (e.g., boundary conditions or material mapping). This represents the corresponding temperature field.
[0111] DeepONet
[27] is a fundamental operator learning framework that employs a branch-backbone structure. The branch network samples and encodes the input function u, generating expansion coefficients; while the backbone network embeds query coordinates and outputs basis functions. The final solution is obtained by combining the coefficients and basis values.
[0112] DeepOHeat
[23] extends DeepONet
[27] for steady-state thermal analysis in 3D ICs by training operators using only controlled PDEs and boundary constraints, thus eliminating the need for real temperature labels.
[0113] PDE residual loss. At a point within the material domain, the predicted temperature field must satisfy equation (2). PDE residual loss penalizes violations:
[0114] (6)
[0115] The driving network satisfies the steady-state thermal equation.
[0116] Boundary condition losses. On exposed surfaces, heat exchange follows Neumann or Robin type conditions. The corresponding losses are defined as:
[0117] (7)
[0118] The first item mandates heat flow (e.g., power distribution map), and the second item mandates convection behavior. These are indicators of the boundary type. Together, these losses ensure that the predicted temperature field is physically consistent with the boundary heat transport.
[0119] 1.3 Problem Statement
[0120] Despite recent advances in operator learning for 3D thermal simulation, no self-supervised framework can yet handle the complex multi-material structures encountered in real-world 3D ICs. In this study, we formulate steady-state multi-material thermal simulation as an operator learning problem. We consider integrated circuits composed of multiple materials, such as silicon, copper, and dielectrics, each exhibiting distinct thermal properties. We define the input as follows:
[0121] Materials Domain The 3D IC stack is broken down into: Each subdomain corresponds to a different material region. With thermal conductivity field Relatedly, the field is a piecewise constant within the material, but discontinuous at the interface.
[0122] Power density distribution : Describes the spatial map (spatial mapping) of volumetric heat generation within the chip, and captures the power consumption of functional modules.
[0123] Boundary condition diagram Thermal boundary conditions specified on the domain boundary include heat flow (Neumann) and convection constraints (Robin).
[0124] The desired output is the steady-state temperature field governed by the heat conduction equations introduced in Section 1.1. Formally, our goal is to learn a neural operator that approximates the fundamental heat map:
[0125] (8)
[0126] in This represents the solver for the Poisson equations for multiple materials.
[0127] The technical concept of the interface constraint expert hybrid model training method for thermal evaluation of 3D integrated circuits provided by this invention is described in detail below:
[0128] 2.1 Overview
[0129] The overall framework of IC-MoE is as follows: Figure 3 As shown. IC-MoE follows the branch-backbone formula commonly used in DeepONet-style operator learning architectures [23, 25, 27]. The branch network encodes the input function, including the power map and boundary conditions, generating a latent representation that captures its global effects. The backbone network takes spatial coordinates and local thermal conductivity ⟨x, κ> as input and outputs the basis function evaluated at that location. The final temperature prediction is obtained as follows:
[0130] (9)
[0131] in These represent the coefficients of the encoded input function from the branch network. This represents the main basis functions.
[0132] The contribution of this invention's technical solution integrates three key components:
[0133] 1) Interface-aware heat flow loss (Section 2.2) separates intra-domain PDE enforcement from inter-domain heat flow matching, enabling physically correct heat transfer between discontinuous thermal conductivities;
[0134] 2) Two-sided heat transfer modeling (Section 2.3): Generate paired interface samples along the normal to calculate the two-sided gradient;
[0135] 3) Domain-aware MoE architecture (Section 2.4) combines a global backbone network with domain-specific residual experts to capture piecewise thermal behavior. These components together enable IC-MoE to provide fast, accurate, and physically consistent temperature predictions in heterogeneous 3D IC stacks.
[0136] 2.2 Interface-Sensitive Heat Flow Loss
[0137] A key challenge in self-supervised thermal learning is handling multi-material domains, where thermal conductivity... It is a piecewise constant but discontinuous at the interface. The conventional PDE loss (Equation (6)) assumes a continuous thermal conductivity field and enforces a single-domain Poisson residual (Equation (2)). However, at the material interface, the temperature derivative is discontinuous, and directly applying the global residual would violate the conservation of heat flux.
[0138] To address this issue, this invention introduces an interface-sensing heat flow loss mechanism that explicitly respects the physical properties of multiple materials. The domain is defined by M sub-regions. Composition, each subregion has a constant thermal conductivity Within each material, the PDE must be independent:
[0139] (10)
[0140] Extending κ(x) to the entire domain produces a distribution divergence:
[0141] (11)
[0142] in Represents the constraints of the heat conduction equation within a material subdomain (volume PDE residual). This indicates the heat flow discontinuity at the material interface. It is the interface between materials. It is a domain The indicator, if ,but Otherwise it is 0, it is only in the first... Activate PDE residuals within each region. The Dirac measure on the interface is used as... The indicator. Therefore, the residual is decomposed into two physically distinct components: the volumetric PDE residual that controls the thermal diffusion within each material, and the interface heat flux residual that enforces conservation across discontinuities. We therefore define interface-sensored heat flux loss:
[0143] (12)
[0144] Furthermore, the temperature must remain continuous at the interface. We enforce this requirement through a value continuity loss:
[0145] (13)
[0146] During training, sampling points are strictly drawn within the subdomain to avoid pseudo-cross-interface gradients. Our final training loss is:
[0147] (14)
[0148] in It is in the The PDE loss within each domain is calculated as shown in equation (6). It is the first The loss due to each boundary condition.
[0149] 2.3 Two-sided heat transfer modeling
[0150] At the material interface, the temperature gradient decomposes into tangential and normal components:
[0151] (15)
[0152] Where n is the interface normal. , Form orthogonal tangential bases.
[0153] To ensure accurate cross-material heat transfer, we generate paired and aligned samples around the interface. For interface points... This generates two offset points, each within the adjacent material:
[0154] (16)
[0155] in It is a small offset distance. Indicates from domain arrive The outward unit normal. These two offset points represent the local neighborhood on either side of the interface. The neural network predicts the temperature field at these two locations. Using automatic differentiation, we obtain the gradient. and Calculate the normal heat flow from it:
[0156] (17)
[0157] This quantifies the mismatch in heat transfer across materials and is minimized to force heat flow continuity.
[0158] To support bilateral heat conduction modeling during training, the model must receive paired and geometrically aligned sampling points at each interface to enable normal offset evaluation and bilateral gradient calculation. To this end, this invention designs a specialized interface-aware sampling strategy that generates normalized sample pairs for each pair of adjacent material regions.
[0159] Algorithm 1 summarizes the interface-aware sampling process. For each pair of adjacent material domains, we first generate a set of hierarchical samples on the 2D parameterization of the shared interface. These uniformly distributed... The samples are then mapped back to 3D Cartesian coordinates and displaced along the interface normal to form a pair of offset points, resulting in two geometrically aligned samples located in adjacent domains.
[0160] By aligning sampling and differentiation along the interface normal, our two-sided scheme directly embeds geometric information into the training process. It acts as a directional physical regularizer, ensuring energy-conserving heat transfer between discontinuous thermal conductivities while preserving the natural tangential variations within each material.
[0161] Algorithm 1: Interface-aware sampling between adjacent domains
[0162] Input: A 3D stack with multiple material domains Number of sampling points per interface
[0163] Output: Set of sample point pairs
[0164] For each located and Interface between ,implement:
[0165] 1.
[0166] / * Create Each hierarchical interval /
[0167] 2.
[0168] 3.
[0169] 4. For each ,implement:
[0170] 5.
[0171] 6.
[0172] 7.
[0173] 8.
[0174] / Map UV parameter domain to Cartesian coordinates* /
[0175] 9. For each ,implement:
[0176] 10.
[0177] 11.
[0178] 12.
[0179] 13.
[0180] Among them, sampling points lie in and The algorithm obtains the result after shuffling. , Indicates from material subdomain point to The interface unit normal vector, This represents a function to randomly shuffle the data. This represents the set of n consecutive sub-intervals that divide the interval [0,1] into equal parts along the direction of the parameterized domain 𝑢. It is specifically used to implement a hierarchical sampling partitioning structure. This represents an ordered list of *k* sampled parameter values generated through hierarchical sampling in the *u* direction of the interface parameterization domain. In the 𝑣 direction of the interface parameterization domain, an ordered list of 𝑁 original sampled parameter values generated through hierarchical sampling is formed. Rand_sample(·) represents the random sampling function, and Rand_shuffle(·) represents the random shuffling function. This represents the shuffled list of parameters for 𝑣. Represents the set of parameter point pairs (𝑢,𝑣). This represents the set of *k* equally wide hierarchical intervals divided along the *k* direction (orthogonal to *k*) in the two-dimensional parameterized domain of the interface. *uv_to_Cartesian(·)* represents the parametric coordinate-Cartesian coordinate mapping function. This indicates a preset small offset distance. This represents the assignment operator, and 𝑁 represents the number of sample point pairs for each interface. , Indicates the coordinates of the paired sampling points;
[0181] 2.4 Proposed Expert Hybrid Architecture
[0182] Modeling heat conduction in multi-material 3D ICs requires capturing two different behaviors:
[0183] (i) Within each material, the temperature field is smooth and controlled by its own thermal conductivity;
[0184] (ii) At material interfaces, temperature gradients may change abruptly, while heat flow must remain continuous. A single, monolithic neural operator implicitly enforces global smoothness. Therefore, it struggles to represent the piecewise structure of real PDE solutions. To address these limitations, we propose a domain-aware MoE architecture specifically tailored to the physical properties of multi-material thermal systems.
[0185] A global backbone network with domain residual experts. Our MoE consists of a shared global backbone network. and a group of domain-specific residual experts The backbone network processes each spatial coordinate, capturing the overall temperature trend that changes smoothly throughout the stack. Each residual expert is then responsible for modeling local corrections within its respective material region, enabling the network to learn thermal conductivity-dependent variations that the backbone network alone could not represent.
[0186] Routing is determined by the material identity of each sampling point. Belongs to domain. The point was only forwarded to experts. The final temperature prediction calculation is as follows:
[0187] (18)
[0188] Interface constraint modeling. The MoE architecture is naturally integrated with our interface processing mechanism (Sec. 3.3). Bilateral interface samples are routed to the corresponding material experts, enabling the model to compute bilateral gradients that capture thermal conductivity discontinuities. These expert outputs are then incorporated into the interface thermal flux consistency and value continuity losses. Simultaneously, the global backbone network maintains the consistency of the overall temperature field. To improve training stability, we employ a warm-up phase during which only the backbone network is optimized, allowing it to learn the global solution structure before domain-specific experts refine material-dependent details.
[0189] 3 Experiments
[0190] 3.1 Experimental Setup
[0191] Training and test cases. We build a dataset based on a three-layer stacked 3D IC structure. The footprint is 1 mm², and the layer thicknesses are 0.2 mm, 0.1 mm and 0.2 mm, with corresponding thermal conductivities of 0.1, 1.0 and 0.1 W / (m·K). A power mapping is applied to the top surface, and the power density is sampled from a Gaussian mixture during training and distributed as a uniform value over the tiling area during testing, with the peak power density set to 2.5 × 10³ W / m²
[23] . The sides are modeled as adiabatic, while the bottom surface follows a convective boundary condition, h = 500 W / (m²·K), and the ambient temperature is 298.15 K. We generate ten test cases with different power mappings to evaluate the generalization performance of the proposed method.
[0192] Evaluation Metrics. We use three metrics to evaluate prediction accuracy: Mean Absolute Error (MAE), Mean Absolute Peak Error (MAPE), and Mean Absolute Relative Error (MARE). MAE measures the mean absolute temperature deviation across all spatial points, while MARE normalizes this deviation using the corresponding true temperatures. For each test case, we calculate the absolute difference between the predicted peak temperature and the true peak temperature. We report the MAPE obtained by averaging this peak temperature error across all test cases.
[0193] Implementation Details. All experiments were conducted on a Linux server equipped with an NVIDIA A100 (80GB) GPU and a 48-core Intel Xeon Silver4130 CPU (running at 3.3 GHz, with 128 GB of memory). The global backbone network and all residual experts in IC-MoE share the same architecture: the input coordinates are first encoded using 64 random Fourier features, followed by a 5-layer feedforward network with 128 hidden dimensions and SiLU activations, and finally a linear output layer. Training used 11,000 samples per iteration, including 4,000 volumetric samples for PDE loss, 3,000 boundary samples for boundary condition loss, and 4,000 interface samples generated by the proposed interface-aware sampling strategy. The loss weights for PDE, boundary, and heat flow losses were all set to 1.0.
[0194] 3.2 Comparison with State-of-the-art Methods
[0195] To evaluate the accuracy of IC-MoE under realistic multimaterial thermal conditions, we compare it with state-of-the-art self-supervised operator learning methods for 3D IC thermal simulation, including DeepOHeat
[23] and DeepOHeat-KAN
[24] . For a fair comparison, we modify both baselines to support the thermal conductivity of heterogeneous materials. In particular, the PDE residuals in each method are calculated using the local thermal conductivity of each domain rather than a single global value, ensuring that all methods operate under the same physical assumptions.
[0196] The quantitative results are summarized in Table 1. IC-MoE achieved significantly lower errors across all metrics, reducing MAE from 0.621K to 0.280K, MAPE from 2.520K to 1.144K, and MARE from 0.203 to 0.092, compared to DeepOHeat
[23] . Compared to DeepOHeat-KAN
[24] , IC-MoE achieved similar or even greater improvements, including a 59.4% reduction in MAE and a 74.7% reduction in MAPE. Figure 5 Further results show that the predicted temperature field matches the real-world situation very well in terms of thermal distribution and absolute value. These results demonstrate that the combination of domain decomposition, interface-aware sampling, and boundary-consistent heat flow modeling enables IC-MoE to capture the piecewise thermal behavior of multi-material conduction more faithfully than existing self-supervised neural operators.
[0197] Table 1: Comparison with state-of-the-art methods. Benchmark: COMSOL
[0198]
[0199] To further evaluate the ability of each model to handle material discontinuities, we constructed an additional case where the top surface is subjected to a spatially uniform heat source. Under this simple excitation, the temperature varies only along the z-direction, making the effect of thermal conductivity jumps particularly pronounced in the three-layer stack. Figure 4 As shown, the true solution exhibits a significant change in temperature gradient at each material interface, consistent with the discontinuous thermal conductivity.
[0200] IC-MoE accurately reproduces these gradient discontinuities, closely matching the true temperature distribution at all layer boundaries. In contrast, DeepOHeat fails to capture abrupt slope changes because it relies solely on the global PDE residuals without explicitly handling the heterogeneous material; its optimization implicitly forces an overly smoothed temperature field and suppresses gradient transitions induced by interfaces. This structural constraint results in a persistent mismatch around the interfaces and explains the significantly larger maximum error observed by DeepOHeat in Table 1.
[0201] 3.3 Efficiency Evaluation
[0202] To evaluate the computational efficiency of our method, the inference time of IC-MoE was compared with several widely used thermal solvers, including HotSpot
[28] , MTA
[29] , and COMSOL
[30] . For each method, we measured the clock time required to compute the same three-layer stacked steady-state temperature field used in our accuracy experiments. For a fair comparison of runtime, IC-MoE was evaluated on 257,091 uniformly sampled spatial points (71×71×51), chosen to closely match the 247,961 degrees of freedom used in COMSOL. For HotSpot, we used its default mesh configuration, resulting in a comparable resolution of 262,144 elements. Due to mesh limitations in MTA, we chose a setting with 151,905 degrees of freedom, which is closest to the available configurations of the other solvers, ensuring a comparison of all methods at similar problem scales.
[0203] Table 2 summarizes the runtime of each solver and the corresponding speedup achieved by IC-MoE. Our method offers a speedup of 269.5 times compared to Hotspot, 60.0 times compared to MTA, and 63.9 times compared to COMSOL. With nearly the same number of degrees of freedom, our method completes inference in 0.47 seconds, while COMSOL requires approximately half a minute. This significant speedup stems from IC-MoE's use of feedforward neural operators instead of iterative numerical solutions, allowing steady-state temperatures to be calculated in a single forward pass without matrix assembly or iterative convergence.
[0204]
[0205] 4.4 Effectiveness of the MoE Architecture
[0206] To validate the necessity of the domain-aware MoE architecture, we compared IC-MoE with a baseline that replaced all domain-specific experts with a single monolithic network. Figure 6 Training curves for the total loss and squared interface heat flux error (heat flux loss) are shown for two configurations. Figure 6 In the left figure, the single-network baseline exhibits a clear conflict between PDE loss and heat flux loss. The interfacial heat flux error decreases only briefly before rising again, eventually converging to a relatively high plateau. As training progresses, the two curves begin to overlap, indicating that the optimizer cannot reduce the heat flux mismatch without increasing the overall PDE residual. This behavior reveals the fundamental limitation of the monolithic architecture: because the model enforces global smoothness, it cannot represent the required temperature gradient discontinuities at material interfaces. The network is forced to compromise between satisfying the PDE within each material and maintaining heat flux continuity between boundaries, ultimately failing at both.
[0207] In comparison, Figure 6 The right figure shows that IC-MoE is stable throughout training and simultaneously reduces both total loss and heat flux loss. Domain-specific experts learn smooth, thermally consistent temperature fields within each material, while the global backbone network maintains cross-domain consistency. This decomposition provides the representational flexibility needed to represent piecewise thermal behavior, including sharp gradient transitions at interfaces, without compromising accuracy in homogeneous regions. Therefore, heat flux loss decreases consistently along with total loss, demonstrating that IC-MoE can faithfully enforce both PDE residuals and heat flux continuity simultaneously.
[0208] 5. Conclusion
[0209] In this paper, we propose IC-MoE, an interface-constrained MoE framework tailored for fast and accurate thermal evaluation in multi-material 3D ICs. By introducing an interface-aware heat flux formula, IC-MoE enforces physically correct heat flux continuity across material boundaries. Through a bilateral, normal-aligned sampling strategy, the model effectively captures temperature gradient discontinuities ignored by conventional sampling. The domain-aware MoE architecture enables IC-MoE to learn piecewise smooth thermal behavior in heterogeneous materials while maintaining global consistency across the stack. Comprehensive experiments demonstrate that IC-MoE is up to 63.9 times faster than commercial solvers inference and reduces temperature error by over 59.4% compared to state-of-the-art operator-learned baseline methods, highlighting its effectiveness and practicality in thermal analysis of multi-material 3D ICs.
[0210] References
[0211] [1] Swapnil S Salvi and Ankur Jain. 2021. A review of recent researchon heat transfer in three dimensional integrated circuits (3-D ICs). IEEETransactions on Components, Packaging and Manufacturing Technology 11, 5(2021), 802–821.
[0212] [2] Tao Li, Jie Hou, Jinli Yan, Rulin Liu, Hui Yang, and Zhigang Sun.2020. Chiplet heterogeneous integration technology—status and challenges.Electronics 9, 4 (2020), 670:1–670:12.
[0213] [3] YoonJo Kim, Yogendra K Joshi, Andrei G Fedorov, Young-Joon Lee,and Sung-Kyu Lim. 2010. Thermal characterization of interlayer microfluidiccooling of three-dimensional integrated circuits with nonuniform heat flux.Journal of Heat Transfer 132, 4 (2010).
[0214] [4] Shunxiang Lan, Min Tang, Liang Chen, and Junfa Mao. 2024. Thermalresistance network derivative (TREND) model for efficient thermal simulationand design of ICs and packages. In ACM / IEEE Design Automation Conference(DAC). 1–6.
[0215] [5] Hao Yu, Joanna Ho, and Lei He. 2009. Allocating power ground viasin 3D ICs for simultaneous power and thermal integrity. ACM Transactions onDesign Automation of Electronic Systems (TODAES) 14, 3 (2009), 1–31.
[0216] [6] Siyuan Miao, Lingkang Zhu, Wenkai Yang, Teng Lu, Yanze Zhou, ChenWu, Zhiping Yu, Ting-Jung Lin, and Lei He. 2025. Electrothermal simulationand vertical interconnect planning for integrated chiplets. In InternationalSymposium of Electronics Design Automation. 705–711.
[0217] [7] Jason Cong, Guojie Luo, and Yiyu Shi. 2011. Thermal-aware celland through-silicon-via co-placement for 3D ICs. In ACM / IEEE DesignAutomation Conference (DAC). 670–675.
[0218] [8] Jaeha Kung, Inhak Han, Sachin Sapatnekar, and Youngsoo Shin.2011. Thermal signature: A simple yet accurate thermal index for floorplanoptimization. In Proceedings of the 48th Design Automation Conference. 108–113.
[0219] [9] Qipan Wang, Xueqing Li, Tianyu Jia, Yibo Lin, Runsheng Wang, andRu Huang. 2024. ATPlace2.5D: Analytical thermal-aware chiplet placementframework for large-scale 2.5D-IC. In IEEE / ACM International Conference onComputer-Aided Design (ICCAD). 1–9.
[0220]
[10] Da-Wei Wang, Bo-Wen Zhang, Le-Tian Wang, Peng Zhang, and Wen-Sheng Zhao. 2025. A Proposal of Fast Thermal Simulation Method For 2.5-DAdvanced Packaging to Enable Efficient Thermal-Aware Placement Optimization.IEEE Transactions on Computer-Aided Design of Integrated Circuits and Systems(2025).
[0221]
[11] Zao Liu, Sahana Swarup, Sheldon X-D Tan, Hai-Bao Chen, and HaiWang. 2014. Compact lateral thermal resistance model of TSVs for fast finite-difference based thermal analysis of 3-D stacked ICs. IEEE Transactions onComputer-Aided Design of Integrated Circuits and Systems 33, 10 (2014), 1490–1502.
[0222]
[12] Zhixuan Dong, Longlong Yang, Cuiyang Ding, Changhao Yan, ZhaoriBi, Sheng-Guo Wang, Dian Zhou, and Xuan Zeng. 2024. pPIRW: An Efficient andAccurate Pre-calculation Path Integral Random Walk Solver for Steady-StateThermal Simulation With Robin Boundary Conditions. IEEE Transactions onComputer-Aided Design of Integrated Circuits and Systems (2024).
[0223]
[13] Peng Li, Lawrence T Pileggi, Mehdi Asheghi, and Rajit Chandra.2006. IC thermal simulation and modeling via efficient multigrid-basedapproaches. IEEE Transactions on Computer-Aided Design of Integrated Circuitsand Systems 25, 9 (2006), 1763–1776.
[0224]
[14] Tianxiang Zhu, Qipan Wang, Yibo Lin, Runsheng Wang, and RuHuang. 2024. FaStTherm: Fast and Stable Full-Chip Transient Thermal PredictorConsidering Nonlinear Effects. In IEEE / ACM International Conference onComputer-Aided Design (ICCAD). 1–9.
[0225]
[15] Haiyang He, Norman Chang, Jie Yang, Akhilesh Kumar, Wenbo Xia,Lang Lin, and Rishikesh Ranade. 2023. Solving fine-grained static 3dicthermal with ml thermal solver enhanced with decay curve characterization. InIEEE / ACM International Conference on Computer Aided Design (ICCAD). IEEE, 1–7.
[0226]
[16] Liang Chen, Jincong Lu, Wentian Jin, and Sheldon X-D Tan. 2023.Fast full-chip parametric thermal analysis based on enhanced physics enforcedneural networks. In IEEE / ACM International Conference on Computer AidedDesign (ICCAD). IEEE, 1–8.
[0227]
[17] Yu-Min Lee, Chi-Wen Pan, Pei-Yu Huang, and Chi-Ping Yang. 2015.LUTSim: A look-up table-based thermal simulator for 3-D ICs. IEEETransactions on Computer-Aided Design of Integrated Circuits and Systems 34,8 (2015), 1250–1263.
[0228]
[18] Nikola B Kovachki, Samuel Lanthaler, and Andrew M Stuart. 2024.Operator learning: Algorithms and analysis. Handbook of Numerical Analysis 25(2024), 419–467.
[0229]
[19] Yue Hua, Zhi-Qiao Wang, Xin-Yi Yuan, Yu Bai Li, Wei-Tao Wu, andNadine Aubry. 2023. Estimation of steady-state temperature field in multichipmodules using deep convolutional neural network. Thermal Science andEngineering Progress 40 (2023), 101755.
[0230]
[20] Zhen Huang, Hong Wang, Wenkai Yang, Muxi Tang, Depeng Xie, Ting-Jung Lin, Yu Zhang, Wei W. Xing, and Lei He. 2025. Self-Attention to OperatorLearning-based 3D-IC Thermal Simulation. In ACM / IEEE Design AutomationConference (DAC). 1–7.
[0231]
[21] Mingyue Wang, Yuanqing Cheng, Weiheng Zeng, Zhenjie Lu, VasilisF Pavlidis, and Wei Xing. 2025. ARO: Autoregressive Operator Learning forTransferable and Multi-fidelity 3D-IC Thermal Analysis With Active Learning.In International Conference on Computer-Aided Design (ICCAD). 1–9.
[0232]
[22] Dong Lu, Fei Ren, Changlong Li, Kailong Liu, and Naxin Cui.2025. An Interpretable Physics Constrained Deep Operator Network for BatteryInternal Temperature Estimation with Limited Data. IEEE Transactions onTransportation Electrification (2025).
[0233]
[23] Ziyue Liu, Yixing Li, Jing Hu, Xinling Yu, Shinyu Shiau, Xin Ai,Zhiyu Zeng, and Zheng Zhang. 2023. DeepOHeat: Operator Learning-based Ultra-fast Thermal Simulation in 3D-IC Design. In ACM / IEEE Design AutomationConference (DAC).
[0234]
[24] Xinling Yu, Ziyue Liu, Hai Li, Yixing Li, Xin Ai, Zhiyu Zeng,Ian Young, and Zheng Zhang. 2025. DeepOHeat-v1: Efficient Operator Learningfor Fast and Trustworthy Thermal Simulation and Optimization in 3D-IC Design.arXiv preprint arXiv:2504.03955 (2025).
[0235]
[25] Zijian Zhou, Min Tang, and Liang Chen. 2025. ASRR-PINN: AdaptiveSub-Regional Random Resampling-Based PINN for Thermal Analysis of 3D-ICs. InACM / IEEE Design Automation Conference (DAC). 1–7.
[0236]
[26] Liang Chen, Wenxing Zhu, Min Tang, Sheldon X-D Tan, Jun-Fa Mao,and Jianhua Zhang. 2024. PISOV: Physics-Informed Separation of VariablesSolvers for Full-Chip Thermal Analysis. IEEE Transactions on Computer-AidedDesign of Integrated Circuits and Systems (2024).
[0237]
[27] Lu Lu, Pengzhan Jin, Guofei Pang, Zhongqiang Zhang, and GeorgeEm Karniadakis. 2021. Learning nonlinear operators via DeepONet based on theuniversal approximation theorem of operators. Nature machine intelligence 3,3 (2021), 218–229.
[0238]
[28] Wei Huang, Shougata Ghosh, Sivakumar Velusamy, KarthikSankaranarayanan, Kevin Skadron, and Mircea R Stan. 2006. HotSpot: A compactthermal modeling methodology for early-stage VLSI design. IEEE Transactionson very large scale integration (VLSI) systems 14, 5 (2006), 501–513.
[0239]
[29] Scott Ladenheim, Yi-Chung Chen, Milan Mihajlović, and Vasilis FPavlidis. 2018. The MTA: An advanced and versatile thermal simulator forintegrated systems. IEEE Transactions on Computer-Aided Design of IntegratedCircuits and Systems 37, 12 (2018), 3123–3136.
[0240]
[30] COMSOL, Inc. [n. d.]. COMSOL Multiphysics. Burlington, MA, USA.
Claims
1. A method for training an interface constraint expert hybrid model for thermal evaluation of 3D integrated circuits, characterized in that, include: Step 1: Obtain material domain information and power density distribution of 3D integrated circuits and boundary condition mapping The material domain information includes M material subdomains. The stacked structure consists of each subdomain Having piecewise constant thermal conductivity ; Step 2: Construct a neural network model, which adopts an expert hybrid architecture, including a shared global backbone network. and multiple residual experts Each residual expert corresponds to a material subdomain in a three-dimensional integrated circuit that is separated by interfaces with discontinuous thermal conductivity. Step 3: Generate paired sampling points at the interface of adjacent material subdomains using a bilateral sampling strategy. The paired sampling points are located within the adjacent material subdomains and aligned along the interface normal direction. Calculate the bilateral normal temperature gradient based on the paired sampling points. Step 4: Calculate the training loss by weighting and summing the various losses according to their respective weights to form the total loss function. The training loss includes: Volume loss Use the corresponding thermal conductivity Independent calculation of volumetric PDE loss This is used to constrain the temperature field within each material subdomain to satisfy the heat conduction equation; Boundary loss Based on the thermal conductivity of the material at the boundary Calculate the corresponding boundary condition loss to constrain the boundaries of the three-dimensional integrated circuit to satisfy the thermal boundary condition; Interfacial heat flow continuity loss The difference in normal heat flow is calculated based on the bilateral normal temperature gradient and the corresponding material thermal conductivity to constrain the continuity of interfacial heat flow and solve the temperature oscillation problem caused by the discontinuity of heat flow at the material interface. Interface temperature continuity loss The temperature prediction difference between the paired sampling points is used to constrain the continuity of the interface temperature. Step 5: Based on the total loss function Optimize the parameters of the neural network model; Step 6: Input the material domain partitioning information, power density distribution data and boundary condition data of the three-dimensional integrated circuit to be evaluated into the trained neural network model, and output the complete thermal field prediction results.
2. The method according to claim 1, characterized in that, Step 1 includes the following steps: Step 1.1: Obtain material domain information, which includes decomposing the 3D integrated circuit stack into... Material subdomains ,in The total number of material subdomains is a positive integer. Indicates the first Each material subdomain Corresponding to different material regions, and with thermal conductivity fields Relatedly, the thermal conductivity field is a piecewise constant within the material subdomain but discontinuous at the material interface; Step 1.2: Obtain power density distribution The power density distribution To describe the spatial distribution of volumetric heat generation within the chip, used to capture the power consumption of functional modules; Step 1.3: Obtain boundary conditions The boundary conditions Thermal boundary conditions are defined on the domain boundary, including Neumann boundary conditions and Robin boundary conditions.
3. The method according to claim 1, characterized in that, Step 2 Includes the following steps: According to the method of claim 1, the step 2 of constructing the neural network model includes the following steps: Step 2.1: Construct a shared global backbone network The shared global backbone network uses spatial coordinates As input, the output basis function The basis functions are used to capture the overall temperature trend that changes smoothly throughout the entire 3D integrated circuit stack; Step 2.2: For each material subdomain Building a residual expert The residual expert uses spatial coordinates and the thermal conductivity of the material subdomain As input, it is used to model the local thermal behavior correction amount within the corresponding material region, where For the first A constant thermal conductivity value within a material subdomain; Step 2.3: Establish a routing mechanism so that for input spatial coordinates ,according to Material subdomain ,Will Routing to the corresponding residual expert ; Step 2.4: Construct a temperature prediction module by adding the output of the shared global backbone network to the output of the routed residual expert to obtain the final temperature prediction value. ,satisfy: in, This represents the temperature field predicted by the neural network model. Represents the set of model parameters. This indicates that the output of the shared global backbone network is being shared. Indicates the first The output of a residual expert For spatial coordinates, For the first Material subdomains The thermal conductivity (constant within the subdomain). The first is separated by a thermal conductivity discontinuity interface. A material subdomain, , This represents the total number of material subdomains.
4. The method according to claim 2, characterized in that, Step 3, which generates paired sampling points at the interface of adjacent material subdomains and calculates the two-sided normal temperature gradient components using a bilateral sampling strategy, includes the following steps: Step 3.1: For each pair of adjacent material subdomains and Determine its shared interface and the shared interface Perform two-dimensional parameterization, in parameter coordinates The upper layer generates stratified sampling points, where and These are the two-dimensional coordinates of the interface. Step 3.2: Set the parameter coordinates Mapping the layered sampling points on the surface to a three-dimensional Cartesian coordinate system yields the interface sampling points. ; Step 3.3: For each interface sampling point Two offset points are generated along the interface normal as paired sampling points, and the two offset points are respectively located in the material subdomain. and Within, and satisfying the following relationship: in, Indicates that it is located in the material subdomain Offset point within, Indicates that it is located in the material subdomain Offset point within, The preset offset distance, To from the material subdomain Pointing to material subdomain The interface unit normal vector; Step 3.4: The paired sampling points and The data is input into a neural network model to obtain the corresponding predicted temperature field value. and ; Step 3.5: Calculate the neural network model at paired sampling points using automatic differentiation. and Temperature gradient vector at and ; Step 3.6: Calculate the two-sided normal temperature gradient components, i.e., calculate... With interface unit normal vector dot product ,as well as With interface unit normal vector dot product .
5. The method according to claim 3, characterized in that, Step 4, calculating the training loss, includes the following steps: Step 4.1: For each material subdomain ( ), calculate volumetric PDE loss The volumetric PDE loss Calculated using the following formula: in, Represents the gradient operator, For material subdomains thermal conductivity, For the neural network model in spatial coordinates The predicted temperature Let the volume heat source function be... The norm squared of the residual vector; Step 4.2: For each boundary condition (total) (each), calculate boundary condition loss. The boundary condition loss Based on the thermal conductivity of the material at the boundary It can be calculated using the following formula: in, This represents the derivative of the temperature field along the outer normal of the boundary. For the specified heat flux density, The heat transfer coefficient is... For ambient temperature, and It is a boundary type indicator and satisfies ; Step 4.3: Calculate the interfacial heat flow continuity loss The interfacial heat flow continuity loss Calculated using the following formula: in, This represents the total area of all material interfaces. Represents a material subdomain and Shared interface between ( ), and For the paired sampling points generated in step 3, and The two-sided normal temperature gradient represents the temperature gradient components along the interface normal at the paired sampling points. The interface area is a micro-element; Step 4.4: Calculate the interface temperature continuity loss The interface temperature continuity loss Calculated using the following formula: in, and These represent the paired sampling points of the neural network model. and The predicted temperature value; Step 4.5: Sum the losses according to their weights to form the total loss function. : in, , , , These are the weighting coefficients for the corresponding loss terms. The total number of material subdomains, This represents the total number of boundary conditions.
6. The method according to claim 5, characterized in that, Step 5 is based on the total loss function. Optimizing the parameters of the neural network model includes the following steps: Step 5.1: Employ a warm-up training strategy, optimizing only the shared global backbone network during the warm-up phase. The parameters are fixed for all domain-specific residual experts. The parameters enable the shared global backbone network to learn the overall temperature trend of smooth changes in the three-dimensional integrated circuit stack; Step 5.2: After the warm-up phase is completed, remove the fixation on the domain-specific residual expert parameters, and simultaneously optimize the shared global backbone network. The parameters and all domain-specific residual experts The parameters are set to make the neural network model converge to the total loss function. The minimum value.
7. The method according to claim 1, characterized in that, Step 6, which inputs the material domain partitioning information, power density distribution data, and boundary condition data of the three-dimensional integrated circuit to be evaluated into the trained neural network model and outputs the complete thermal field prediction result, includes the following sub-steps: Step 6.1: Obtain the material domain partitioning information of the three-dimensional integrated circuit to be evaluated. The material domain partitioning information includes dividing the three-dimensional integrated circuit into... Material subdomains ,in The total number of material subdomains is a positive integer. Indicates the first Each material subdomain It has thermal conductivity ; Step 6.2: Obtain the power density distribution of the 3D integrated circuit to be evaluated. and boundary conditions ; Step 6.3: Generate a uniformly distributed set of sampling points within the spatial domain of the three-dimensional integrated circuit. ,in The total number of sampling points. For the first The three-dimensional spatial coordinates of each sampling point; Step 6.4: Power density distribution and boundary conditions The input is fed into the branch network of the trained neural network model to generate the encoded latent representation; Step 6.5: For each sampling point The material subdomain to which it belongs is determined based on the material domain partitioning information. , sampling points and local thermal conductivity Input to shared global backbone network Obtaining basis functions and sampling points Thermal conductivity of its material subdomain Input to the corresponding residual expert ; Step 6.6: Perform a single forward pass using the trained neural network model to calculate each sampling point. Temperature prediction at the location Among them, for those belonging to the material subdomain sampling points The predicted temperature value is calculated using the following formula: in, Indicates at the sampling point Predicted temperature at the location, Indicates sharing the global backbone network in The output at that location, Indicates the first A residual expert inputs spatial coordinates. and thermal conductivity Output at time; Step 6.7: Compile all sampling points The temperature prediction value is output, providing a complete thermal field prediction result covering the entire three-dimensional integrated circuit stack.