A thermal-mechanical combined simulation data processing method and system for ceramic packaging

Abstract

The application provides a kind of ceramic package-oriented heat-force combined simulation data processing method and system, belong to packaging simulation analysis and optimization technical field, the application can significantly improve the calculation efficiency and solution accuracy of ceramic package heat-force combined simulation, by optimizing the layout direction of heat conduction channel, effectively solve the heat dissipation bottleneck and stress concentration problem caused by the mismatch of heat flow direction and channel trend in traditional design, finally the heat conduction channel structure not only realizes the shortest heat dissipation path, effectively suppresses the chip hot spot temperature, at the same time, the maximum thermal stress is strictly controlled below the material yield strength, so as to significantly shorten the design iteration cycle under the premise of ensuring the reliability of packaging, improve the design efficiency of multi-physical field collaborative optimization.

Description

Technical Field

[0001] This invention belongs to the field of packaging simulation analysis and optimization technology, specifically relating to a thermal-mechanical co-simulation data processing method and system for ceramic packaging. Background Technology

[0002] Ceramic packaging, due to its excellent thermal stability, hermeticity, and insulation properties, has been widely used in aerospace, military electronics, high-power LEDs, and RF microwave devices. As electronic devices continue to evolve towards higher integration and higher power density, the heat flux density inside ceramic packages is increasing dramatically, leading to increasingly significant thermo-mechanical coupling effects. Localized overheating not only causes material performance degradation but also induces thermal stress concentration due to thermal expansion coefficient mismatch, potentially resulting in package cracking or solder joint fatigue failure. Therefore, accurately predicting the temperature distribution and thermal stress state of ceramic packages under thermo-mechanical synergy has become a critical aspect of reliability design. Currently, multiphysics simulations of ceramic packages typically employ the finite element method (FEM), which involves constructing a three-dimensional geometric model, applying thermal loads and mechanical boundary conditions, and solving the heat conduction and thermoelasticity equations. To further improve heat dissipation performance while maintaining structural strength, topology optimization techniques have been introduced into the packaging design field. Using material density as a design variable, this method automatically generates the optimal heat conduction channel layout under given volume and stress constraints.

[0003] However, existing thermal-mechanical co-simulation and optimization methods for ceramic packaging still have many shortcomings. First, in the mesh generation stage, traditional methods typically employ a uniform mesh size control strategy for the entire model. Ceramic packaging contains multiple materials; while a globally fine mesh can ensure the solution accuracy in high thermal conductivity regions, it leads to excessively dense meshes in low thermal conductivity regions, causing a sharp increase in computational scale. Conversely, a globally coarse mesh cannot accurately capture the temperature gradient in high thermal conductivity regions, resulting in distorted simulation results. Second, in the topology optimization process, existing methods often use heat dissipation weakness or minimum temperature as the sole optimization objective, neglecting the synergistic relationship between the direction of heat conduction channels and the direction of heat flow. In traditional optimization frameworks, stress constraints are usually used as a post-hoc verification index; if the optimization result exceeds the stress limit, parameters need to be manually adjusted and recalculated, resulting in low design efficiency and difficulty in guaranteeing convergence. Summary of the Invention

[0004] The purpose of this section is to outline some aspects of the embodiments of the present invention and to briefly introduce some preferred embodiments. Some simplifications or omissions may be made in this section, as well as in the abstract and title of the present application, to avoid obscuring the purpose of this section, the abstract and title of the invention. Such simplifications or omissions shall not be used to limit the scope of the present invention.

[0005] In view of the aforementioned existing problems, the present invention is proposed.

[0006] Therefore, the technical problem solved by this invention is the stress concentration caused by inefficient mesh generation and the mismatch between the direction of heat conduction channels and the direction of heat flow in existing ceramic packaging thermo-mechanical co-simulation.

[0007] To address the aforementioned technical problems, the present invention provides the following technical solution:

[0008] A thermal-mechanical co-simulation data processing method for ceramic packaging includes: acquiring a three-dimensional model and material parameters of the ceramic package, wherein the material parameters include the thermal conductivity of each material, and dividing the three-dimensional model into high heat flux data slices and low heat flux data slices according to the thermal conductivity threshold; generating first grid data by fine meshing of the high heat flux data slices, generating second grid data by coarse meshing of the low heat flux data slices, and combining the first grid data and the second grid data to form a joint mesh model; performing thermal simulation based on the joint mesh model to obtain temperature field distribution data, identifying the heat flow direction based on the temperature field distribution data, and chaining the heat flow direction with the preset heat conduction channel direction to generate a direction correction factor; and iteratively adjusting the layout direction of the heat conduction channel based on the direction correction factor until optimized channel structure data that satisfies the shortest heat dissipation path and stress value is lower than the stress limit threshold is obtained.

[0009] In a preferred embodiment of the present invention, the following steps are taken: reading the material identifier of each entity unit in the three-dimensional model of the ceramic package; retrieving the corresponding thermal conductivity data from the material parameter library based on the material identifier; assigning the thermal conductivity data to the corresponding entity unit; calculating the median of the thermal conductivity data of all entity units; setting the median as the thermal conductivity threshold; and comparing the thermal conductivity data of each entity unit with the thermal conductivity threshold; marking entity units with thermal conductivity data greater than the thermal conductivity threshold as high heat flux units; and marking entity units with thermal conductivity data less than or equal to the thermal conductivity threshold as low heat flux units.

[0010] In a preferred embodiment of the present invention, the division of the high heat flux data slice and the low heat flux data slice includes: grouping the same marked entity units sharing a face or edge into continuous regions through neighborhood search to form a high heat flux unit group and a low heat flux unit group; for each high heat flux unit group and each low heat flux unit group, traversing the faces of all entity units within the group, extracting the face owned by only one unit as the boundary triangle facet of the corresponding group, recording the three vertex numbers and vertex coordinates of the boundary triangle facet to form the boundary surface data of the corresponding group; and simultaneously recording the index list of all entity units within the group to form data slice identification data.

[0011] In a preferred embodiment of the present invention, the following steps are taken: The data slice identification data and the boundary surface data are read to identify the shared boundary triangle facets between different data slices, and the regions where the boundary triangle facets are located are marked as conformal boundaries; simultaneously, the spatial range of each data slice is determined based on the boundary surface data; the unit index list in the data slice identification data is read to locate the original solid units of each data slice in the initial 3D model, and the geometric dimension data of the original solid units is extracted: for high heat flux data slices, the smallest side length of the contained original solid units is taken. One-tenth is used as the target side length for the fine mesh; for low heat flux data patches, the average feature size of the original solid cells is taken. One-tenth of the target side length is used as the coarse-grained mesh target side length; this target side length is only used for mesh generation in this preprocessing stage and remains unchanged in subsequent optimization iterations; for each conformal boundary, the shared boundary triangle facet is used as a geometric reference, and the shared boundary triangle facet is re-divided using the fine mesh target side length to generate a conformal surface mesh; the nodes on the conformal surface mesh completely overlap between adjacent data slices; for each high heat flux data slice and each low heat flux data slice, the boundary triangle facet on the non-conformal boundary is used as a geometric reference, and the target side length corresponding to the data slice is used to re-divide the data slice surface mesh to generate a data slice surface mesh; the data slice surface mesh belongs to only a single data slice.

[0012] As a preferred embodiment of the present invention, the formation of the joint mesh model includes: for each high heat flux data slice, splicing the conformal surface mesh portion and the data slice surface mesh portion to form a complete boundary mesh for the corresponding data slice; using the complete boundary mesh as the boundary and the fine mesh target edge length as the internal mesh size control, generating a tetrahedral mesh to obtain the first mesh data; for each low heat flux data slice, splicing the conformal surface mesh portion and the data slice surface mesh portion to form a complete boundary mesh for the corresponding data slice; using the complete boundary mesh as the boundary and the coarse-grained mesh target edge length as the internal mesh size control, generating a tetrahedral mesh to obtain the second mesh data; integrating the first mesh data and the second mesh data into the same mesh file, recording the node correspondence on the conformal boundary, and forming a joint mesh model.

[0013] As a preferred embodiment of the present invention, the identification of heat flow direction includes: inputting the joint mesh model and the node temperature data obtained from thermal simulation, and inputting a preset heat conduction channel direction vector; the heat conduction channel direction vector is a global constant or locally specified according to the spatial location; for each mesh cell in the joint mesh model, the temperature gradient at the center of the cell is calculated using the node temperature data and the derivative of the cell shape function, and the opposite direction of the temperature gradient is used as the cell heat flow direction vector.

[0014] As a preferred embodiment of the present invention, the step of performing chained data linkage to generate direction correction factors includes: for each node, collecting all mesh cells containing the node, weighting the heat flow direction vector of the cell by the cell volume to obtain the local heat flow direction vector of the node, and performing normalization processing; for each node, calculating the cosine value of the angle between the normalized local heat flow direction vector of the node and the heat conduction channel direction vector, using the cosine value of the angle as the direction deviation coefficient, and then subtracting the direction deviation coefficient from 1 as the node direction correction factor; and outputting an array of direction correction factors for all nodes.

[0015] As a preferred embodiment of the present invention, the iterative adjustment of the layout direction of the heat conduction channel includes: assigning a density design variable to each element in the joint mesh model; applying convective heat transfer boundary conditions and mechanical boundary conditions to the joint mesh model, performing thermo-mechanical coupled finite element solution to obtain nodal temperature field data and nodal stress field data; using the nodal temperature field data of the current iteration as input, calculating the nodal orientation correction factor array of the current iteration; interpolating the nodal orientation correction factor array to the mesh element, and for each mesh element, taking the average value of the orientation correction factors of the contained nodes to obtain the element orientation correction factor; projecting the density design variable to the nodes, calculating the nodal density gradient, and then interpolating it back to the element center to obtain the element density gradient vector, and normalizing it to obtain the element density gradient direction vector.

[0016] As a preferred embodiment of the present invention, the iterative adjustment of the layout direction of the heat conduction channel further includes: constructing a comprehensive objective function containing a shortest heat dissipation path term, a direction penalty term, and a stress constraint penalty term; calculating the sensitivity of the comprehensive objective function to the density design variable using the adjoint method; updating the density design variable using the optimization criterion method or the moving asymptote method based on the sensitivity to obtain a new density design variable; determining whether the maximum change between the new density design variable and the density design variable in the previous iteration step is less than the optimization convergence tolerance; if so, terminating the iteration, and determining that the mesh element with a density value greater than a preset density threshold is a material element, and outputting the heat conduction channel structure data composed of the material element; if not, returning to the thermal-mechanical coupled finite element solution and continuing the iterative adjustment.

[0017] On the other hand, the present invention also provides a thermal-mechanical co-simulation data processing system for ceramic packaging, comprising: a data slice partitioning module, which acquires a three-dimensional model and material parameters of the ceramic package, wherein the material parameters include the thermal conductivity of each material, and partitions the three-dimensional model into high heat flux data slices and low heat flux data slices according to the thermal conductivity threshold; a joint mesh generation module, which generates first mesh data by fine mesh partitioning of the high heat flux data slices and generates second mesh data by coarse mesh partitioning of the low heat flux data slices, and combines the first mesh data and the second mesh data to form a joint mesh model; a correction factor calculation module, which performs thermal simulation based on the joint mesh model to obtain temperature field distribution data, identifies the heat flow direction based on the temperature field distribution data, and performs chain-like data linkage between the heat flow direction and a preset heat conduction channel direction to generate a direction correction factor; and a channel structure optimization module, which iteratively adjusts the layout direction of the heat conduction channel based on the direction correction factor until optimized channel structure data that satisfies the shortest heat dissipation path and the stress value is lower than the stress limit threshold is obtained.

[0018] The beneficial effects of this invention are as follows: Compared with the prior art, the technical effects of this invention are as follows: This invention can significantly improve the computational efficiency and solution accuracy of thermo-mechanical co-simulation of ceramic packaging. By optimizing the layout direction of the heat conduction channel, it effectively solves the heat dissipation bottleneck and stress concentration problems caused by the mismatch between the heat flow direction and the channel direction in traditional design. The final heat conduction channel structure not only achieves the shortest heat dissipation path and effectively suppresses the chip hot spot temperature, but also strictly controls the maximum thermal stress below the material yield strength. Thus, while ensuring the reliability of packaging, it significantly shortens the design iteration cycle and improves the design efficiency of multi-physics field co-optimization. Attached Figure Description

[0019] Figure 1 This is a flowchart of a thermal-mechanical co-simulation data processing method for ceramic packaging as described in this invention.

[0020] Figure 2 This is a cloud map showing the chip temperature distribution in an embodiment of the present invention.

[0021] Figure 3 This is a cloud map showing the overall temperature distribution of the package in an embodiment of the present invention.

[0022] Figure 4 This is a side temperature distribution cloud map in an embodiment of the present invention.

[0023] Figure 5 This is a schematic diagram of the simulation environment in an embodiment of the present invention.

[0024] Figure 6This is a structural diagram of a thermal-mechanical co-simulation data processing system for ceramic packaging as described in this invention. Detailed Implementation

[0025] To make the objectives, technical solutions, and advantages of this invention clearer, the technical solutions of this invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments of this invention. The embodiments described in this application are merely some embodiments of this invention, and not all embodiments. Other embodiments obtained by those skilled in the art based on the spirit of this invention without creative effort are all within the protection scope of this invention.

[0026] like Figure 1 As shown, the thermal-mechanical co-simulation data processing method for ceramic packaging according to the present invention includes:

[0027] S1: Obtain the three-dimensional model and material parameters of the ceramic package. The material parameters include the thermal conductivity of each material. Based on the thermal conductivity threshold, the three-dimensional model is divided into high heat flux data slices and low heat flux data slices.

[0028] S1.1: Read the material identifiers of each solid unit in the 3D model of the ceramic package, retrieve the corresponding thermal conductivity data from the material parameter library based on the material identifiers, and assign the thermal conductivity data to the corresponding solid unit. The 3D model of the ceramic package includes a multi-layered structure such as a chip and a ceramic casing.

[0029] S1.2: Calculate the median of the thermal conductivity data of all solid elements, set the median as the thermal conductivity threshold, and compare the thermal conductivity data of each solid element with the thermal conductivity threshold: mark solid elements with thermal conductivity data greater than the thermal conductivity threshold as high heat flux elements, and mark solid elements with thermal conductivity data less than or equal to the thermal conductivity threshold as low heat flux elements.

[0030] If the total number of units is odd, the median is the value that is in the middle after sorting; if the total number of units is even, the median is the arithmetic mean of the two middle values ​​after sorting.

[0031] It should be noted that the median is not sensitive to extreme values. When there are a very small number of materials with extremely high or low thermal conductivity in the model, the median can stably reflect the overall distribution center of the material's thermal conductivity, avoiding the threshold from deviating from the actual physical meaning due to individual extreme values, thereby ensuring that the division ratio of high / low thermal flux units is balanced and reasonable.

[0032] S1.3: Through neighborhood search, the same labeled entity units sharing a face or edge are grouped into continuous regions to form high heat flux unit groups and low heat flux unit groups. For each high heat flux unit group and each low heat flux unit group, the faces of all entity units in the group are traversed, and the faces owned by only one unit are extracted as the boundary triangle facets of the corresponding group. The three vertex numbers and vertex coordinates of the boundary triangle facets are recorded to form the boundary surface data of the corresponding group. At the same time, the index list of all entity units in the group is recorded to form the data patch identification data.

[0033] It should be noted that after unit labeling, the entity units in the model are divided into two categories of labels, but units with the same label may appear as multiple isolated regions that are not connected to each other in space. To identify physically independent heat conduction channel regions, spatial connectivity merging is required. In this embodiment of the invention, the processor uses a neighborhood search algorithm. First, an empty queue is established, and a high heat flux unit that has not yet been visited is randomly selected as a seed unit, pushed into the queue, and marked as visited. When the queue is not empty, a unit is popped from the head of the queue, and all adjacent units that share a face or edge with this unit are retrieved. It is determined whether the adjacent units also carry the high heat flux label and have not been visited. If so, the adjacent unit is pushed into the tail of the queue and marked as visited.

[0034] Repeat the above process until the queue is empty. At this point, all visited cells form a high heat flux cell group. Then, continue searching for the next unvisited high heat flux cell as a new seed, repeating the search process until all high heat flux cells are assigned to a cell group. Use the same method to process low heat flux cells, obtaining several low heat flux cell groups.

[0035] Furthermore, for each unit group obtained through the above steps, the processor traverses all entity units within the unit group, accessing the constituent faces of each unit. For each face, the processor obtains the numbers of the three vertices constituting the face, and uses the combination of these three vertex numbers as the unique identifier of the face. The processor maintains a face counter. When a face combination is encountered for the first time, the face is recorded as a face of the current unit group, and the count is set to 1; when the same face combination is encountered again, the count of the face is incremented by 1. After traversal, all faces with a count of 1 are faces owned by only one unit in the unit group, i.e., the boundary triangle facets of the unit group. For each such boundary triangle facet, the processor records the global numbers of the three vertices and the coordinate values ​​of each vertex in three-dimensional space, forming the boundary surface data of the unit group. At the same time, the processor records the index list of all entity units within the unit group, forming the data patch identifier data.

[0036] It should be noted that this invention integrates spatially dispersed but materially similar units into physically continuous independent regions through connectivity merging, providing a clear processing object for subsequent differentiated meshing of each independent region; and by extracting boundary triangular patches, the discrete unit mesh information is transformed into a continuous boundary surface representation, providing an accurate geometric reference for generating a consistent surface mesh on the conformal boundary in subsequent steps, avoiding mesh penetration or gap problems caused by unit boundary mismatch.

[0037] S2: The high heat flux data slice is divided into a fine mesh to generate the first mesh data, and the low heat flux data slice is divided into a coarse mesh to generate the second mesh data. The first mesh data and the second mesh data are then combined to form a joint mesh model.

[0038] S2.1: Read the data slice identifier data and boundary surface data, identify the boundary triangle facets shared between different data slices, and mark the area where the boundary triangle facets are located as conformal boundaries; at the same time, determine the spatial range of each data slice based on the boundary surface data.

[0039] In practice, the data slice identifier data contains a list of entity unit indexes corresponding to each data slice, and the boundary surface data contains information about the triangle faces on the boundary of each data slice. Each triangle face records the global ID and spatial coordinates of its three vertices. The processor traverses the boundary surface data of all data slices and constructs a hash mapping table with the combination of triangle face vertex combinations as keys and the data slice ID as values.

[0040] For each boundary triangle facet, the global numbers of its three vertices are taken and sorted in ascending order to form a unique triangle identifier. This identifier is used as the key, and the current data slice number is used as the value to insert it into a hash map table. If, during insertion, the triangle identifier is found to already exist, and the corresponding data slice number is different from the current data slice number, it indicates that the triangle facet belongs to two different data slices simultaneously, meaning that the triangle facet is a shared boundary triangle facet.

[0041] The processor marks the region containing such shared boundary triangle faces as a conformal boundary and records the numbers of the two data patches involved in the conformal boundary, as well as the vertex information of the shared triangle faces. Simultaneously, based on the vertex coordinate range of all triangle faces in the boundary surface data of each data patch, the processor determines the minimum and maximum coordinate values ​​in the X, Y, and Z directions for each data patch, which serve as the spatial bounding box of that data patch.

[0042] S2.2: Read the cell index list in the data slice identifier data, locate the original solid cell of each data slice in the initial 3D model, calculate the Euclidean distance of all side lengths for each original solid cell, and take the minimum value as the minimum side length of the cell; at the same time, calculate the cube root of the volume as the average feature size of the cell.

[0043] For high heat flux data slices, the processor collects the minimum side length of all original physical units within the data slice. One-tenth of the target side length for the fine mesh is used in this embodiment of the invention. Set to 5; for low heat flux data slices, the processor collects the average feature size of all raw solid units within the data slice, calculates the arithmetic mean of these average feature sizes, and then... In this embodiment of the invention, one-tenth is used as the target side length of the coarse-grained grid. Set to 2; the target side length is only used for this preprocessing mesh generation, and the mesh remains unchanged in subsequent optimization iterations.

[0044] It is important to emphasize that the target side length is only used for this preprocessing mesh generation. In the subsequent topology optimization iteration process, the mesh will remain fixed and will not be regenerated.

[0045] It can be seen that using one-fifth of the minimum side length as the target side length for high heat flux data patches ensures the generation of a sufficiently fine mesh in high thermal conductivity regions to accurately capture heat flux density distribution and temperature gradient changes, such as... Figure 2 The chip temperature distribution cloud map shown reveals a clear temperature gradient on the chip surface. A fine mesh helps to accurately characterize this gradient. Furthermore, using half the average feature size as the target side length for low heat flux data chips minimizes the number of meshes while maintaining basic solution accuracy. Figure 3 In the overall temperature distribution cloud diagram of the package shown, the temperature change of the ceramic tube shell is relatively gradual, and the coarse-grained grid is sufficient to meet the solution requirements.

[0046] S2.3: For each conformal boundary, using the shared boundary triangle facet as a geometric reference, the shared boundary triangle facet is re-divided using the fine mesh target side length to generate a conformal surface mesh; the nodes on the conformal surface mesh completely overlap between adjacent data patches.

[0047] The process of re-division includes:

[0048] The processor first preserves the original nodes at the vertices of the shared boundary triangular facets, then inserts new nodes along the edges of the shared boundary according to the target edge length of the fine mesh, and finally generates new triangular elements within the region enclosed by these nodes. Since the target edge length of the fine mesh is consistent for the data patches on both sides of the conformal boundary, and the re-meshing process is performed with the same geometric reference and the same dimensional parameters, the generated conformal surface mesh has completely overlapping nodes and element edges between adjacent data patches.

[0049] The processor records the node coordinates, node numbers, and cell connection relationships of the newly generated conformal surface mesh, and establishes a mapping relationship between the new nodes and the vertices of the original shared boundary triangle facets.

[0050] It should be noted that this invention fundamentally eliminates the problem of mismatched interface nodes caused by independent meshing in traditional methods by generating a consistent fine surface mesh on the conformal boundary.

[0051] S2.4: For each high heat flux data patch and each low heat flux data patch, the boundary triangle facets on the non-conformal boundary are used as geometric references, and the data patch surface mesh is generated by re-dividing the data patch using the target side length corresponding to the data patch; the data patch surface mesh belongs to only a single data patch.

[0052] Non-conformal boundaries include the outer surface of the data chip that contacts the external environment, and the internal boundaries of the data chip that are not adjacent to other data chips.

[0053] Specifically, for each high heat flux data slice, the processor extracts all boundary triangle faces of the data slice from the boundary surface data generated in step S1.3, and removes the shared boundary triangle faces that have been processed in step S2.3. The remaining part is the non-conformal boundary triangle faces of the data slice.

[0054] Using these non-conformal boundary triangular patches as geometric references, and employing the fine-grained target side length of the high-heat-flux data patch as the size control parameter, the non-conformal boundaries are re-meshed using the same mesh generation algorithm as in step S2.3 to generate the surface mesh of the data patch. For each low-heat-flux data patch, the same processing method is used, but the mesh size control parameter adopts the coarse-grained target side length of the low-heat-flux data patch.

[0055] Since these surface meshes serve only a single data patch, the nodes and cells do not need to match other data patches, so they can be generated independently.

[0056] The processor records the node coordinates, node numbers, and cell connection relationships of the surface mesh for each data slice, and establishes a boundary correspondence with the conformal surface mesh generated in step S2.3.

[0057] It should be noted that decoupling the meshing of non-conformal boundaries from that of conformal boundaries allows each data patch to independently control the mesh density of its internal region based on its own physical characteristics, while ensuring interface matching, thus achieving a balance between local fineness and global efficiency.

[0058] S2.5: For each high heat flux data slice, the conformal surface mesh portion and the data slice surface mesh portion are spliced ​​together to form the complete boundary mesh of the corresponding data slice. Using the complete boundary mesh as the boundary, the fine mesh target side length is used as the internal mesh size control to generate a tetrahedral mesh and obtain the first mesh data.

[0059] During the stitching process, since the nodes of the conformal surface mesh have already been matched with the adjacent data slices, while the nodes of the non-conformal surface mesh only belong to the current data slice, no additional processing is required during the stitching process. The nodes and cells of the two meshes can be merged to form the complete boundary mesh of the data slice.

[0060] The processor uses a complete boundary mesh as the constraint boundary and a fine mesh target edge length as the internal mesh size control parameter. It applies a 3D Delaunay triangulation algorithm or a leading-edge method to generate tetrahedral elements within the spatial region enclosed by the boundary mesh. The generated tetrahedral mesh satisfies the node position constraints on the boundary mesh, has uniform internal node distribution, and its element quality meets the requirements for finite element solution.

[0061] The processor stores the generated node coordinates, node numbers, cell connection relationships, and the data slice identifier to which the cell belongs, forming the first grid data. For example... Figure 2 and Figure 3 The temperature distribution cloud map shown is supported by a mesh model that needs to be obtained through this kind of fine mesh generation method to ensure accurate capture of the temperature gradient in the chip area.

[0062] S2.6: For each low heat flux data slice, the conformal surface mesh portion and the data slice surface mesh portion are spliced ​​together to form the complete boundary mesh of the corresponding data slice. Using the complete boundary mesh as the boundary, the target side length of the coarse-grained mesh is used as the internal mesh size control to generate a tetrahedral mesh and obtain the second mesh data.

[0063] For each low heat flux data slice, the same processing procedure as step S2.5 is used.

[0064] Because the target side length of the low heat flux data patch is relatively large, the number of generated cells is significantly less than that of the high heat flux data patch of the same volume, thereby greatly reducing the overall mesh size while ensuring the solution accuracy.

[0065] The processor stores the generated node coordinates, node numbers, cell connection relationships, and data slice identifiers to which the cells belong, forming the second grid data.

[0066] like Figure 4 As shown in the side temperature distribution cloud map, the temperature change of the ceramic tube shell is gradual, and a coarse-grained grid is sufficient to meet the solution requirements. The second grid data generated by this invention can reduce the computational cost while ensuring the computational accuracy.

[0067] S2.7: Integrate the first grid data and the second grid data into the same grid file, record the node correspondence on the conformal boundary, and form a joint grid model.

[0068] During integration, the processor globally renumbers the node numbers independently generated for each data slice to avoid node number conflicts; simultaneously, it updates the node connection relationships of each unit to ensure that the references are globally unique node numbers. For nodes on conformal boundaries, since node overlap was already ensured in step S2.3, these nodes retain the same global number between adjacent data slices after global renumbering.

[0069] The processor writes the integrated array of node coordinates, element connectivity, and element material properties into the same mesh file. The file format can be a common finite element mesh format such as INP, MSH, or CGNS. Simultaneously, the processor maintains a separate table of conformal boundary node correspondences. Each row in this table records the global ID of a conformal boundary node and its associated data patch ID, facilitating rapid identification of cross-boundary nodes in subsequent thermo-mechanical co-simulations.

[0070] The resulting joint mesh model has the following characteristics: fine mesh in high heat flux regions, coarse mesh in low heat flux regions, and conformal continuous mesh at data slice interfaces.

[0071] S3: Perform thermal simulation based on the joint mesh model to obtain temperature field distribution data, identify the heat flow direction based on the temperature field distribution data, and link the heat flow direction with the preset heat conduction channel direction in a chain to generate a direction correction factor.

[0072] It should be noted that in the thermo-mechanical co-simulation of ceramic packaging, accurate identification of the heat flow direction is crucial for optimizing the layout of heat conduction channels. Traditional topology optimization methods typically only aim at minimizing the heat dissipation path or the temperature, neglecting the synergistic relationship between the heat flow direction and the material distribution direction. This leads to a mismatch between the optimized heat conduction channel structure and the actual heat flow direction. To address this issue, this invention, based on obtaining a joint mesh model and performing thermal simulation, further identifies the local heat flow direction based on temperature field distribution data and generates a direction correction factor. Specifically:

[0073] S3.1: Input the joint mesh model and the nodal temperature data obtained from thermal simulation, and input the preset heat conduction channel direction vector. The heat conduction channel direction vector can be a global constant or locally specified based on spatial location.

[0074] The node temperature data is obtained as follows: The processor reads the node temperature data calculated by the thermal simulation solver based on the joint mesh model. The node temperature data is stored in array form, with the array length equal to the total number of nodes. Each array element corresponds to the temperature value of one node. The thermal simulation solution process uses the steady-state heat conduction equation for solution. Figure 5 The simulation environment shown applies JEDEC standard static air convection heat transfer boundary conditions, with the ambient temperature set at 25°C. The PCB uses a JEDEC standard 4-layer (2s2p) structure, resulting in the following... Figure 2 The chip temperature distribution cloud map shown Figure 3 The overall temperature distribution cloud map of the package shown is as follows: Figure 4 The image shows the side temperature distribution cloud map.

[0075] In addition, the direction vector of the heat conduction channel It can be a global constant; for example, for a package structure with a main heat dissipation direction, a global direction vector can be set to point towards the main heat dissipation path. Alternatively, it can be locally specified based on spatial location; for example, a vertically downward direction vector can be set in the area directly below the chip, and a horizontal outward direction vector can be set in the package edge area. The setting of the direction vector is based on factors including the thermal boundary conditions of the package, the distribution of heat sources on the chip, and the external heat dissipation environment. The processor will assign the heat conduction channel direction vector... It is stored as vector field data associated with node coordinates. For the case of global constant, all nodes share the same direction vector; for the case of local specification, each node or each cell region has its own direction vector.

[0076] S3.2: For each grid cell in the joint mesh model, use the nodal temperature data and the derivative of the cell shape function to calculate the temperature gradient at the cell center, and use the opposite direction of the temperature gradient as the cell heat flow direction vector.

[0077] For each grid cell, the processor performs the following sub-steps to calculate the temperature gradient at the cell center:

[0078] Obtain the list of node numbers for this unit, and extract the temperature value of the corresponding node from the node temperature data array based on the node number, denoted as . (For tetrahedral elements); obtain the spatial coordinates of the four nodes of the element, calculate the three edge vectors pointing from node 1 to nodes 2, 3, and 4, and construct the Jacobian matrix. Inverting the Jacobian matrix yields the Jacobian inverse matrix. ; Obtain the derivative matrix of the tetrahedral element shape function with respect to local coordinates. (This is a constant matrix, determined by the element type); The Jacobian inverse matrix... With the derivative matrix of the shape function Multiplying them together yields the derivative matrix of the shape function with respect to the global coordinates (x, y, z). . Transform the node temperature numerical vector With the global derivative matrix of shape functions Multiplying these two components yields the temperature gradient vector at the center of the element: According to Fourier's law of heat conduction, the direction of heat flow is opposite to the direction of the temperature gradient. Therefore, the opposite direction of the temperature gradient vector is taken as the unit heat flow direction vector. The unit vector of heat flow direction of the unit is obtained by normalization. .

[0079] S3.3: For each node, collect all mesh elements containing the node, weight the element heat flow direction vectors by element volume to obtain the local heat flow direction vector of the node, and then normalize it.

[0080] The element volume can be calculated using the coordinates of the four vertices of the tetrahedral element. The weighted average vector is the nodal local heat flow direction vector, which is not yet normalized and its magnitude may not be 1. The processor further normalizes this vector to obtain the nodal local heat flow direction unit vector. .

[0081] It should be noted that by using volume-weighted averaging, the influence of large-volume elements on the nodal orientation is greater, which is consistent with the spatial distribution characteristics of physical field quantities; at the same time, the central quantity of the element is smoothly transitioned to the nodal quantity, avoiding numerical instability caused by abrupt changes in orientation between elements.

[0082] S3.4: For each node, calculate the cosine of the angle between the normalized node local heat flow direction vector and the heat conduction channel direction vector, use the cosine of the angle as the direction deviation coefficient, and use the value obtained by subtracting the direction deviation coefficient from 1 as the node direction correction factor. Output the direction correction factor array for all nodes.

[0083] For each node, the processor reads the heat conduction channel direction vector corresponding to that node. If it is a global constant, then all nodes share the same one. If specified locally, the corresponding direction vector is retrieved based on the node coordinates. Assume... Having been normalized to unit vectors, the processor calculates the dot product of the two vectors, i.e.:

[0084] ;

[0085] The dot product value is the cosine of the angle between the two direction vectors, and its value ranges from [-1, 1]. When the heat flow direction is exactly the same as the preset direction, When the two are in opposite directions, = -1; when the two are perpendicular = 0. The processor defines the cosine of the included angle as the direction deviation coefficient, which reflects the degree of consistency between the heat flow direction at the current node and the design target direction.

[0086] Furthermore, the processor generates a node orientation correction factor based on the orientation deviation coefficient. The calculation formula is: The correction factor ranges from [0, 2]: when the heat flow direction is consistent with the preset direction, =1, = 0 indicates no correction is needed; when the heat flow direction is perpendicular to the preset direction, = 0, = 1 indicates that a moderate level of correction is needed; when the heat flow direction is opposite to the preset direction, = -1, = 2 indicates that the maximum correction is required. The processor organizes the orientation correction factors of all nodes into a one-dimensional array according to the node number order, and outputs the node orientation correction factor array.

[0087] It should be noted that the deviation between the heat flow direction and the channel direction not only affects heat dissipation efficiency but also leads to significant stress concentration. The node orientation correction factor generated in this invention is a quantitative description of this deviation. Through this chain-like data linkage mechanism, this invention transforms the temperature field information obtained from thermal simulation into an orientation correction factor that can directly guide topology optimization, achieving closed-loop feedback from thermal analysis results to structural optimization design.

[0088] S4: Based on the direction correction factor, the layout direction of the heat conduction channel is iteratively adjusted until optimized channel structure data that satisfies the requirement of the shortest heat dissipation path and the stress value being lower than the stress limit threshold is obtained.

[0089] First, a density design variable is assigned to each element in the joint mesh model. Specifically, the processor iterates through all mesh elements in the generated joint mesh model and initializes a density design variable for each mesh element. Among them, the range of values ​​for the density design variable is: ,in A value of 0.001 is used to avoid numerical singularities. In the initial iteration step, the density design variable for all elements is set to 0.5, indicating that the material distribution is in a completely grayscale state. A solid isotropic material penalty (SIMP) model is used to map the density design variable to the material properties of the elements, such as element thermal conductivity.

[0090] ;

[0091] Elemental elastic modulus:

[0092] ;

[0093] in, The penalty factor is 3 in this embodiment. and These represent the thermal conductivity and elastic modulus of the solid material, respectively. The density design variable will be continuously updated in subsequent iterations, eventually converging to 0 or 1, forming a clear 0-1 material distribution.

[0094] S4.1: Apply convective heat transfer boundary conditions and mechanical boundary conditions to the joint mesh model, perform thermo-mechanical coupled finite element solution, and obtain nodal temperature field data and nodal stress field data.

[0095] First, based on the density design variables of each element, the thermal conductivity and elastic modulus of each element are calculated using the SIMP model, and the material property array in the finite element model is updated.

[0096] Subsequently, boundary conditions are applied. The thermal boundary conditions utilize a JEDEC standard static air environment, applying convective heat transfer boundaries to all outer surfaces of the model. The heat transfer coefficient is set according to the JEDEC standard (e.g., 5~10 W / m²K under natural convection conditions), and the ambient temperature is set to 25℃. For the bottom of the PCB or fixed areas, fixed temperature boundaries or adiabatic boundaries can be applied. Mechanical boundary conditions are set according to the actual assembly of the package, such as applying fixed constraints to the bottom of the package and applying uniform pressure or free boundaries to the chip surface. Furthermore, thermal strain loads are applied, using the temperature field obtained in the current iteration as the thermal load, and calculating the thermal strain using the coefficient of thermal expansion. The processor calls the finite element solver to simultaneously solve the steady-state heat conduction equation and the linear elasticity equation, obtaining nodal temperature field data (a one-dimensional array with a length equal to the total number of nodes) and nodal stress field data (each node contains six stress components or von Mises equivalent stress).

[0097] It should be noted that this invention uses a thermo-mechanical coupling solution, which allows subsequent optimization to simultaneously consider heat dissipation performance and structural strength, avoiding the separate process of optimization followed by verification in traditional methods.

[0098] S4.2: Using the node temperature field data of the current iteration as input, calculate the node orientation correction factor array for the current iteration.

[0099] That is, after obtaining the node temperature field data of the current iteration step, the processor takes it as input, executes steps S3.2 to S3.4 in their entirety, outputs the node orientation correction factor array of the current iteration step, and simultaneously obtains the element heat flow direction vector of each element.

[0100] S4.3: Interpolate the nodal orientation correction factor array to the grid cells. For each grid cell, take the average value of the orientation correction factors of the nodes it contains to obtain the cell orientation correction factor.

[0101] It should be noted here that, since the direction penalty term in the subsequent objective function needs to be defined on the cell, the processor needs to convert the direction correction factor on the node to the cell.

[0102] For each element in the joint mesh model, the processor obtains the global IDs of all nodes constituting that element, reads the corresponding correction factor values ​​from the node orientation correction factor array of the current iteration step based on the node IDs, calculates the arithmetic mean of these node correction factors, and uses it as the element orientation correction factor for that element. .

[0103] For tetrahedral elements, the average value of four nodes is used; for hexahedral elements, the average value of eight nodes is used. This element orientation correction factor reflects the degree of deviation between the heat flow direction in the element's region and the preset direction: if... If the value is close to 0, it indicates that the heat flow direction in this region is consistent with the preset direction; if... A value close to 2 indicates that the two directions are opposite.

[0104] S4.4: Project the density design variables onto the nodes, calculate the node density gradient, and then interpolate it back to the element center to obtain the element density gradient vector. Normalize the vector to obtain the element density gradient direction vector.

[0105] Due to density design variables Since the density is defined at the cell center (i.e., a constant value for each cell), the density gradient inside the cell cannot be directly calculated. Therefore, the processor first projects the cell density onto the nodes: for each node, it collects all cells containing that node and sets the density design variables for those cells. The density value of the node is obtained by weighting the density by the element volume. Furthermore, for the nodal density values, the density gradient at the center of each element is calculated using the finite element shape function derivative: for each element, the density gradient vector at the element center is calculated using the nodal coordinates and nodal density values ​​through the shape function derivative matrix. The specific calculation method is similar to that in S3.2 for calculating the temperature gradient; simply replace the node temperature with the node density.

[0106] After obtaining the density gradient vector, normalization is performed to obtain the unit vector of the cell density gradient direction:

[0107] ;

[0108] in This is the density gradient magnitude. This direction vector represents the direction of the fastest change in material density, i.e., the orientation of the material distribution. In topology optimization, the density gradient direction can be considered as the axial direction of the heat conduction channel.

[0109] S4.5: Construct a comprehensive objective function that includes the shortest heat dissipation path term, the direction penalty term, and the stress constraint penalty term.

[0110] Preferably, the processor constructs a synthesis objective function of the following form:

[0111] ;

[0112] in, The option with the shortest heat dissipation path. For directional penalty items, This is a stress constraint penalty term. and These are the directional penalty weighting coefficient and the stress penalty coefficient (both are preset normal values ​​that can be adjusted through numerical experiments).

[0113] Among them, the shortest heat dissipation path item The integral of the square of the nodal temperature is used as a measure of heat dissipation weakness, and its discrete form is: ,in Let n be the temperature of node n. The node control volume is obtained by distributing the volume among adjacent cells. Minimizing this objective corresponds to reducing the overall temperature level and shortening the heat flow path.

[0114] Directional penalty item This is used to penalize regions where the direction of heat flow does not align with the direction of material distribution. First, define the element direction penalty coefficient: ,in The unit vector of the element heat flux direction calculated in S3.2. The unit vector of the cell density gradient direction calculated in S4.4, dot product , representing the cosine of the angle between the heat flow direction and the material gradient direction, when the two directions are consistent. =1, =0; when the two directions are perpendicular =0, = When the two are in opposite directions =-1, e=2 .at the same time, The deviation between the heat flow and the preset direction has already been reflected, therefore It penalizes two types of deviations: the degree to which heat flow deviates from the preset direction and the degree to which material distribution deviates from the direction of heat flow.

[0115] Directional penalty term is defined as ,in For the unit density gradient mode, This refers to the unit volume. This is to make the penalty related to the interface area, avoiding unnecessary penalties in regions with uniform density. It should be noted that the guided optimization process ensures that the direction of the heat conduction channel (characterized by the density gradient) is consistent with the direction of heat flow, while the direction of heat flow is constrained by a preset direction, thereby achieving synergy between heat flow and structure while meeting the design intent.

[0116] Stress constraint penalty term This is used to limit the maximum stress in the optimized structure to not exceed the material's yield strength. The processor extracts the maximum von Mises stress value for each element from the nodal stress field solved in S4.1. (The maximum value of the stress at the element integration point or nodal point can be taken). When the direction of the heat conduction channel is opposite to the direction of heat flow, the stress at the root of the fin increases significantly, therefore stress constraint is necessary. The penalty function method is used to transform the stress constraint into an objective term: ,in This is the stress limit threshold, i.e., the yield strength of the material. This term applies when the element stress approaches or exceeds the stress limit threshold. The stress increases dramatically, forcing optimization algorithms to adjust material distribution to reduce stress.

[0117] S4.6: The sensitivity of the integrated objective function to the density design variables is calculated using the adjoint method. Based on the sensitivity, the density design variables are updated using the optimization criterion method or the moving asymptote method to obtain new density design variables.

[0118] It should be noted that, due to It implicitly depends on the temperature field and stress field, which in turn are related to the density design variables through the finite element equation. The direct differential calculation is computationally intensive. Therefore, this invention employs the adjoint method for sensitivity analysis. First, the adjoint equations for the heat conduction problem and the linear elasticity problem are constructed. The contributions of the temperature field and stress field to the objective function are obtained through a first adjoint solution, and then analytically calculated. .

[0119] The specific implementation of the adjoint method is well known to those skilled in the art and will not be elaborated here.

[0120] After obtaining the sensitivity, the processor updates the design variables using either the Optimization Criterion Method (OC) or the Moving Asymptote Method (MMA). The Optimization Criterion Method is suitable for problems such as minimum compliance under volume constraints, and the update format is as follows:

[0121] ;

[0122] in Using Lagrange multipliers, the volume constraint is satisfied by solving the bisection method. This is the damping coefficient (usually taken as 0.5). This represents the element density design variable value for the current iteration step (i.e., after the previous iteration step has been updated). The MMA method is applicable to more general objective functions and constraints, solving them iteratively by constructing convex approximate subproblems. In this embodiment, the MMA method can be used to ensure convergence stability. The updated density design variable is then obtained. Density filtering and projection techniques are applied to eliminate grayscale cells and checkerboard patterns.

[0123] Finally, determine whether the maximum change of the new density design variable from the density design variable in the previous iteration step is less than the optimization convergence tolerance (e.g., 0.01). If so, terminate the iteration and identify the mesh elements with a density value greater than the preset density threshold (e.g., 0.5) as material elements (i.e., heat conduction channels), and output the heat conduction channel structure data composed of material elements. If not, return to perform thermal-mechanical coupled finite element solution, i.e., return to S4.1, and continue iterative adjustment.

[0124] In summary, this invention achieves topology optimization through multi-physics coupling of heat, force, and power, automatically aligning the direction of the heat conduction channel with the direction of heat flow, and strictly controlling the stress level, ultimately obtaining an optimized structure that combines efficient heat dissipation with high reliability.

[0125] like Figure 6 As shown, the present invention also provides a thermal-mechanical co-simulation data processing system for ceramic packaging, including: a data slice division module, which acquires a three-dimensional model and material parameters of the ceramic package, wherein the material parameters include the thermal conductivity of each material, and divides the three-dimensional model into high heat flux data slices and low heat flux data slices according to the thermal conductivity threshold;

[0126] The joint mesh generation module generates first mesh data by using fine mesh division for the high heat flux data patch and second mesh data by using coarse mesh division for the low heat flux data patch, and combines the first mesh data and the second mesh data to form a joint mesh model;

[0127] The correction factor calculation module performs thermal simulation based on the joint mesh model to obtain temperature field distribution data, identifies the heat flow direction based on the temperature field distribution data, and links the heat flow direction with the preset heat conduction channel direction in a chain to generate a direction correction factor.

[0128] The channel structure optimization module iteratively adjusts the layout direction of the heat conduction channel based on the direction correction factor until optimized channel structure data that satisfies the requirement of the shortest heat dissipation path and stress value below the stress limit threshold is obtained.

[0129] The system also includes one or more processors and memory.

[0130] The memory is used to store operable instructions that, when executed by the one or more processors, cause the one or more processors to perform operations, including the flow of the personalized nutrition management method based on multimodal data fusion and deep learning described in the foregoing embodiments, especially... Figure 1 The flowchart of the method is shown.

[0131] Other aspects disclosed in the embodiments of the present invention also propose a computer-readable medium for storing software including instructions executable by one or more computers, which, upon execution, cause the one or more computers to perform operations including the flow of the personalized nutrition management method based on multimodal data fusion and deep learning of the foregoing embodiments, particularly... Figure 1 The flowchart of the method is shown.

[0132] It should be recognized that embodiments of the present invention may be implemented or carried out by computer hardware, a combination of hardware and software, or by computer instructions stored in a non-transitory computer-readable storage medium.

[0133] The method can be implemented using standard programming techniques, including a non-transitory computer-readable storage medium configured with a computer program in the computer program, wherein the storage medium is configured such that the computer operates in a specific and predefined manner.

[0134] Each program can be implemented in a high-level procedural or object-oriented programming language to communicate with the computer system; however, if required, the program can be implemented in assembly or machine language.

[0135] In any case, the language can be either compiled or interpreted.

[0136] Furthermore, for this purpose, the program can run on a programmed application-specific integrated circuit.

[0137] The processes described herein (or variations and / or combinations thereof) can be executed under the control of one or more computer systems configured with executable instructions, and can be implemented by hardware or a combination thereof as code (e.g., executable instructions, one or more computer programs, or one or more applications) that commonly executes on one or more processors. The computer program includes a plurality of instructions executable by one or more processors.

[0138] Furthermore, the method can be implemented in any suitable computing platform, including but not limited to personal computers, minicomputers, mainframes, workstations, networked or distributed computing environments, standalone or integrated computer platforms, or in communication with charged particle tools or other imaging devices.

[0139] Various aspects of the present invention can be implemented in machine-readable code stored on a non-transitory storage medium or device, whether portable or integrated into a computing platform, such as a hard disk, optical read and / or write storage medium, RAM, ROM, etc., such that it can be read by a programmable computer, and when the storage medium or device is read by the computer, it can be used to configure and operate the computer to perform the processes described herein.

[0140] Furthermore, machine-readable code, or parts thereof, can be transmitted via wired or wireless networks.

[0141] When such media includes instructions or programs that combine with a microprocessor or other data processor to implement the steps described above, the invention described herein includes these and other different types of non-transitory computer-readable storage media.

[0142] It should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit it. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical solutions of the present invention without departing from the spirit and scope of the technical solutions of the present invention, and all such modifications or substitutions should be covered within the scope of the claims of the present invention.

Claims

1. A method for thermal-mechanical co-simulation data processing for ceramic packaging, characterized in that, include: A three-dimensional model and material parameters of the ceramic package are obtained, including the thermal conductivity of each material. The three-dimensional model is then divided into high heat flux data slices and low heat flux data slices based on the thermal conductivity threshold. The material identifiers of each entity unit in the 3D model of the ceramic package are read. Based on these material identifiers, the corresponding thermal conductivity data is retrieved from the material parameter library and assigned to the corresponding entity unit. The median of the thermal conductivity data for all entity units is calculated and set as the thermal conductivity threshold. The thermal conductivity data of each entity unit is compared with the thermal conductivity threshold. Entity units with thermal conductivity data greater than the thermal conductivity threshold are marked as high heat flux units, and entity units with thermal conductivity data less than or equal to the thermal conductivity threshold are marked as low heat flux units. The division of high-heat-flux data slices and low-heat-flux data slices includes: grouping identically labeled entity units sharing a face or edge into continuous regions through neighborhood search to form high-heat-flux unit groups and low-heat-flux unit groups; for each high-heat-flux unit group and each low-heat-flux unit group, traversing the faces of all entity units within the group, extracting the face owned by only one unit as the boundary triangle facet of the corresponding group, recording the three vertex numbers and vertex coordinates of the boundary triangle facet to form the boundary surface data of the corresponding group; and simultaneously recording the index list of all entity units within the group to form data slice identification data. Read the data slice identification data and the boundary surface data to identify the shared boundary triangle facets between different data slices, and mark the regions where the boundary triangle facets are located as conformal boundaries; simultaneously, determine the spatial range of each data slice based on the boundary surface data; read the unit index list in the data slice identification data to locate the original solid units of each data slice in the initial 3D model, and extract the geometric dimension data of the original solid units: for high heat flux data slices, take the smallest side length of the original solid units contained therein. One-tenth is used as the target side length for the fine mesh; for low heat flux data patches, the average feature size of the original solid cells is taken. One-tenth of the target side length is used as the coarse-grained mesh target side length; this target side length is only used for mesh generation in this preprocessing stage and remains unchanged in subsequent optimization iterations; for each conformal boundary, the shared boundary triangle facet is used as a geometric reference, and the shared boundary triangle facet is re-divided using the fine mesh target side length to generate a conformal surface mesh; the nodes on the conformal surface mesh completely overlap between adjacent data slices; for each high heat flux data slice and each low heat flux data slice, the boundary triangle facet on the non-conformal boundary is used as a geometric reference, and the data slice surface mesh is re-divided using the target side length corresponding to the data slice to generate a data slice surface mesh; the data slice surface mesh belongs to only a single data slice; The high heat flux data slice is divided into a fine grid to generate the first grid data, and the low heat flux data slice is divided into a coarse grid to generate the second grid data. The first grid data and the second grid data are then combined to form a joint grid model. Thermal simulation is performed based on the joint mesh model to obtain temperature field distribution data. The heat flow direction is identified based on the temperature field distribution data. The heat flow direction is then linked with the preset heat conduction channel direction in a chain to generate a direction correction factor. Based on the aforementioned direction correction factor, the layout direction of the heat conduction channel is iteratively adjusted until optimized channel structure data is obtained that satisfies the requirement of the shortest heat dissipation path and a stress value lower than the stress limit threshold.

2. The thermal-mechanical co-simulation data processing method for ceramic packaging according to claim 1, characterized in that, The formation of the joint mesh model includes: For each high heat flux data slice, the conformal surface mesh portion and the data slice surface mesh portion are spliced ​​together to form the complete boundary mesh of the corresponding data slice. Using the complete boundary mesh as the boundary, the fine mesh target side length is used as the internal mesh size control to generate a tetrahedral mesh and obtain the first mesh data. For each low heat flux data slice, the conformal surface mesh portion and the data slice surface mesh portion are spliced ​​together to form the complete boundary mesh of the corresponding data slice. Using the complete boundary mesh as the boundary, the coarse-grained mesh target side length is used as the internal mesh size control to generate a tetrahedral mesh and obtain the second mesh data. The first grid data and the second grid data are integrated into the same grid file, and the node correspondence on the conformal boundary is recorded to form a joint grid model.

3. The thermal-mechanical co-simulation data processing method for ceramic packaging according to claim 2, characterized in that, The identification of heat flow direction includes: Input the joint mesh model and the node temperature data obtained from thermal simulation, and input the preset heat conduction channel direction vector; The direction vector of the heat conduction channel is a global constant or locally specified according to the spatial location; For each grid cell in the joint grid model, the temperature gradient at the cell center is calculated using the node temperature data and the derivative of the cell shape function, and the opposite direction of the temperature gradient is used as the cell heat flow direction vector.

4. The thermal-mechanical co-simulation data processing method for ceramic packaging according to claim 3, characterized in that, The process of chaining data and generating a direction correction factor includes: For each node, all mesh cells containing the node are collected, and the heat flow direction vectors of the cells are weighted and averaged according to the cell volume to obtain the local heat flow direction vector of the node, and then normalized. For each node, calculate the cosine of the angle between the normalized node local heat flow direction vector and the heat conduction channel direction vector, use the cosine of the angle as the direction deviation coefficient, and then subtract the direction deviation coefficient from 1 as the node direction correction factor. Output the orientation correction factor array for all nodes.

5. The thermal-mechanical co-simulation data processing method for ceramic packaging according to claim 4, characterized in that, Iterative adjustments to the layout orientation of the heat conduction channels include: Assign a density design variable to each element in the joint mesh model; Convective heat transfer boundary conditions and mechanical boundary conditions are applied to the joint mesh model, and thermo-mechanical coupled finite element solution is performed to obtain nodal temperature field data and nodal stress field data. Using the node temperature field data of the current iteration as input, calculate the node orientation correction factor array for the current iteration; The node orientation correction factor array is interpolated to the grid cell, and for each grid cell, the average value of the orientation correction factors of the nodes contained therein is taken to obtain the cell orientation correction factor; The density design variables are projected onto the nodes, the node density gradient is calculated, and then interpolated back to the cell center to obtain the cell density gradient vector. The normalized vector is then used to obtain the cell density gradient direction vector.

6. The thermal-mechanical co-simulation data processing method for ceramic packaging according to claim 5, characterized in that, Iterative adjustments to the layout orientation of the heat conduction channels also include: Construct a comprehensive objective function that includes a shortest heat dissipation path term, a direction penalty term, and a stress constraint penalty term; wherein, the direction penalty term is constructed based on the element direction correction factor, the element heat flow direction vector, and the element density gradient direction vector; The sensitivity of the integrated objective function to the density design variable is calculated using the adjoint method. Based on the sensitivity, the density design variable is updated using the optimization criterion method or the moving asymptote method to obtain the new density design variable. Determine whether the maximum change between the new density design variable and the density design variable in the previous iteration step is less than the optimization convergence tolerance. If so, terminate the iteration and determine the grid cells with density values ​​greater than the preset density threshold as material cells, and output the heat conduction channel structure data formed by the material cells. If not, return to the thermo-mechanical coupled finite element solution and continue iterative adjustments.

7. A thermal-mechanical co-simulation data processing system for ceramic packaging, based on the thermal-mechanical co-simulation data processing method for ceramic packaging according to any one of claims 1 to 6, characterized in that: The data chip segmentation module acquires the three-dimensional model and material parameters of the ceramic package, including the thermal conductivity of each material, and divides the three-dimensional model into high heat flux data chips and low heat flux data chips according to the thermal conductivity threshold. The joint mesh generation module generates first mesh data by using fine mesh division for the high heat flux data patch and second mesh data by using coarse mesh division for the low heat flux data patch, and combines the first mesh data and the second mesh data to form a joint mesh model; The correction factor calculation module performs thermal simulation based on the joint mesh model to obtain temperature field distribution data, identifies the heat flow direction based on the temperature field distribution data, and links the heat flow direction with the preset heat conduction channel direction in a chain to generate a direction correction factor. The channel structure optimization module iteratively adjusts the layout direction of the heat conduction channel based on the direction correction factor until optimized channel structure data that satisfies the requirement of the shortest heat dissipation path and stress value below the stress limit threshold is obtained.

Images