A calculation method for transient thermal simulation of a chiplet packaging structure

Abstract

The present invention provides a calculation method for transient thermal simulation of a chiplet packaging structure. The method includes: processing the time-domain heat conduction equation based on Laguerre polynomials to establish an equation in the form of order stepping; discretizing the chiplet packaging structure using structured grids in the spatial domain; setting the material properties, heat sources, and boundary conditions of the chiplet packaging structure to obtain a step-by-step equation in matrix form; solving the obtained step-by-step equation order by order through the operator splitting method to sequentially obtain the Laguerre coefficients of each order; and reconstructing the time-domain temperature result of the chiplet packaging structure based on the obtained Laguerre coefficients of each order. The present invention has the advantages of unconditional stability, high calculation efficiency, and low consumption of calculation resources.

Description

Technical Field

[0001] The present invention relates to the technical field of numerical simulation calculation, and specifically, to a calculation method for transient thermal simulation of a chiplet packaging structure. Background Art

[0002] Chiplet-based packaging integration has many advantages, including improving chip yield, reducing manufacturing costs, and shortening the R & D cycle. Different chiplets can be separately manufactured according to the required appropriate process technology, and then assembled through advanced packaging technology, with extremely high flexibility. However, the continuous increase in power density poses challenges to the thermal management of chiplet-based packaging integration. High temperatures may seriously reduce the performance of the system and shorten its lifespan. To solve this problem in the early design stage, accurate and effective transient thermal simulation is very important for evaluating the thermal behavior of chiplet-based systems.

[0003] The numerical solution methods of the transient heat conduction equation can be divided into two categories: time-stepping based and order-stepping based. The order-stepping based method proposed in recent years uses weighted Laguerre polynomials as the basis functions in the time domain and has the advantage of unconditional stability. However, in the process of solving step by step in the existing method, it is necessary to solve a large sparse matrix, facing the problems of long calculation time and high consumption of calculation resources. Summary of the Invention

[0004] Aiming at the defects in the prior art, the purpose of the present invention is to provide a calculation method for transient thermal simulation of a chiplet packaging structure. The method of the present invention has the advantages of unconditional stability, high calculation efficiency, and low consumption of calculation resources.

[0005] To solve the above problems, the technical solution of the present invention is as follows:

[0006] A calculation method for transient thermal simulation of a chiplet packaging structure, comprising the following steps:

[0007] Process the time-domain heat conduction equation based on Laguerre polynomials to establish an equation in the form of order stepping;

[0008] Discretize the chiplet packaging structure using structured grids in the spatial domain;

[0009] Set the material properties, heat sources, and boundary conditions of the chiplet packaging structure to obtain the step-by-step equation in matrix form;

[0010] Solve the step-by-step equation obtained step by step through the operator splitting method to obtain the Laguerre coefficients of each order in turn;

[0011] Reconstruct the time-domain temperature result of the chiplet packaging structure based on the obtained Laguerre coefficients of each order.

[0012] Preferably, in the step of processing the time-domain heat conduction equation based on Laguerre polynomials to establish an equation in a step-by-step form of order, the time-domain heat conduction equation satisfied by the die package structure is:

[0013]

[0014] where r represents the position vector, T is the temperature value, κ is the thermal conductivity of the material, ρ is the density, c is the heat capacity, g represents the time-varying heat source, and the temperature T is expanded in the time domain using Laguerre basis functions as:

[0015]

[0016] where T p is the p-th order Laguerre coefficient, s is the time scale factor, is the weighted Laguerre polynomial, which is defined as: where L p (st) is the p-th order Laguerre polynomial.

[0017] Preferably, in the step of processing the time-domain heat conduction equation based on Laguerre polynomials to establish an equation in a step-by-step form of order, the first derivative of the temperature T with respect to time is expressed in the Laguerre domain as:

[0018]

[0019] Substituting equations and into the heat conduction equation, we get:

[0020]

[0021] Using the orthogonality of the Laguerre basis functions, i.e.:

[0022]

[0023] Using as the test function for Galerkin testing, we get:

[0024]

[0025] where,

[0026]

[0027] where T f is the time when the actual heat source decays to zero.

[0028] Preferably, in the step of discretizing the die package structure using a structured grid in the spatial domain, the three-dimensional die package structure is discretized based on a structured grid of a Cartesian rectangular coordinate system, and the temperature unknowns to be solved are located at the grid nodes. The formula is discretized as:

[0029]

[0030] where T q i,j,k represents the q-th order Laguerre coefficient at the grid point index (i, j, k).

[0031]

[0032] Preferably, in the step of setting the material properties, heat sources, and boundary conditions of the die package structure to obtain the step-by-step equation in matrix form, the material properties of the die package structure include specific heat capacity, density, and thermal conductivity; the heat source is the volumetric heat source density loaded onto the die; the boundary conditions include constant temperature boundary conditions, heat flux boundary conditions, and convective boundary conditions. The discrete equations of the internal nodes and boundary nodes are assembled and combined to obtain the step-by-step equation in matrix form as:

[0033] AT q = g q - β q-1

[0034] where A is the coefficient matrix, g q is the vector related to the heat source and boundary conditions, and β q－1 is the sum of the Laguerre coefficients of the previous q - 1 terms.

[0035] Preferably, in the step of successively obtaining the Laguerre coefficients of each order by solving the step-by-step equation obtained by the operator splitting method, the operator splitting process is divided into the following three steps:

[0036] Step I:

[0037]

[0038] Step II:

[0039]

[0040] Step III:

[0041]

[0042] where,

[0043]

[0044] ΔT q+1 / 3= T q+1 / 3 -T q

[0045] ΔT q+2 / 3 = T q+2 / 3 -T q

[0046] ΔT q+1 = T q+1 -T q

[0047] The matrix forms are respectively:

[0048] Step I:

[0049]

[0050] Step II:

[0051]

[0052] Step III:

[0053]

[0054] where K x , K y and K z are the heat conduction matrices in the x, y, and z directions respectively. The above matrix equations are efficiently solved using the chase method with linear complexity to obtain the Laguerre coefficients of each order.

[0055] Preferably, in the step of reconstructing the time-domain temperature result of the die package structure based on the obtained Laguerre coefficients of each order, the formula for reconstructing the temperature result of the die package structure in the time domain according to the obtained Laguerre coefficients of each order is:

[0056]

[0057] where N is the total number of Laguerre orders in the actual calculation.

[0058] Compared with the prior art, the present invention has the following beneficial effects:

[0059] 1. When solving the transient heat conduction equation, by using the weighted Laguerre polynomial as the time-domain basis function, it is not necessary to consider the size of the time step during the calculation, so it has the advantage of unconditional stability;

[0060] 2. When solving the step-by-step equation of each order in the Laguerre domain, the large sparse matrix is split into three tridiagonal matrices, which can be solved by the chase method with linear complexity, so it has the advantage of high calculation efficiency;

[0061] 3. When solving the step-by-step equation in the Laguerre domain, since only a tridiagonal matrix needs to be stored, the memory requirement is small, and it can be used for the simulation of large-scale complex packaging structures. BRIEF DESCRIPTION OF THE DRAWINGS

[0062] Other features, objects, and advantages of the present invention will become more apparent by reading the following detailed description of non-limiting embodiments with reference to the accompanying drawings:

[0063] Figure 1 It is a flowchart of a calculation method for transient thermal simulation of a chiplet packaging structure;

[0064] Figure 2 It is a schematic diagram of the operator splitting solution process in the Laguerre domain;

[0065] Figure 3 It is an overall geometric model diagram of a specific embodiment;

[0066] Figure 4 It is a side view of the geometric model of a specific embodiment;

[0067] Figure 5 It is a top view of the geometric model of a specific embodiment;

[0068] Figure 6 It is a temperature response curve graph of an observation point in a specific embodiment;

[0069] Figure 7 It is a simulation time graph of a specific embodiment under different numbers of unknowns. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0070] The present invention will be described in detail below with reference to specific embodiments. The following embodiments will help those skilled in the art to further understand the present invention, but do not limit the present invention in any form. It should be noted that those of ordinary skill in the art can make several changes and improvements without departing from the concept of the present invention. These all belong to the protection scope of the present invention.

[0071] Specifically, the present invention provides a calculation method for transient thermal simulation of a chiplet packaging structure, as Figure 1 shown, the method includes the following steps:

[0072] S1: Process the time-domain heat conduction equation based on the Laguerre polynomial to establish an equation in the form of step-by-step order;

[0073] Specifically, in the step S1, the time-domain heat conduction equation satisfied by the chiplet packaging structure is:

[0074]

[0075] Among them, r represents the position vector, T is the temperature value, κ is the thermal conductivity of the material, ρ is the density, c is the heat capacity, and g represents the time-varying heat source.

[0076] The temperature T is expanded in the time domain using Laguerre basis functions as:

[0077]

[0078] Among them, T p is the Laguerre coefficient of order p, s is the time scale factor, is the weighted Laguerre polynomial, which is defined as:

[0079]

[0080] Among them, L p (st) is the Laguerre polynomial of order p.

[0081] The first derivative of the temperature T with respect to time is expressed in the Laguerre domain as:

[0082]

[0083] Substituting Eqs. (2) and (4) into Eq. (1), we get:

[0084]

[0085] Using the orthogonality of the Laguerre basis functions, that is:

[0086]

[0087] Using as the test function to perform the Galerkin test on (5), we get:

[0088]

[0089] Among them,

[0090]

[0091] In the above equation, T f is the time when the actual heat source decays to zero.

[0092] S2: Discretize the die package structure using a structured grid in the spatial domain;

[0093] Specifically, in the step S2, the three-dimensional die package structure is discretized based on the structured grid of the Cartesian rectangular coordinate system, and the temperature unknowns to be solved are located at the grid nodes. Equation (7) is discretized as:

[0094]

[0095] Among them, T q i,j,k represents the q-th order Laguerre coefficient at the grid point index (i, j, k).

[0096]

[0097] S3: Set the material properties, heat source, and boundary conditions of the die package structure to obtain the step-by-step equation in matrix form;

[0098] Specifically, the material properties of the die package structure include specific heat capacity, density, and thermal conductivity; the heat source is the volumetric heat source density loaded on the die; the boundary conditions include constant temperature boundary conditions, heat flux boundary conditions, and convective boundary conditions.

[0099] Assemble and combine the discrete equations of the internal nodes and boundary nodes to obtain the step-by-step equation in matrix form as:

[0100] AT q = g q -β q-1 (11)

[0101] Among them, A is the coefficient matrix, g q is the vector related to the heat source and boundary conditions, and β q－1 is the sum of the Laguerre coefficients of the previous q - 1 terms.

[0102] S4: Solve the step-by-step equation obtained by the operator splitting method order by order to obtain the Laguerre coefficients of each order;

[0103] Specifically, as Figure 2 shown, the operator splitting process is divided into the following three steps:

[0104] Step I:

[0105]

[0106] Step II:

[0107]

[0108] Step III:

[0109]

[0110] Among them,

[0111]

[0112] The matrix form of Eqs. (12)-(14) is:

[0113] Step I:

[0114]

[0115] Step II:

[0116]

[0117] Step III:

[0118]

[0119] where K x , K y and K z are the heat conduction matrices in the x, y, and z directions respectively. The above matrix equation is efficiently solved using the chase method with linear complexity to obtain the Laguerre coefficients of each order.

[0120] S5: Based on the obtained Laguerre coefficients of each order, reconstruct the time-domain temperature result of the die package structure.

[0121] Specifically, the formula for reconstructing the time-domain temperature result of the die package structure according to the Laguerre coefficients of each order obtained in step S4 is:

[0122]

[0123] where N is the total order of Laguerre in the actual calculation.

[0124] According to the above simulation method, a specific embodiment is calculated.

[0125] In this embodiment, the geometric model of the die package structure is as Figure 3 shown. Figure 4 is the stacked structure of the specific embodiment, Figure 5 is the dimension diagram of the specific embodiment.

[0126] The heat sources of the three dies are set as:

[0127]

[0128] where,

[0129] The top surface of the die is a constant temperature boundary of 300K, representing the function of the heat sink, and all other outer surfaces are set as convective boundaries with a convective coefficient of 5 W / (m2·K).

[0130] Two observation points (P1, P2) are set at the centers of die 1 and die 3 respectively. Thermal simulations are carried out by the traditional Laguerre method and the method of the present invention, and the order is 300 for both. Their temperature responses are compared in Figure 6 , and the mean absolute errors are 0.077K and 0.028K respectively, which fully demonstrates the accuracy of the method of the present invention.

[0131] Next, for cases with different numbers of unknowns, calculations were respectively carried out on a server with Intel Xeon Gold 5218 CPUs (2.30 GHz) and 256 GB of memory using the traditional Laguerre method and the method of the present invention. The simulation time consumption is as Figure 7 shown. It can be seen that compared with the traditional Laguerre method, the method of the present invention has an order of magnitude improvement in efficiency, and as the number of unknowns increases, the required calculation time increases linearly, fully demonstrating the efficiency advantage of the method of the present invention.

[0132] The specific embodiments of the present invention have been described above. It should be understood that the present invention is not limited to the above specific embodiments, and those skilled in the art can make various changes or modifications within the scope of the claims, which does not affect the essence of the present invention. Without conflict, the embodiments of the present application and the features in the embodiments can be combined with each other arbitrarily.

Claims

1. A calculation method for transient thermal simulation of a chiplet packaging structure, characterized in that, the method comprises the following steps: Processing the time-domain heat conduction equation based on Laguerre polynomials to establish an equation in the form of order stepping; Discretizing the chiplet packaging structure using structured grids in the spatial domain; Setting the material properties, heat sources, and boundary conditions of the chiplet packaging structure to obtain a step-by-step equation in matrix form. The material properties of the chiplet packaging structure include specific heat capacity, density, and thermal conductivity; the heat source is the volumetric heat source density loaded onto the chiplet; the boundary conditions include constant temperature boundary conditions, heat flux boundary conditions, and convective boundary conditions. Assembling and combining the discrete equations of internal nodes and boundary nodes, the step-by-step equation in matrix form is: AT q = g q - β q-1 where A is the coefficient matrix, g q is a vector related to the heat source and boundary conditions, and β q－1 is the sum of the Laguerre coefficients of the previous q - 1 terms; Solving the step-by-step equation obtained by the operator splitting method order by order to obtain the Laguerre coefficients of each order in sequence; Based on the obtained Laguerre coefficients of each order, reconstructing the time-domain temperature result of the chiplet packaging structure.

2. The calculation method for transient thermal simulation of a chiplet packaging structure according to claim 1, characterized in that, in the step of processing the time-domain heat conduction equation based on Laguerre polynomials to establish an equation in the form of order stepping, the time-domain heat conduction equation satisfied by the chiplet packaging structure is: where r represents the position vector, T is the temperature value, κ is the thermal conductivity of the material, ρ is the density, c is the heat capacity, g represents the time-varying heat source, and the temperature T is expanded in the time domain using Laguerre basis functions as: Among them, T p is the Laguerre coefficient of order p, and s is the time scale factor. is the weighted Laguerre polynomial, which is defined as: Among them, L p (st) is the Laguerre polynomial of order p.

3. The calculation method for transient thermal simulation of a chiplet packaging structure according to claim 2, characterized in that, in the step of processing the time-domain heat conduction equation based on Laguerre polynomials to establish an equation in the form of order stepping, the first derivative of the temperature T with respect to time is expressed in the Laguerre domain as: Substitute the expressions and into the heat conduction equation, and we get: Utilizing the orthogonality of Laguerre basis functions, that is: Use Perform the Galerkin test using it as the test function, and obtain: where, Among them, T f is the time when the actual heat source decays to zero.

4. The calculation method for transient thermal simulation of a chiplet packaging structure according to claim 3, characterized in that, In the step of discretizing the die packaging structure using a structured grid in the spatial domain, the three-dimensional die packaging structure is discretized based on a structured grid of a Cartesian rectangular coordinate system, and the temperature unknowns to be solved are located at the grid nodes. Equation is discretized as: where T q i,j,k represents the q-th order Laguerre coefficient at the grid point index (i, j, k).

5. The calculation method for transient thermal simulation of a chiplet packaging structure according to claim 1, characterized in that, in the step of solving the step-by-step equation obtained by the operator splitting method order by order to obtain the Laguerre coefficients of each order in sequence, the operator splitting process is divided into the following three steps: Step I: Step II: Step III: where, ΔT q+1 / 3 = T q+1 / 3 - T q ΔT q+2 / 3 = T q+2 / 3 - T q ΔT q+1 = T q+1 - T q The matrix forms are respectively: Step I: Step II: Step III: where K x , K y and K z are the heat conduction matrices in the x, y, and z directions respectively. The above matrix equation is efficiently solved by using the chase method with linear complexity to obtain the Laguerre coefficients of each order.

6. The calculation method for transient thermal simulation of a chiplet packaging structure according to claim 5, characterized in that, in the step of reconstructing the time-domain temperature result of the chiplet packaging structure based on the obtained Laguerre coefficients of each order, the formula for reconstructing the temperature result of the chiplet packaging structure in the time domain according to the obtained Laguerre coefficients of each order is: where N is the total order of Laguerre in actual calculation.

Images