Integrated circuit packaging quality assessment method and system based on nondestructive testing technology

Abstract

The present invention provides an integrated circuit packaging quality assessment method and system based on nondestructive testing technology, which relates to the field of integrated circuit technology, including using a two-photon excitation microscope system to obtain a three-dimensional fluorescence response signal of a packaging structure, and establishing a stress distribution characteristic map and a stress gradient classification map by analyzing the spatial distribution of the fluorescence signal; arranging a quantum nitrogen vacancy probe array in a stress abnormality area, collecting its electromagnetic field response signal under stress, and realizing accurate mapping of the stress field and the electromagnetic field of the packaging structure; based on the acquired multi-physical field data, establishing a group of thermo-electromagnetic coupling equations, and calculating the interaction relationship between each physical field; by analyzing the evolution characteristics of each detection area under the multi-field coupling effect, establishing a stability evaluation system including mechanical failure risk and electromagnetic failure risk, and realizing accurate assessment of the quality level of the packaging structure.

Description

Technical Field

[0001] The present invention relates to the field of integrated circuit technology, and in particular to an integrated circuit packaging quality assessment method and system based on non-destructive testing technology. Background Art

[0002] With the rapid development of integrated circuit technology, the integration of chips is constantly improving, and the packaging structure is becoming more and more complex. As an important interface between the chip and the external environment, the quality of the package is directly related to the performance and reliability of the integrated circuit. In practical applications, the packaging structure is subjected to complex thermal stress, mechanical stress and electromagnetic field effects. The distribution characteristics and mutual coupling effects of these physical fields have an important impact on the stability of chip functions.

[0003] There are some problems with the current packaging quality assessment methods. The detection means are single, and it is difficult to simultaneously obtain multi-physical field information such as stress field and electromagnetic field in the packaging structure; the detection accuracy is limited, and it is impossible to accurately identify abnormal stress areas at the micron and sub-micron scales; there is a lack of systematic analysis of the coupling effects of multiple physical fields, making it difficult to accurately assess the reliability risks of the packaging structure in the actual working environment; the quality assessment standards are not perfect, and it is difficult to establish a comprehensive evaluation system covering mechanical failure and electromagnetic failure.

[0004] In summary, there is an urgent need for an integrated circuit packaging quality assessment method to solve the technical problems of the existing packaging quality assessment methods, such as the single detection means, insufficient accuracy, and lack of multi-field coupling analysis; to achieve high-precision synchronous measurement of the stress field and electromagnetic field of the packaging structure; to accurately assess the failure risk of each detection area; and finally to establish a stability evaluation system to provide a comprehensive and accurate technical solution for the quality assessment of the integrated circuit packaging structure. The present invention can solve the problems in the prior art. Summary of the invention

[0005] The embodiments of the present invention provide a method and system for evaluating the quality of integrated circuit packaging based on nondestructive testing technology, which can solve the problems in the prior art.

[0006] According to a first aspect of the embodiments of the present invention,

[0007] A method for evaluating the quality of integrated circuit packaging based on nondestructive testing technology is provided, comprising:

[0008] The package structure of the integrated circuit to be inspected is placed on the sample stage of the two-photon excitation microscope system, a femtosecond laser is used to generate a two-photon excitation pulse, and the energy density of the two-photon excitation pulse is controlled in real time by an acousto-optic modulator, so that the package structure generates a stable fluorescence response signal; the confocal scanning unit is controlled to scan the package structure layer by layer along the Z-axis direction to obtain the three-dimensional distribution data of the fluorescence response signal; a stress distribution characteristic map of the package structure is established according to the three-dimensional distribution data, and the area where the stress gradient exceeds the preset threshold is marked in the stress distribution characteristic map to generate a stress anomaly location map; the stress gradient value is calculated for each area in the stress anomaly location map to establish a stress gradient grading map;

[0009] Determine a detection area according to the stress gradient grading map, arrange a quantum nitrogen vacancy probe array in the detection area, and set radio frequency pulse sequence parameters according to the spatial distribution of the quantum nitrogen vacancy probe array so that the quantum nitrogen vacancy probe array achieves a corresponding detection sensitivity in each detection area; collect the electromagnetic field response signal generated by the quantum nitrogen vacancy probe array under stress, reconstruct the electromagnetic field distribution map of the packaging structure, and spatially correspond it with the stress gradient grading map to form a stress-electromagnetic field coupling model;

[0010] Based on the stress-electromagnetic field coupling model, a group of thermomechanical and electromagnetic coupling equations for the packaging structure is established; the evolution characteristics of the group of thermomechanical and electromagnetic coupling equations in each detection area are calculated to obtain the interaction relationship between the physical fields; based on the interaction relationship, the failure probability of each detection area under the multi-field coupling is calculated, and a stability evaluation system including mechanical failure risk and electromagnetic failure risk is established; the quality grade of the packaging structure is determined by comparing with the stability evaluation system, and a quality assessment report of the packaging structure is generated.

[0011] In an optional embodiment, establishing a stress distribution characteristic diagram of the packaging structure according to the three-dimensional distribution data includes:

[0012] Inputting the three-dimensional distribution data of the fluorescence response signal into an optical system model, wherein the optical system model includes a lateral resolution parameter determined by the numerical aperture of the objective lens and a longitudinal resolution parameter determined by the excitation wavelength;

[0013] Generate a Gaussian point spread function according to the optical system model, perform a deconvolution operation on the three-dimensional distribution data and the Gaussian point spread function, and perform Wiener filtering on the result of the deconvolution operation to obtain three-dimensional fluorescence distribution data after eliminating the influence of the optical system;

[0014] Measure the photoelastic coefficient tensor and refractive index anisotropy coefficient of the material corresponding to the packaging structure, and establish the photoelastic effect response equation of the material stress field and refractive index change;

[0015] Establishing a nonlinear relationship equation between the refractive index and the fluorescence intensity during two-photon excitation, combining the nonlinear relationship equation with the photoelastic effect response equation to obtain a quantitative mapping model between the stress field and the fluorescence intensity;

[0016] Selecting an initial stress field, calculating a theoretical fluorescence intensity according to the quantitative mapping model, comparing the theoretical fluorescence intensity with the measured fluorescence intensity of the three-dimensional fluorescence distribution data, and calculating a mean square error as an optimization objective function; performing a gradient calculation, calculating a gradient vector of the optimization objective function, determining a conjugate search direction, performing a conjugate search, and iteratively searching for an optimal step size;

[0017] The stress field distribution is updated using the optimal step size. When the convergence value of the optimization objective function is greater than a preset convergence threshold, the gradient calculation and conjugate search are continued. When the convergence value of the optimization objective function is less than or equal to the convergence threshold, the current stress field distribution is used as a stress distribution characteristic diagram of the packaging structure.

[0018] In an optional embodiment, marking an area where the stress gradient exceeds a preset threshold in the stress distribution characteristic map, and generating a stress anomaly location map includes:

[0019] Using the Sobel operator to calculate the gradient of the stress distribution characteristic diagram in the horizontal direction to obtain a first gradient component matrix, and using the Sobel operator to calculate the gradient of the stress distribution characteristic diagram in the vertical direction to obtain a second gradient component matrix;

[0020] Calculating the square sum of the first gradient component matrix and the second gradient component matrix, and performing a square root operation to obtain a gradient matrix of the stress distribution characteristic diagram; performing a power-of-2 multi-scale decomposition on the gradient matrix, wherein the multi-scale decomposition ranges from 2×2 to 16×16, and obtaining sub-matrices of multiple scales;

[0021] Calculating the grayscale mean and grayscale standard deviation of each local area in the sub-matrices of the multiple scales respectively, determining the segmentation threshold of the corresponding local area according to the grayscale mean and the grayscale standard deviation, and determining the adaptive segmentation threshold;

[0022] Based on an adaptive segmentation threshold, the sub-matrices of the multiple scales are binarized respectively to obtain binarized sub-matrices of the multiple scales; the grayscale variance of each binarized sub-matrix of the multiple scales is calculated respectively, and the grayscale variance is normalized to obtain the fusion weight of each scale; the binarized sub-matrices of the multiple scales are multiplied and superimposed with the corresponding fusion weights respectively to obtain a fused segmentation matrix; the fused segmentation matrix is ​​binarized to obtain a mask matrix, and in the mask matrix, the stress abnormality area is determined based on the area with a value of 1; the stress abnormality area is marked in the stress distribution characteristic map to generate a stress abnormality location map.

[0023] In an optional embodiment, determining a detection area according to the stress gradient grading map, arranging a quantum nitrogen vacancy probe array in the detection area, and setting radio frequency pulse sequence parameters according to the spatial distribution of the quantum nitrogen vacancy probe array so that the quantum nitrogen vacancy probe array achieves a corresponding detection sensitivity in each detection area includes:

[0024] Constructing a quantum state evolution measurement system, applying an initialization radio frequency pulse to the quantum state nitrogen vacancy probe array, and initializing the quantum state nitrogen vacancy probe array to a spin ground state;

[0025] Using a quantum state readout sequence to obtain the spin state information of the quantum state nitrogen vacancy probe array, record the spin state evolution time, measure the quantum coherence time of the quantum state nitrogen vacancy probe array in each detection area, and establish a quantum state response model;

[0026] Using a quantum state subspace decomposition method, the quantum state response model is characteristically decomposed to obtain a quantum state characteristic vector of each detection area; a quantum decoherence noise estimator is constructed based on the quantum state characteristic vector, the influence of environmental noise on the quantum state nitrogen vacancy probe of each detection area is measured, and the quantum state decoherence rate of the detection area is calculated;

[0027] According to the quantum state decoherence rate, the minimum detection sensitivity corresponding to each detection area is determined, and a multi-channel quantum state manipulation timing scheme is established to achieve synchronous quantum state manipulation between detection areas;

[0028] For each detection area, a quantum state characteristic transformation matrix is ​​constructed to extract quantum coherence information of the quantum state nitrogen vacancy probe to obtain the quantum state response characteristics of the detection area;

[0029] Inputting the quantum state response characteristics into the quantum information fusion framework, and calculating the sensitivity parameters of each detection area based on the quantum state evolution analysis algorithm;

[0030] Based on the sensitivity parameters, the radio frequency pulse phase accumulation time corresponding to each detection area is calculated, the radio frequency pulse repetition period is determined, and radio frequency pulse sequence parameters that meet the detection sensitivity requirements of each detection area are generated.

[0031] In an optional embodiment, collecting the electromagnetic field response signal generated by the quantum nitrogen vacancy probe array under stress, reconstructing the electromagnetic field distribution diagram of the packaging structure, and spatially corresponding it with the stress gradient classification diagram to form a stress-electromagnetic field coupling model includes:

[0032] The electromagnetic field response signal of the quantum nitrogen vacancy probe array is collected in each detection area to construct an electromagnetic field response time series; the electromagnetic field response time series is discretely sampled, a signal analytical expression is obtained by Hilbert transform, and the instantaneous value of the electromagnetic field intensity in each detection area is extracted; the instantaneous value of the electromagnetic field intensity is constructed as an electromagnetic field response matrix, and the main eigenvectors of the electromagnetic field distribution are extracted by principal component analysis method;

[0033] Using continuous wavelet transform to perform time-frequency decomposition on the electromagnetic field response matrix to obtain electromagnetic field spectrum components at different scales; performing adaptive threshold filtering on the electromagnetic field spectrum components to eliminate measurement noise and obtain frequency domain characteristics of the electromagnetic field signal;

[0034] A compressed sensing reconstruction algorithm is used to sparsely reconstruct the frequency domain features, a mathematical model of electromagnetic field spatial distribution is established, an electromagnetic field inversion equation is constructed based on the mathematical model, and a conjugate gradient iteration method is used to solve the electromagnetic field distribution to obtain an electromagnetic field distribution map of the packaging structure; the electromagnetic field distribution map is affinely transformed and registered with the stress gradient classification map to establish a spatial coordinate correspondence relationship;

[0035] Based on the spatial coordinate correspondence, a kernel function regression method is used to extract the mapping relationship between the electromagnetic field intensity and the stress gradient, and the electromagnetic field-stress response data is obtained; based on the electromagnetic field-stress response data, a nonlinear coupling equation group is established, and the coupling coefficient is solved by the Newton iteration method;

[0036] Substituting the coupling coefficient into the Maxwell stress-electromagnetic field tensor equation, a stress-electromagnetic field coupling model of the packaging structure is formed.

[0037] In an optional embodiment, a group of thermomechanical and electromagnetic coupling equations of the packaging structure is established based on the stress-electromagnetic field coupling model; the evolution characteristics of the group of thermomechanical and electromagnetic coupling equations in each detection area are calculated to obtain the interaction relationship between the physical fields, and the failure probability of each detection area under the multi-field coupling is calculated according to the interaction relationship, and a stability evaluation system including mechanical failure risk and electromagnetic failure risk is established; the quality level of the packaging structure is determined by comparing with the stability evaluation system, and the quality assessment report of the packaging structure is generated, including:

[0038] The temperature field is introduced into the stress-electromagnetic field coupling model to construct a stress distribution function including a thermal expansion coefficient; a corresponding relationship between temperature and conductivity is established based on the stress distribution function, and an exponential function including a parameter of the ratio of activation energy to the Boltzmann constant is used to determine the change characteristics of conductivity with temperature;

[0039] The stress distribution function is combined with a predetermined electromagnetic field basic equation group, a trial function is selected to construct a weight function, a discretization equation is established by piecewise linear interpolation, and a thermal electromagnetic coupling equation group of the packaging structure is obtained by minimizing the residual error;

[0040] The solution time interval is divided into multiple subintervals, and four iterative calculations including midpoint estimation and full interval correction are performed in each subinterval. The calculation step size is dynamically adjusted based on the local error of the subinterval to obtain the numerical solution of the thermo-electromagnetic coupling equations.

[0041] Performing signal decomposition operation on the numerical solution to extract characteristic mode functions that characterize the inherent change characteristics of the physical field of the packaging structure; calculating the ratio of the joint probability density to the edge probability density of the characteristic mode functions to obtain the interaction strength between the physical fields in the packaging structure;

[0042] By calculating the predicted force of the known physical field on the target physical field, a quantitative characterization model of the interaction between the physical fields of the packaging structure is established, and based on the quantitative characterization model, a failure judgment function including the physical field intensity threshold and the load response is constructed;

[0043] Based on the expectation value minimization criterion, an optimized probability density function is constructed, and multiple random sampling calculations are performed according to the failure judgment function to obtain the failure probability of the packaging structure;

[0044] Establishing a stability evaluation index system based on the failure probability of the packaging structure, and calculating a membership relationship matrix of each stability evaluation index to the quality grade;

[0045] The weight coefficient vector of the stability evaluation index is calculated by pairwise comparison, and the weight coefficient vector is multiplied by the membership matrix to obtain the stability evaluation result of the packaging structure.

[0046] In an optional embodiment, the solution time interval is divided into multiple subintervals, four iterative calculations including midpoint estimation and full interval correction are performed in each subinterval, and the calculation step size is dynamically adjusted based on the local error of the subinterval to obtain the numerical solution of the thermo-electromagnetic coupling equations, including:

[0047] Divide the solution time interval into a plurality of subintervals, set an initial time step at the initial time of the subinterval; at the initial time of the subinterval, use an explicit Euler format to calculate the initial value of the physical field of the thermo-electromagnetic coupling equation group;

[0048] Calculate the predicted value of the physical field at the intermediate time, the formula is as follows:

[0049] ;

[0050] in, t i Indicates the current moment, t i+1 / 2 represents the intermediate moment of advancing half a step, u * (·) represents the predicted value of the physical field, u (·) represents the physical field value, Δ t represents the time step, f (·) represents the value of the physical field function;

[0051] Calculate the predicted value of the physical field at the end time, the formula is as follows:

[0052] ;

[0053] in, u * ( t i+1 ) represents the predicted value of the physical field at the termination time;

[0054] Calculate the physical field correction value at the end time, the formula is as follows:

[0055] ;

[0056] in, u ' (·) represents the physical field correction value;

[0057] Calculate the local error estimate as follows:

[0058] ;

[0059] in, E represents the local error estimate, ||·|| represents the norm;

[0060] Compare the local error estimate with a preset error threshold, and when the local error estimate is greater than the error threshold, halve the initial time step and return to execute the calculation step;

[0061] When the local error estimate is less than one tenth of the error threshold, the initial time step is doubled and the physical field correction value at the termination time is used as the physical field initial value of the next subinterval;

[0062] When the local error estimate is between one tenth of the error threshold and the error threshold, the initial time step is kept unchanged and the physical field correction value at the termination time is used as the physical field initial value of the next sub-interval;

[0063] The calculation steps are repeatedly performed on the subintervals until the calculation of the entire solution time interval is completed to obtain the numerical solution of the thermal electromagnetic coupling equations.

[0064] According to a second aspect of the embodiments of the present invention,

[0065] Provided is an integrated circuit packaging quality assessment system based on non-destructive testing technology, comprising:

[0066] The first unit is used to place the package structure of the integrated circuit to be detected on the sample stage of the two-photon excitation microscope system, use a femtosecond laser to generate a two-photon excitation pulse, and use an acousto-optic modulator to control the energy density of the two-photon excitation pulse in real time, so that the package structure generates a stable fluorescence response signal; control the confocal scanning unit to scan the package structure layer by layer along the Z-axis direction to obtain the three-dimensional distribution data of the fluorescence response signal; establish a stress distribution characteristic map of the package structure according to the three-dimensional distribution data, and mark the area where the stress gradient exceeds the preset threshold in the stress distribution characteristic map to generate a stress anomaly location map; calculate the stress gradient value for each area in the stress anomaly location map, and establish a stress gradient classification map;

[0067] The second unit is used to determine the detection area according to the stress gradient grading map, arrange the quantum nitrogen vacancy probe array in the detection area, and set the radio frequency pulse sequence parameters according to the spatial distribution of the quantum nitrogen vacancy probe array, so that the quantum nitrogen vacancy probe array can achieve the corresponding detection sensitivity in each detection area; collect the electromagnetic field response signal generated by the quantum nitrogen vacancy probe array under the action of stress, reconstruct the electromagnetic field distribution map of the packaging structure, and spatially correspond it with the stress gradient grading map to form a stress-electromagnetic field coupling model;

[0068] The third unit is used to establish a group of thermomechanical and electromagnetic coupling equations for the packaging structure based on the stress-electromagnetic field coupling model; calculate the evolution characteristics of the group of thermomechanical and electromagnetic coupling equations in each detection area, obtain the interaction relationship between each physical field, and calculate the failure probability of each detection area under the multi-field coupling according to the interaction relationship, and establish a stability evaluation system including mechanical failure risk and electromagnetic failure risk; determine the quality level of the packaging structure by comparing with the stability evaluation system, and generate a quality assessment report for the packaging structure.

[0069] According to a third aspect of the embodiments of the present invention,

[0070] An electronic device is provided, comprising:

[0071] processor;

[0072] a memory for storing processor-executable instructions;

[0073] The processor is configured to call the instructions stored in the memory to execute the aforementioned method.

[0074] A fourth aspect of the embodiments of the present invention is:

[0075] A computer-readable storage medium is provided, on which computer program instructions are stored. When the computer program instructions are executed by a processor, the aforementioned method is implemented.

[0076] In an embodiment of the present invention, two-photon excitation microscopy is combined with a quantum nitrogen vacancy probe array, and the combination realizes precise positioning and quantitative analysis of the stress distribution inside the packaging structure; two-photon excitation microscopy can obtain three-dimensional stress distribution characteristics, while the quantum nitrogen vacancy probe array can detect electromagnetic field changes at the nanoscale, thereby improving the accuracy and sensitivity of detection; by setting the detection area and radio frequency pulse sequence parameters according to the stress gradient grading diagram, the detection sensitivity of the quantum nitrogen vacancy probe can be optimized, so that it can achieve the best performance in different stress gradient areas, further improving the detection effect; a stress-electromagnetic field coupling model and a group of thermal electromagnetic coupling equations are established, which not only pay attention to stress distribution, but also consider the influence of electromagnetic field and temperature field. By establishing a coupling model and a group of equations, the failure mechanism of the packaging structure can be more comprehensively analyzed; the interaction relationship between each physical field is calculated, and by analyzing the interaction between each physical field, the failure risk of the packaging structure in the actual working environment can be more accurately predicted, avoiding the limitations brought by a single physical field analysis; a stability evaluation system including mechanical failure risk and electromagnetic failure risk is established, considering the electromagnetic failure risk, and establishing a more comprehensive stability evaluation system. BRIEF DESCRIPTION OF THE DRAWINGS

[0077] Figure 1It is a flow chart of an integrated circuit packaging quality assessment method based on nondestructive testing technology according to an embodiment of the present invention;

[0078] Figure 2 It is a schematic diagram of the structure of an integrated circuit packaging quality assessment system based on nondestructive testing technology according to an embodiment of the present invention. DETAILED DESCRIPTION

[0079] In order to make the purpose, technical solution and advantages of the embodiments of the present invention clearer, the technical solution in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0080] The technical solution of the present invention is described in detail with specific embodiments below. The following specific embodiments can be combined with each other, and the same or similar concepts or processes may not be described in detail in some embodiments.

[0081] Figure 1 FIG. 1 is a flow chart of an integrated circuit packaging quality assessment method based on nondestructive testing technology according to an embodiment of the present invention. Figure 1 As shown, the method includes:

[0082] The package structure of the integrated circuit to be inspected is placed on the sample stage of the two-photon excitation microscope system, a femtosecond laser is used to generate a two-photon excitation pulse, and the energy density of the two-photon excitation pulse is controlled in real time by an acousto-optic modulator, so that the package structure generates a stable fluorescence response signal; the confocal scanning unit is controlled to scan the package structure layer by layer along the Z-axis direction to obtain the three-dimensional distribution data of the fluorescence response signal; a stress distribution characteristic map of the package structure is established according to the three-dimensional distribution data, and the area where the stress gradient exceeds the preset threshold is marked in the stress distribution characteristic map to generate a stress anomaly location map; the stress gradient value is calculated for each area in the stress anomaly location map to establish a stress gradient grading map;

[0083] Determine a detection area according to the stress gradient grading map, arrange a quantum nitrogen vacancy probe array in the detection area, and set radio frequency pulse sequence parameters according to the spatial distribution of the quantum nitrogen vacancy probe array so that the quantum nitrogen vacancy probe array achieves a corresponding detection sensitivity in each detection area; collect the electromagnetic field response signal generated by the quantum nitrogen vacancy probe array under stress, reconstruct the electromagnetic field distribution map of the packaging structure, and spatially correspond it with the stress gradient grading map to form a stress-electromagnetic field coupling model;

[0084] Based on the stress-electromagnetic field coupling model, a group of thermomechanical and electromagnetic coupling equations for the packaging structure is established; the evolution characteristics of the group of thermomechanical and electromagnetic coupling equations in each detection area are calculated to obtain the interaction relationship between the physical fields; based on the interaction relationship, the failure probability of each detection area under the multi-field coupling is calculated, and a stability evaluation system including mechanical failure risk and electromagnetic failure risk is established; the quality grade of the packaging structure is determined by comparing with the stability evaluation system, and a quality assessment report of the packaging structure is generated.

[0085] In a specific implementation, an integrated circuit package structure sample to be inspected is prepared and fixed on a sample stage of a two-photon excitation microscope system.

[0086] A femtosecond laser is used as the excitation light source to generate two-photon excitation pulses. In order to ensure that the package structure produces a stable fluorescence response signal, the energy density of the excitation pulse needs to be controlled in real time. This can be achieved through an acousto-optic modulator, which can adjust the parameters of the acousto-optic modulator, such as the modulation depth or frequency, according to the real-time monitoring of the fluorescence signal intensity feedback, so as to finely control the excitation pulse energy.

[0087] The confocal scanning unit of the microscope is controlled to scan the package structure layer by layer along the Z-axis direction. The scanning step is set to, for example, 0.5 microns to ensure that sufficiently fine three-dimensional structural information is obtained. During the scanning process, the fluorescence response signal intensity of each position is recorded, and finally a complete data set containing three-dimensional spatial coordinates and corresponding fluorescence intensity is obtained.

[0088] The collected three-dimensional fluorescence data is imported into the data processing software, and a specific algorithm (such as finite element analysis) is used to establish a stress distribution characteristic map of the packaging structure. The characteristic map represents the magnitude of stress by the change of color or grayscale value. A stress gradient threshold is set, such as 10MPa / μm, and the area where the stress gradient exceeds the threshold in the characteristic map is marked to generate a stress anomaly location map.

[0089] The stress gradient value of each marked area in the stress anomaly location map is calculated and graded according to the size of the gradient value, for example, the gradient value is divided into three levels: low, medium and high, to generate a stress gradient grading map. Different levels are represented by different colors or gray values ​​for easy distinction.

[0090] According to the stress gradient classification map, the areas that need further fine detection are selected as detection areas. In these detection areas, quantum nitrogen vacancy probe arrays are arranged. The spatial distribution density of the probe array is adjusted according to the stress gradient level, for example, more dense probes are arranged in areas with high stress gradients to improve spatial resolution.

[0091] For each nitrogen vacancy probe in the detection area, the corresponding radio frequency pulse sequence parameters are set to achieve the best detection sensitivity. For example, the probe's response to the local electromagnetic field can be optimized by adjusting parameters such as the frequency, width and repetition rate of the radio frequency pulse.

[0092] Collect the electromagnetic field response signal generated by the probe array under stress. Since stress affects the spin state of nitrogen vacancy defects and thus changes their response to radio frequency pulses, electromagnetic field information can be obtained by measuring the fluorescence intensity or other optical properties of the probe.

[0093] The collected electromagnetic field response signals are processed and reconstructed to generate the electromagnetic field distribution map of the packaging structure. Then the distribution map is spatially matched with the previously obtained stress gradient classification map to establish a stress-electromagnetic field coupling model. This model reflects the spatial correlation between stress and electromagnetic field.

[0094] Based on the established stress-electromagnetic field coupling model, combined with the physical properties of the material and related physical laws, the thermal-electromagnetic coupling equations of the packaging structure are established. For example, the heat conduction equation, the electromagnetic field equation, and the stress-strain relationship can be considered.

[0095] The thermo-electromagnetic coupling equations are numerically solved in each detection area to analyze their evolution characteristics, such as temperature changes, changes in electromagnetic field intensity, and stress changes. By analyzing the interaction between these physical fields, such as temperature increases leading to stress changes, changes in electromagnetic field intensity affecting current distribution, etc., we can have a deeper understanding of the physical processes inside the packaging structure.

[0096] According to the multi-field coupling evolution characteristics of each detection area, the failure probability of each area under the multi-field coupling is calculated. For example, the possibility of mechanical failure or electromagnetic failure in each area can be evaluated based on indicators such as the strength limit of the material and the electromagnetic interference threshold.

[0097] Based on the failure probability of each detection area, a stability evaluation system including mechanical failure risk and electromagnetic failure risk is established. This system can classify the quality of the packaging structure according to the size of the failure probability, for example, into qualified, good, and excellent.

[0098] Finally, according to the evaluation results of the stability evaluation system, the quality level of the packaging structure is determined, and a quality evaluation report of the packaging structure is generated. The report can include detailed information such as the failure probability, stress distribution, and electromagnetic field distribution of each detection area, providing a reference for subsequent packaging process improvements.

[0099] For example, in a specific case, the stress gradient of a certain area in the chip packaging structure was observed to be 15MPa / μm through a two-photon microscope, which is a high stress gradient area. A high-density nitrogen vacancy probe array is arranged in this area, and the electromagnetic field response signal of the probe is collected by optimizing the radio frequency pulse parameters. Through data processing and reconstruction, it is found that the electromagnetic field intensity in this area is high, and there is an obvious spatial correlation with the stress distribution. By solving the thermo-electromagnetic coupling equations, it is found that the temperature in this area will increase under working conditions, further aggravating the stress concentration and resulting in a high probability of mechanical failure in this area. Finally, according to the stability evaluation system, the quality grade of the packaging structure is rated as "good".

[0100] In this embodiment, through multi-field coupling analysis, the influence of various factors on the stability of the packaging structure can be considered more comprehensively, so as to more accurately evaluate its quality level; by calculating the failure probability of each detection area, the potential risks of the packaging structure can be quantitatively evaluated, providing more precise guidance for subsequent process improvements; by analyzing the multi-field coupling mechanism, the weak links in the packaging process can be identified and optimized in a targeted manner, so as to improve the reliability and stability of the packaging structure.

[0101] In an optional implementation, establishing a stress distribution characteristic diagram of a packaging structure according to the three-dimensional distribution data includes:

[0102] Inputting the three-dimensional distribution data of the fluorescence response signal into an optical system model, wherein the optical system model includes a lateral resolution parameter determined by the numerical aperture of the objective lens and a longitudinal resolution parameter determined by the excitation wavelength;

[0103] Generate a Gaussian point spread function according to the optical system model, perform a deconvolution operation on the three-dimensional distribution data and the Gaussian point spread function, and perform Wiener filtering on the result of the deconvolution operation to obtain three-dimensional fluorescence distribution data after eliminating the influence of the optical system;

[0104] Measure the photoelastic coefficient tensor and refractive index anisotropy coefficient of the material corresponding to the packaging structure, and establish the photoelastic effect response equation of the material stress field and refractive index change;

[0105] Establishing a nonlinear relationship equation between the refractive index and the fluorescence intensity during two-photon excitation, combining the nonlinear relationship equation with the photoelastic effect response equation to obtain a quantitative mapping model between the stress field and the fluorescence intensity;

[0106] Selecting an initial stress field, calculating a theoretical fluorescence intensity according to the quantitative mapping model, comparing the theoretical fluorescence intensity with the measured fluorescence intensity of the three-dimensional fluorescence distribution data, and calculating a mean square error as an optimization objective function; performing a gradient calculation, calculating a gradient vector of the optimization objective function, determining a conjugate search direction, performing a conjugate search, and iteratively searching for an optimal step size;

[0107] The stress field distribution is updated using the optimal step size. When the convergence value of the optimization objective function is greater than a preset convergence threshold, the gradient calculation and conjugate search are continued. When the convergence value of the optimization objective function is less than or equal to the convergence threshold, the current stress field distribution is used as a stress distribution characteristic diagram of the packaging structure.

[0108] The stress distribution inside the package structure is crucial to its performance and reliability. In order to accurately measure and reconstruct the three-dimensional stress distribution inside the package structure, a measurement method based on two-photon fluorescence and photoelastic effect is proposed. This method uses the refractive index change of the material under stress and combines the fluorescence intensity of two-photon excitation to quantitatively characterize the stress.

[0109] In a specific embodiment, a confocal microscope is used to collect the three-dimensional fluorescence response signal inside the package structure. Taking a cubic package structure with a side length of 100 microns as an example, the scanning step size is set to 1 micron, and a 100×100×100 three-dimensional data matrix can be obtained, and each element in the matrix represents the fluorescence intensity at that position.

[0110] Then, considering the influence of the optical system on the measurement results, the collected fluorescence signal needs to be corrected. It is known that the numerical aperture of the objective lens used is 0.8 and the excitation wavelength is 800 nanometers. According to these parameters, the lateral resolution of the optical system can be determined to be 220 nanometers and the longitudinal resolution is 500 nanometers. A Gaussian point spread function is constructed using these resolution parameters. The collected three-dimensional fluorescence data is deconvolved with the point spread function, and the deconvolution result is denoised using the Wiener filter to obtain the three-dimensional fluorescence distribution data after eliminating the influence of the optical system. The parameters of the Wiener filter are set according to the noise level. Assuming the signal-to-noise ratio is 10, the corresponding Wiener filter parameters are set to 0.1.

[0111] Next, the stress-optical properties of the material need to be determined. The photoelastic coefficient tensor and refractive index anisotropy coefficient of the corresponding material of the packaging structure are measured experimentally. Assuming that the photoelastic coefficient of the material is C 11 =0.2, C 12 =0.1, C 44 =0.15, and the refractive index anisotropy coefficient is 0.05. These parameters are used to establish the photoelastic effect response equation of material stress field and refractive index change.

[0112] At the same time, the nonlinear relationship equation between the refractive index and the fluorescence intensity during two-photon excitation is established. Assuming that the fluorescence intensity is proportional to the square of the refractive index and the proportionality coefficient is 0.5, the corresponding equation can be established. Combining this equation with the photoelastic effect response equation, a quantitative mapping model between stress field and fluorescence intensity is obtained.

[0113] In order to reconstruct the stress distribution inside the packaging structure, an iterative optimization method is used. First, an initial stress field is set, for example, assuming that the stress is uniformly distributed and the value is 10MPa. According to the established quantitative mapping model, the theoretical fluorescence intensity distribution corresponding to the stress field is calculated. The calculated theoretical fluorescence intensity is compared with the experimentally measured three-dimensional fluorescence distribution data, and the mean square error between the two is calculated as the optimization objective function.

[0114] Then, perform gradient calculation, calculate the gradient vector of the optimization objective function, and determine the conjugate search direction. Perform conjugate search along this direction and iterate to find the optimal step size. Assuming the initial step size is 1, the optimal step size is determined to be 0.5 through linear search.

[0115] Update the stress field distribution using the optimal step size. Repeat the steps of gradient calculation, conjugate search, and stress field update until the convergence value of the optimization objective function is less than the preset convergence threshold. Assuming the convergence threshold is 0.01, when the optimization objective function is less than the threshold, stop the iteration and use the current stress field distribution as the final packaging structure stress distribution characteristic map.

[0116] In this embodiment, the three-dimensional stress distribution inside the packaging structure can be accurately reconstructed, overcoming the limitation that the traditional method can only measure surface stress; by using optical measurement means, there is no need to physically destroy the packaging structure, and the integrity of the sample can be maintained; a quantitative mapping model of stress field and fluorescence intensity is established, which can achieve quantitative measurement of stress field rather than just qualitative analysis.

[0117] In an optional implementation, marking an area in the stress distribution characteristic map where the stress gradient exceeds a preset threshold, and generating a stress anomaly location map includes:

[0118] Using the Sobel operator to calculate the gradient of the stress distribution characteristic diagram in the horizontal direction to obtain a first gradient component matrix, and using the Sobel operator to calculate the gradient of the stress distribution characteristic diagram in the vertical direction to obtain a second gradient component matrix;

[0119] Calculating the square sum of the first gradient component matrix and the second gradient component matrix, and performing a square root operation to obtain a gradient matrix of the stress distribution characteristic diagram; performing a power-of-2 multi-scale decomposition on the gradient matrix, wherein the multi-scale decomposition ranges from 2×2 to 16×16, and obtaining sub-matrices of multiple scales;

[0120] Calculating the grayscale mean and grayscale standard deviation of each local area in the sub-matrices of the multiple scales respectively, determining the segmentation threshold of the corresponding local area according to the grayscale mean and the grayscale standard deviation, and determining the adaptive segmentation threshold;

[0121] Based on an adaptive segmentation threshold, the sub-matrices of the multiple scales are binarized respectively to obtain binarized sub-matrices of the multiple scales; the grayscale variance of each binarized sub-matrix of the multiple scales is calculated respectively, and the grayscale variance is normalized to obtain the fusion weight of each scale; the binarized sub-matrices of the multiple scales are multiplied and superimposed with the corresponding fusion weights respectively to obtain a fused segmentation matrix; the fused segmentation matrix is ​​binarized to obtain a mask matrix, and in the mask matrix, the stress abnormality area is determined based on the area with a value of 1; the stress abnormality area is marked in the stress distribution characteristic map to generate a stress abnormality location map.

[0122] A stress distribution feature map is obtained. The feature map represents the stress magnitude with pixel grayscale values. For example, a 512×512 pixel grayscale image has a grayscale value range of 0 to 255. A larger value indicates a greater stress.

[0123] First, the Sobel operator is used to calculate the horizontal gradient of the stress distribution feature map. The Sobel operator is a 3×3 matrix used to calculate the horizontal gradient of each pixel in the image. For example, for a pixel with coordinates (x, y) in the feature map, its horizontal gradient can be obtained by performing a specific calculation with the grayscale values ​​of the 8 pixels around it. The horizontal gradient values ​​of all pixels are combined into the first gradient component matrix.

[0124] Then, the same Sobel operator is used to calculate the vertical gradient of the stress distribution characteristic map to obtain the second gradient component matrix. The calculation method is similar to the horizontal gradient calculation, but a different Sobel operator is used.

[0125] Next, the stress gradient amplitude matrix is ​​calculated. For each pixel in the feature map, the square root of the sum of the squares of the first gradient component value and the second gradient component value corresponding to the pixel is calculated to obtain the gradient amplitude of the pixel. The gradient amplitudes of all pixels form the stress gradient amplitude matrix.

[0126] The gradient magnitude matrix is ​​decomposed into a power of 2 multi-scale, ranging from 2 × 2 to 16 × 16. For example, a 512 × 512 gradient magnitude matrix is ​​decomposed into multiple 2 × 2, 4 × 4, 8 × 8, and 16 × 16 sub-matrices.

[0127] For each scale submatrix, the grayscale mean and grayscale standard deviation of each local region are calculated. For example, for a 4×4 submatrix, it can be divided into four 2×2 local regions, and the grayscale mean and grayscale standard deviation of each 2×2 region are calculated respectively.

[0128] According to the grayscale mean and grayscale standard deviation of each local area, the adaptive segmentation threshold of the area is determined. For example, the grayscale mean of the local area plus a certain multiple of the grayscale standard deviation can be used as the segmentation threshold.

[0129] The calculated adaptive segmentation threshold is used to binarize the submatrix of each scale. For example, if the gradient magnitude of a pixel is greater than or equal to the segmentation threshold of the region, its value is set to 1, otherwise it is set to 0. In this way, binarized submatrices of multiple scales are obtained.

[0130] Calculate the grayscale variance of the binary submatrix at each scale. Then, normalize the grayscale variance at each scale to obtain the fusion weight at each scale. For example, the grayscale variance at each scale can be divided by the sum of the grayscale variances at all scales to obtain the normalized fusion weight.

[0131] The binarized submatrix of each scale is multiplied by its corresponding fusion weight, and then the results of all scales are superimposed to obtain the fused segmentation matrix.

[0132] The fused segmentation matrix is ​​binarized to obtain a mask matrix. In the mask matrix, the area with a value of 1 is the stress abnormality area.

[0133] Finally, the stress anomaly area is marked on the original stress distribution feature map. For example, the pixels in the stress anomaly area can be marked with a specific color or symbol to generate a stress anomaly location map. Assume a 4×4 gradient amplitude matrix as follows:

[0134] ;

[0135] After 2×2 scale decomposition, four 2×2 sub-matrices are obtained. Assume that one of the sub-matrices is:

[0136] ;

[0137] The grayscale mean of this submatrix is ​​19.25, and the grayscale standard deviation is 9.43. Assume that the segmentation threshold is the mean plus the standard deviation, which is 28.68. Then the submatrix after binarization is:

[0138] ;

[0139] In this embodiment, through multi-scale decomposition and adaptive threshold segmentation, stress anomaly areas of different sizes can be identified more accurately, avoiding the limitations of the single-scale method, thereby improving the positioning accuracy; fusing the segmentation results of multiple scales can effectively reduce the impact of noise on the segmentation results and improve the robustness of the algorithm; the adaptive threshold segmentation method can automatically adjust the segmentation threshold according to the stress distribution in different areas without manual intervention, thereby improving the adaptability of the algorithm.

[0140] In an optional embodiment, determining a detection area according to the stress gradient grading map, arranging a quantum nitrogen vacancy probe array in the detection area, and setting radio frequency pulse sequence parameters according to the spatial distribution of the quantum nitrogen vacancy probe array so that the quantum nitrogen vacancy probe array achieves a corresponding detection sensitivity in each detection area includes:

[0141] Constructing a quantum state evolution measurement system, applying an initialization radio frequency pulse to the quantum state nitrogen vacancy probe array, and initializing the quantum state nitrogen vacancy probe array to a spin ground state;

[0142] Using a quantum state readout sequence to obtain the spin state information of the quantum state nitrogen vacancy probe array, record the spin state evolution time, measure the quantum coherence time of the quantum state nitrogen vacancy probe array in each detection area, and establish a quantum state response model;

[0143] Using a quantum state subspace decomposition method, the quantum state response model is characteristically decomposed to obtain a quantum state characteristic vector of each detection area; a quantum decoherence noise estimator is constructed based on the quantum state characteristic vector, the influence of environmental noise on the quantum state nitrogen vacancy probe of each detection area is measured, and the quantum state decoherence rate of the detection area is calculated;

[0144] According to the quantum state decoherence rate, the minimum detection sensitivity corresponding to each detection area is determined, and a multi-channel quantum state manipulation timing scheme is established to achieve synchronous quantum state manipulation between detection areas;

[0145] For each detection area, a quantum state characteristic transformation matrix is ​​constructed to extract quantum coherence information of the quantum state nitrogen vacancy probe to obtain the quantum state response characteristics of the detection area;

[0146] Inputting the quantum state response characteristics into the quantum information fusion framework, and calculating the sensitivity parameters of each detection area based on the quantum state evolution analysis algorithm;

[0147] Based on the sensitivity parameters, the radio frequency pulse phase accumulation time corresponding to each detection area is calculated, the radio frequency pulse repetition period is determined, and radio frequency pulse sequence parameters that meet the detection sensitivity requirements of each detection area are generated.

[0148] The quantum coherence time refers to the length of time a quantum system can maintain coherence. Coherence refers to the phase relationship between states in a quantum system, which is crucial to the effectiveness of quantum computing and quantum sensing. Within the quantum coherence time of the system, the quantum state can maintain a stable quantum superposition state and will not be significantly affected by environmental noise or other interference. Once the quantum coherence time is exceeded, the coherence of the quantum state will gradually disappear, manifesting as quantum decoherence. The longer the quantum coherence time, the better the information transmission and processing effect of the quantum state during the evolution process. Therefore, the quantum coherence time is an important performance indicator in quantum state measurement.

[0149] The quantum state decoherence rate refers to the rate at which the coherence of a quantum system disappears, and measures the extent to which the quantum system loses coherence due to environmental influences. The decoherence rate is related to factors such as environmental noise and the system's spin interactions. A higher decoherence rate means that the coherence of the system will decay significantly in a shorter period of time, resulting in the loss of quantum state information. Therefore, in quantum information processing applications, a lower decoherence rate is usually more ideal. By calculating the quantum state decoherence rate, the impact of the environment on the quantum state nitrogen vacancy probe can be quantified, thereby evaluating and improving the detection sensitivity of the detection area.

[0150] First, prepare a diamond sample containing nitrogen vacancy color centers. As quantum probes, the quantum state of nitrogen vacancy color centers is extremely sensitive to stress changes in the surrounding environment. Assume that the stress distribution on the surface of the sample shows a gradient change, for example, the stress gradually increases from left to right.

[0151] Next, the surface of the diamond sample is divided into different detection areas according to the stress gradient map. For example, the sample is divided into three areas: low stress area, medium stress area and high stress area. In each detection area, a probe array consisting of nitrogen vacancy color centers is arranged. The spatial distribution of nitrogen vacancy color centers in the array can be adjusted according to actual needs, such as uniform distribution or non-uniform distribution according to stress gradient changes.

[0152] Then, a quantum state evolution measurement system was built. The system includes a microwave source and radio frequency coil for controlling the spin state of nitrogen vacancy color centers, as well as an optical system and detector for collecting the fluorescence signal of nitrogen vacancy color centers. Using this system, an initialization radio frequency pulse is applied to the nitrogen vacancy probe array to initialize the spin state of all nitrogen vacancy color centers to the ground state.

[0153] Use the quantum state readout sequence to obtain the spin state information of the nitrogen vacancy probe array. For example, the Ramsey pulse sequence is used to measure the changes in the fluorescence signal of the nitrogen vacancy color center to obtain the information of its spin state at different evolution times. The spin state evolution time is recorded, and the quantum coherence time of the nitrogen vacancy probe array in each detection area is measured, that is, the time the spin state remains in the quantum superposition state. By analyzing the law of the change of quantum coherence time with the evolution time, a quantum state response model is established, which describes the evolution law of the spin state of the nitrogen vacancy color center under the action of external stress.

[0154] The quantum state subspace decomposition method is used to perform feature decomposition on the established quantum state response model. For example, the principal component analysis method is used to decompose the quantum state response model into multiple mutually orthogonal quantum state feature vectors. These feature vectors represent the main evolution modes of the nitrogen vacancy color center spin state in different detection regions.

[0155] Based on the acquired quantum state eigenvector, a quantum decoherence noise estimator is constructed. By measuring the influence of environmental noise (such as temperature fluctuations and magnetic field noise) on the spin state of the nitrogen vacancy color center in each detection area, the quantum state decoherence rate of the detection area is calculated. The decoherence rate reflects the speed at which the nitrogen vacancy color center spin state degenerates from the quantum superposition state to the classical state due to environmental noise.

[0156] According to the calculated quantum state decoherence rate, the minimum detection sensitivity corresponding to each detection area is determined. For example, the higher the decoherence rate, the lower the minimum detection sensitivity. A multi-channel quantum state manipulation timing scheme is established to achieve synchronous manipulation of the spin states of nitrogen vacancy color centers in different detection areas.

[0157] For each detection area, a quantum state characteristic transformation matrix is ​​constructed. This matrix is ​​used to transform the quantum state information of the nitrogen vacancy color center from the original space to the characteristic space, extract the quantum coherence information of the nitrogen vacancy color center, and obtain the quantum state response characteristics of the detection area.

[0158] The extracted quantum state response characteristics are input into the quantum information fusion framework. Based on the quantum state evolution analysis algorithm, such as the Bayesian estimation method, the sensitivity parameter of each detection area is calculated. This parameter reflects the sensitivity of the nitrogen vacancy color center to the stress change in the area.

[0159] Based on the calculated sensitivity parameters, the RF pulse phase accumulation time corresponding to each detection area is calculated. The RF pulse repetition period is determined, and RF pulse sequence parameters that meet the detection sensitivity requirements of each detection area are generated. For example, for areas with high sensitivity requirements, a shorter RF pulse repetition period and a more accurate phase accumulation time are required.

[0160] Assuming that the sample is divided into three regions, the measured quantum state decoherence rates of the three regions are 0.1ms⁻¹, 0.2ms⁻¹ and 0.3ms⁻¹ respectively. Based on these decoherence rates, the minimum detection sensitivities of the three regions are calculated to be 10nN / μm², 5nN / μm² and 3nN / μm² respectively. The final generated RF pulse sequence parameters can ensure that each detection area reaches its corresponding minimum detection sensitivity.

[0161] In this embodiment, by dividing the detection area according to the stress gradient and optimizing the RF pulse parameters for each area, refined measurement of different stress areas can be achieved, thereby improving the spatial resolution of stress measurement; by considering the influence of environmental noise on the quantum state and adjusting the RF pulse parameters according to the quantum state decoherence rate, the influence of noise can be effectively suppressed, thereby improving the sensitivity of stress measurement; it can adapt to the measurement requirements of different stress gradients and optimize parameters for different areas, thereby achieving accurate measurement of complex stress distribution.

[0162] In an optional embodiment, collecting the electromagnetic field response signal generated by the quantum nitrogen vacancy probe array under stress, reconstructing the electromagnetic field distribution map of the packaging structure, and spatially corresponding it with the stress gradient classification map to form a stress-electromagnetic field coupling model includes:

[0163] The electromagnetic field response signal of the quantum nitrogen vacancy probe array is collected in each detection area to construct an electromagnetic field response time series; the electromagnetic field response time series is discretely sampled, a signal analytical expression is obtained by Hilbert transform, and the instantaneous value of the electromagnetic field intensity in each detection area is extracted; the instantaneous value of the electromagnetic field intensity is constructed as an electromagnetic field response matrix, and the main eigenvectors of the electromagnetic field distribution are extracted by principal component analysis method;

[0164] Using continuous wavelet transform to perform time-frequency decomposition on the electromagnetic field response matrix to obtain electromagnetic field spectrum components at different scales; performing adaptive threshold filtering on the electromagnetic field spectrum components to eliminate measurement noise and obtain frequency domain characteristics of the electromagnetic field signal;

[0165] A compressed sensing reconstruction algorithm is used to sparsely reconstruct the frequency domain features, a mathematical model of electromagnetic field spatial distribution is established, an electromagnetic field inversion equation is constructed based on the mathematical model, and a conjugate gradient iteration method is used to solve the electromagnetic field distribution to obtain an electromagnetic field distribution map of the packaging structure; the electromagnetic field distribution map is affinely transformed and registered with the stress gradient classification map to establish a spatial coordinate correspondence relationship;

[0166] Based on the spatial coordinate correspondence, a kernel function regression method is used to extract the mapping relationship between the electromagnetic field intensity and the stress gradient, and the electromagnetic field-stress response data is obtained; based on the electromagnetic field-stress response data, a nonlinear coupling equation group is established, and the coupling coefficient is solved by the Newton iteration method;

[0167] Substituting the coupling coefficient into the Maxwell stress-electromagnetic field tensor equation, a stress-electromagnetic field coupling model of the packaging structure is formed.

[0168] The compressed sensing is a signal processing technology that can reconstruct sparse or compressible signals with fewer sampling points than required by the traditional sampling theorem. The compressed sensing reconstruction algorithm uses the sparsity of the signal to obtain approximate sampling data of the high-dimensional signal through a small number of linear measurements, and then restores the original signal through an optimization algorithm (such as minimizing a certain loss function). The compressed sensing reconstruction algorithm is used to sparsely reconstruct the frequency domain characteristics of the electromagnetic field signal, thereby effectively restoring the spatial distribution characteristics of the electromagnetic field, reducing the amount of calculation, and improving the measurement accuracy.

[0169] Electromagnetic field inversion in the electromagnetic field inversion equation refers to the process of deriving or reconstructing the unknown electromagnetic field distribution in the physical system through the electromagnetic field response. The electromagnetic field inversion equation is a mathematical model that establishes the relationship between the observed data and the electromagnetic field, and infers the spatial distribution of the electromagnetic field through the known response signal (such as the instantaneous value of the electromagnetic field intensity). This process usually requires the use of numerical methods (such as the finite element method or the conjugate gradient method) to solve. By establishing the inversion equation of the electromagnetic field, the measured data is used to infer the electromagnetic field distribution inside the packaging structure.

[0170] The affine transformation registration specifically refers to an image or data set alignment technology, the purpose of which is to align or match point sets from different sources in two spaces (such as electromagnetic field distribution maps and stress gradient maps) to a common coordinate system. Affine transformation is a linear transformation that includes not only translation, rotation, and scaling, but may also include shear transformations, which can be used to perform geometric transformations on data or images. Specifically, affine transformation registration will accurately align the electromagnetic field distribution map and the stress gradient map based on the correspondence between spatial coordinates, so that the two can correctly correspond in spatial coordinates, providing an accurate basis for subsequent analysis and coupled modeling.

[0171] First, a package structure sample containing a quantum nitrogen vacancy probe array is prepared. The sample can be a chip, in which the nitrogen vacancy probe array is embedded in a specific area of ​​the chip. A known stress gradient is applied to the sample, for example, by bending or stretching the sample. The stress gradient can be set in stages, for example, divided into several levels, so as to subsequently analyze the relationship between stress and electromagnetic field response.

[0172] Next, a specific device, such as a light detection magnetic resonance device, is used to collect the electromagnetic field response signal of the nitrogen vacancy probe array in each detection area of ​​the sample. Each detection area corresponds to a specific stress level in the stress gradient grading diagram. The electromagnetic field response signal is reflected in the spin state change of the nitrogen vacancy probe, which can be converted into an electromagnetic field intensity value that changes with time, thereby constructing an electromagnetic field response time series.

[0173] In order to facilitate subsequent processing, the collected electromagnetic field response time series is discretely sampled. For example, the electromagnetic field intensity value is recorded at a certain sampling frequency. Then, the discretely sampled signal is Hilbert transformed to obtain the analytical expression of the signal. From the analytical expression of the signal, the instantaneous value of the electromagnetic field intensity in each detection area can be extracted.

[0174] The instantaneous values ​​of the electromagnetic field intensity in all detection areas are combined into a matrix, called the electromagnetic field response matrix. The rows of the matrix represent different detection areas, and the columns represent different sampling time points. The electromagnetic field response matrix is ​​subjected to principal component analysis to extract the main eigenvectors of the electromagnetic field distribution. These eigenvectors represent the main change patterns of the electromagnetic field distribution.

[0175] In order to analyze the frequency domain characteristics of the electromagnetic field signal, the electromagnetic field response matrix is ​​decomposed in time and frequency using continuous wavelet transform. This will decompose the electromagnetic field signal into spectral components at different scales. The spectral components at different scales are adaptively threshold filtered to eliminate measurement noise. The filtered spectral components represent the true frequency domain characteristics of the electromagnetic field signal.

[0176] In order to reconstruct the electromagnetic field distribution map of the packaging structure, the compressed sensing reconstruction algorithm is used to sparsely reconstruct the filtered frequency domain features. This will establish a mathematical model of the electromagnetic field spatial distribution. Based on this mathematical model, the electromagnetic field inversion equation is constructed, and the conjugate gradient iteration method is used to solve the electromagnetic field distribution. Finally, the electromagnetic field distribution map of the packaging structure is obtained.

[0177] In order to establish the spatial correspondence between stress and electromagnetic field, the obtained electromagnetic field distribution map is registered with the pre-set stress gradient classification map by affine transformation, which will establish the spatial coordinate correspondence between the two images.

[0178] Based on the established spatial coordinate correspondence, the kernel function regression method is used to extract the mapping relationship between electromagnetic field intensity and stress gradient. This will obtain a series of electromagnetic field-stress response data.

[0179] Based on the obtained electromagnetic field-stress response data, a set of nonlinear coupling equations is established, and the coupling coefficients are solved using the Newton iteration method. These coupling coefficients describe the interaction between stress and electromagnetic field.

[0180] Finally, the solved coupling coefficient is substituted into the Maxwell stress-electromagnetic field tensor equation to form the stress-electromagnetic field coupling model of the packaging structure.

[0181] For example, suppose a chip sample is subjected to three levels of stress gradient. At each stress level, the electromagnetic field intensity values ​​at 10 time points are collected from three detection areas. Through the above steps, a model describing the coupling relationship between stress and electromagnetic field in the chip sample can be finally obtained.

[0182] In this embodiment, by utilizing the sensitivity of the quantum nitrogen vacancy probe array to the electromagnetic field, accurate measurement of tiny stress changes can be achieved; by reconstructing the electromagnetic field distribution map of the packaging structure, the distribution of the stress field can be intuitively displayed; by constructing a stress-electromagnetic field coupling model, the interaction relationship between stress and electromagnetic field can be quantitatively described, providing a theoretical basis for a deep understanding of the mechanical properties of the packaging structure.

[0183] In an optional implementation, a group of thermomechanical and electromagnetic coupling equations for the packaging structure is established based on the stress-electromagnetic field coupling model; the evolution characteristics of the group of thermomechanical and electromagnetic coupling equations in each detection area are calculated to obtain the interaction relationship between the physical fields, and the failure probability of each detection area under the multi-field coupling is calculated according to the interaction relationship, and a stability evaluation system including mechanical failure risk and electromagnetic failure risk is established; the quality level of the packaging structure is determined by comparing with the stability evaluation system, and the quality assessment report of the packaging structure is generated, including:

[0184] The temperature field is introduced into the stress-electromagnetic field coupling model to construct a stress distribution function including a thermal expansion coefficient; a corresponding relationship between temperature and conductivity is established based on the stress distribution function, and an exponential function including a parameter of the ratio of activation energy to the Boltzmann constant is used to determine the change characteristics of conductivity with temperature;

[0185] The stress distribution function is combined with a predetermined electromagnetic field basic equation group, a trial function is selected to construct a weight function, a discretization equation is established by piecewise linear interpolation, and a thermal electromagnetic coupling equation group of the packaging structure is obtained by minimizing the residual error;

[0186] The solution time interval is divided into multiple subintervals, and four iterative calculations including midpoint estimation and full interval correction are performed in each subinterval. The calculation step size is dynamically adjusted based on the local error of the subinterval to obtain the numerical solution of the thermo-electromagnetic coupling equations.

[0187] Performing signal decomposition operation on the numerical solution to extract characteristic mode functions that characterize the inherent change characteristics of the physical field of the packaging structure; calculating the ratio of the joint probability density to the edge probability density of the characteristic mode functions to obtain the interaction strength between the physical fields in the packaging structure;

[0188] By calculating the predicted force of the known physical field on the target physical field, a quantitative characterization model of the interaction between the physical fields of the packaging structure is established, and based on the quantitative characterization model, a failure judgment function including the physical field intensity threshold and the load response is constructed;

[0189] Based on the expectation value minimization criterion, an optimized probability density function is constructed, and multiple random sampling calculations are performed according to the failure judgment function to obtain the failure probability of the packaging structure;

[0190] Establishing a stability evaluation index system based on the failure probability of the packaging structure, and calculating a membership relationship matrix of each stability evaluation index to the quality grade;

[0191] The weight coefficient vector of the stability evaluation index is calculated by pairwise comparison, and the weight coefficient vector is multiplied by the membership matrix to obtain the stability evaluation result of the packaging structure.

[0192] The activation energy specifically refers to the energy barrier that a substance needs to overcome to change from one state to another. The activation energy is related to the property of conductivity changing with temperature. Conductivity usually increases with increasing temperature, especially in semiconductor or conductor materials, where electrons or ions are more likely to "transition" to a conductive state at high temperatures. The minimum energy required for this transition is the activation energy. In the conductivity model, the ratio of the activation energy to the Boltzmann constant is used to construct an exponential function relationship between conductivity and temperature to characterize the effect of temperature on conductivity.

[0193] The joint probability density specifically refers to the probability density that describes two or more variables (such as different physical field variables) taking certain specific values ​​at the same time. In the stress-electromagnetic coupling model, there is interaction between different physical fields (such as stress field and electromagnetic field), so the joint probability density is needed to characterize the distribution of these fields in the packaging structure. The joint probability density function can reflect the interdependence between multiple physical field variables and their common distribution characteristics.

[0194] The marginal probability density specifically refers to the integral of the joint probability density over one or more variables, which is used to represent the probability distribution of one of the variables while ignoring the influence of other variables. In this model, the marginal probability density is used to characterize the independent distribution characteristics of a single physical field (such as a single temperature field or a single stress field). The ratio of the joint probability density to the marginal probability density can be used to evaluate the interaction strength between the physical fields, that is, to determine whether the two physical fields are independent of each other or have a strong coupling relationship.

[0195] In a specific implementation, a stress-electromagnetic field-temperature field coupling model is established. The temperature field is introduced into the existing stress-electromagnetic field coupling model, and the influence of the thermal expansion coefficient is added to the stress distribution function. For example, assuming that the thermal expansion coefficient of the material is 2.5 × 10⁻ 6 / °C, the temperature change is 50°C, then the corresponding thermal stress term is added to the stress distribution function. At the same time, the corresponding relationship between temperature and conductivity is established. An exponential function is used to describe the change characteristics of conductivity with temperature. The exponential function contains the ratio of activation energy to Boltzmann constant. For example, assuming that the activation energy is 0.5eV and the Boltzmann constant is 8.617 × 10⁻ 5 eV / K, the specific functional relationship of conductivity changing with temperature can be determined.

[0196] Next, establish the thermo-electromagnetic coupling equations of the packaging structure. Combine the stress distribution function with the predetermined electromagnetic field basic equations, and select a suitable trial function to construct the weight function. Use the piecewise linear interpolation method to discretize the continuous physical field, establish the discretized equations, and obtain the thermo-electromagnetic coupling equations of the packaging structure by minimizing the residual. For example, divide the packaging structure into 100 grids, approximate the physical quantities in each grid using linear interpolation, and then solve the discretized equations using the least squares method.

[0197] Then, solve the thermo-electromagnetic coupling equations. Divide the solution time interval into multiple subintervals. For example, divide the time interval from 0 to 1 second into 100 subintervals. In each subinterval, use a four-iteration calculation method including midpoint estimation and full-interval correction to solve the equations. For example, use the fourth-order Runge-Kutta method for numerical calculations. At the same time, dynamically adjust the calculation step size according to the local error of the subinterval to ensure the calculation accuracy. For example, set the error tolerance to 10⁻ 6 , when the local error exceeds the tolerance, reduce the step size; when the local error is less than the tolerance, increase the step size.

[0198] After obtaining the numerical solution, a signal decomposition operation is performed to extract the characteristic mode function that characterizes the inherent change characteristics of the physical field of the packaging structure. For example, the empirical mode decomposition method is used to extract the characteristic mode function. Then, the ratio of the joint probability density and the edge probability density of the characteristic mode function is calculated to obtain the interaction strength between the physical fields in the packaging structure. For example, the interaction strength between the temperature field and the electromagnetic field is calculated to be 0.8, indicating that there is a strong coupling relationship between the two.

[0199] Next, a quantitative characterization model of the interaction between the physical fields of the packaging structure is established. By calculating the predicted power of the known physical field on the target physical field, for example, calculating the predicted power of the temperature field on the stress field. Based on this quantitative characterization model, a failure judgment function including the physical field strength threshold and load response is constructed. For example, the stress threshold is set to 100MPa, and when the stress exceeds this threshold, it is judged as a failure.

[0200] Then, the failure probability of the package structure is calculated. Based on the expected value minimization criterion, an optimized probability density function is constructed. Multiple random sampling calculations are performed according to the failure judgment function, for example, 10,000 Monte Carlo simulations are performed to obtain the failure probability of the package structure. For example, the failure probability of the package structure is calculated to be 0.01.

[0201] Finally, the stability of the package structure is evaluated. A stability evaluation index system is established based on the failure probability. For example, a failure probability less than 0.01 is defined as a high quality grade, 0.01 to 0.1 is defined as a medium quality grade, and greater than 0.1 is defined as a low quality grade. The membership matrix of each stability evaluation index to the quality grade is calculated. For example, the membership of the high quality grade to the failure probability less than 0.01 is 1, the membership to the failure probability between 0.01 and 0.1 is 0.5, and the membership to the failure probability greater than 0.1 is 0. The weight coefficient vector of the stability evaluation index is calculated by pairwise comparison. For example, the weight coefficient of the failure probability is 0.8, and the weight coefficients of other indicators are 0.2. The weight coefficient vector is multiplied by the membership matrix to obtain the stability evaluation result of the package structure. For example, the quality grade of the package structure finally calculated is high quality. A quality evaluation report of the package structure is generated, which contains a detailed evaluation process and results.

[0202] In this embodiment, by considering the coupling effect of multiple physical fields, the failure risk of the packaging structure can be more accurately evaluated, thereby guiding the design and optimization of the packaging structure and improving its reliability; the performance prediction of the packaging structure under different working conditions can be provided, providing a theoretical basis for the design of the packaging structure, thereby optimizing the design of the packaging structure so that it can meet the performance requirements while having lower cost and smaller size; the use of numerical calculation methods can quickly and accurately evaluate the stability of the packaging structure, thereby improving the efficiency of quality assessment and shortening the product development cycle.

[0203] In an optional embodiment, the solution time interval is divided into multiple subintervals, four iterative calculations including midpoint estimation and full interval correction are performed in each subinterval, and the calculation step size is dynamically adjusted based on the local error of the subinterval to obtain the numerical solution of the thermo-electromagnetic coupling equations, including:

[0204] Divide the solution time interval into a plurality of subintervals, set an initial time step at the initial time of the subinterval; at the initial time of the subinterval, use an explicit Euler format to calculate the initial value of the physical field of the thermo-electromagnetic coupling equation group;

[0205] Calculate the predicted value of the physical field at the intermediate time, the formula is as follows:

[0206] ;

[0207] in, t i Indicates the current moment, t i+1 / 2 represents the intermediate moment of advancing half a step, u * (·) represents the predicted value of the physical field, u (·) represents the physical field value, Δ t represents the time step, f (·) represents the value of the physical field function;

[0208] Calculate the predicted value of the physical field at the end time, the formula is as follows:

[0209] ;

[0210] in, u * ( t i+1 ) represents the predicted value of the physical field at the termination time;

[0211] Calculate the physical field correction value at the end time, the formula is as follows:

[0212] ;

[0213] in, u ' (·) represents the physical field correction value;

[0214] Calculate the local error estimate as follows:

[0215] ;

[0216] in, E represents the local error estimate, ||·|| represents the norm;

[0217] Compare the local error estimate with a preset error threshold, and when the local error estimate is greater than the error threshold, halve the initial time step and return to execute the calculation step;

[0218] When the local error estimate is less than one tenth of the error threshold, the initial time step is doubled and the physical field correction value at the termination time is used as the physical field initial value of the next subinterval;

[0219] When the local error estimate is between one tenth of the error threshold and the error threshold, the initial time step is kept unchanged and the physical field correction value at the termination time is used as the physical field initial value of the next sub-interval;

[0220] The calculation steps are repeatedly performed on the subintervals until the calculation of the entire solution time interval is completed to obtain the numerical solution of the thermal electromagnetic coupling equations.

[0221] A method for numerically solving the thermo-electromagnetic coupling equations divides the entire time solution interval into multiple smaller subintervals. In each subinterval, a four-iteration calculation method based on midpoint estimation and full interval correction is used, and the time step is dynamically adjusted according to the local error.

[0222] First, set the initial time step. At the beginning of the subinterval, use the explicit Euler format to calculate the initial values ​​of each physical field of the thermo-electromagnetic coupling equations. For example, the initial value of the temperature field is 293K, the initial value of the electric field is 0V / m, and the initial value of the magnetic field is 0T.

[0223] Next, the first estimate is made. Based on the current physical field value and the time step, the predicted value of the physical field at the intermediate moment after advancing half a time step is calculated. This is achieved by adding the current physical field value to the product of the physical field function value and half a time step. For example, if the temperature field function value is 10K / s and the time step is 0.1s, then the predicted value of the temperature field at the intermediate moment is 293K+10K / s*0.05s=293.5K.

[0224] Then, a second estimate is made. The predicted value of the physical field at the intermediate moment is used to calculate the predicted value of the physical field at the end moment after the entire time step is advanced. The calculation method is similar to the first estimate, using the predicted value of the physical field at the intermediate moment plus the product of the physical field function value (calculated based on the predicted value at the intermediate moment) and half the time step. For example, assuming that the function value of the temperature field at the intermediate moment is 11K / s, the predicted value of the temperature field at the end moment is 293.5K+11K / s*0.05s=294.05K.

[0225] After that, a correction calculation is performed. The physical field value at the current moment and the predicted value of the physical field at the end moment are used to calculate the physical field correction value at the end moment. The calculation method of the correction value can adopt, for example, the central difference format. Assuming that the central difference is used and the temperature field function value at the end moment is 12K / s, the corrected temperature value at the end moment is 293K+12K / s*0.1s=294.2K.

[0226] Next, perform a local error estimate. Calculate the norm of the difference between the corrected value of the physical field at the end time and the predicted value, and use this as the local error estimate. For example, if the corrected temperature at the end time is 294.2K and the predicted value is 294.05K, the local error is |294.2K-294.05K|=0.15K.

[0227] Compare the calculated local error with the preset error threshold. If the local error is greater than the preset threshold (for example, 0.2K), the time step is halved (for example, changed to 0.05s), and the physical field values ​​at the intermediate and terminal moments are returned to be recalculated. If the local error is less than one tenth of the error threshold (for example, 0.02K), the time step is doubled (for example, changed to 0.2s), and the corrected physical field value at the terminal moment is used as the initial value of the next sub-interval. If the local error is between one tenth of the error threshold and the error threshold (for example, between 0.02K and 0.2K), the time step is kept unchanged, and the corrected physical field value at the terminal moment is used as the initial value of the next sub-interval.

[0228] Repeat the above steps until the calculation of all sub-intervals is completed, and finally obtain the numerical solution of the thermo-electromagnetic coupling equations in the entire time solution interval.

[0229] In this embodiment, by midpoint estimation and full interval correction, as well as step control based on local errors, the numerical calculation error can be effectively reduced and the accuracy of the numerical solution can be improved; by dynamically adjusting the time step, the time step can be increased as much as possible while ensuring the calculation accuracy, thereby reducing the amount of calculation and improving the calculation efficiency; the local error control mechanism can effectively avoid the divergence of the numerical solution and enhance the stability of the algorithm, especially when dealing with complex physical field problems, it can better capture the changing laws of the physical field.

[0230] Figure 2 FIG. 1 is a schematic diagram of the structure of an integrated circuit packaging quality assessment system based on nondestructive testing technology according to an embodiment of the present invention. Figure 2 As shown, the system comprises:

[0231] The first unit is used to place the package structure of the integrated circuit to be detected on the sample stage of the two-photon excitation microscope system, use a femtosecond laser to generate a two-photon excitation pulse, and use an acousto-optic modulator to control the energy density of the two-photon excitation pulse in real time, so that the package structure generates a stable fluorescence response signal; control the confocal scanning unit to scan the package structure layer by layer along the Z-axis direction to obtain the three-dimensional distribution data of the fluorescence response signal; establish a stress distribution characteristic map of the package structure according to the three-dimensional distribution data, and mark the area where the stress gradient exceeds the preset threshold in the stress distribution characteristic map to generate a stress anomaly location map; calculate the stress gradient value for each area in the stress anomaly location map, and establish a stress gradient classification map;

[0232] The second unit is used to determine the detection area according to the stress gradient grading map, arrange the quantum nitrogen vacancy probe array in the detection area, and set the radio frequency pulse sequence parameters according to the spatial distribution of the quantum nitrogen vacancy probe array, so that the quantum nitrogen vacancy probe array can achieve the corresponding detection sensitivity in each detection area; collect the electromagnetic field response signal generated by the quantum nitrogen vacancy probe array under the action of stress, reconstruct the electromagnetic field distribution map of the packaging structure, and spatially correspond it with the stress gradient grading map to form a stress-electromagnetic field coupling model;

[0233] The third unit is used to establish a group of thermomechanical and electromagnetic coupling equations for the packaging structure based on the stress-electromagnetic field coupling model; calculate the evolution characteristics of the group of thermomechanical and electromagnetic coupling equations in each detection area, obtain the interaction relationship between each physical field, and calculate the failure probability of each detection area under the multi-field coupling according to the interaction relationship, and establish a stability evaluation system including mechanical failure risk and electromagnetic failure risk; determine the quality level of the packaging structure by comparing with the stability evaluation system, and generate a quality assessment report for the packaging structure.

[0234] According to a third aspect of the embodiments of the present invention,

[0235] An electronic device is provided, comprising:

[0236] processor;

[0237] a memory for storing processor-executable instructions;

[0238] The processor is configured to call the instructions stored in the memory to execute the aforementioned method.

[0239] A fourth aspect of the embodiments of the present invention is:

[0240] A computer-readable storage medium is provided, on which computer program instructions are stored. When the computer program instructions are executed by a processor, the aforementioned method is implemented.

[0241] The present invention may be a method, an apparatus, a system and / or a computer program product. The computer program product may include a computer-readable storage medium carrying computer-readable program instructions for executing various aspects of the present invention.

[0242] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, rather than to limit it. Although the present invention has been described in detail with reference to the aforementioned embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the aforementioned embodiments, or replace some or all of the technical features therein with equivalents. However, these modifications or replacements do not cause the essence of the corresponding technical solutions to deviate from the scope of the technical solutions of the embodiments of the present invention.

Claims

1. An integrated circuit packaging quality assessment method based on nondestructive testing technology, characterized in that: include: The package structure of the integrated circuit to be inspected is placed on the sample stage of the two-photon excitation microscope system, a femtosecond laser is used to generate a two-photon excitation pulse, and the energy density of the two-photon excitation pulse is controlled in real time by an acousto-optic modulator, so that the package structure generates a stable fluorescence response signal; the confocal scanning unit is controlled to scan the package structure layer by layer along the Z-axis direction to obtain the three-dimensional distribution data of the fluorescence response signal; a stress distribution characteristic map of the package structure is established according to the three-dimensional distribution data, and the area where the stress gradient exceeds the preset threshold is marked in the stress distribution characteristic map to generate a stress anomaly location map; Calculate the stress gradient value for each area in the stress anomaly location map and establish a stress gradient classification map; Determine a detection area according to the stress gradient grading map, arrange a quantum nitrogen vacancy probe array in the detection area, and set radio frequency pulse sequence parameters according to the spatial distribution of the quantum nitrogen vacancy probe array so that the quantum nitrogen vacancy probe array achieves a corresponding detection sensitivity in each detection area; collect the electromagnetic field response signal generated by the quantum nitrogen vacancy probe array under stress, reconstruct the electromagnetic field distribution map of the packaging structure, and spatially correspond it with the stress gradient grading map to form a stress-electromagnetic field coupling model; Establishing a group of thermal-electromagnetic coupling equations for the packaging structure based on the stress-electromagnetic field coupling model; Calculate the evolution characteristics of the thermo-electromagnetic coupling equations in each detection area, obtain the interaction relationship between the physical fields, calculate the failure probability of each detection area under the multi-field coupling effect according to the interaction relationship, and establish a stability evaluation system including mechanical failure risk and electromagnetic failure risk; The quality grade of the packaging structure is determined by comparing the stability evaluation system, and a quality evaluation report of the packaging structure is generated.

2. The method according to claim 1, characterized in that Establishing a stress distribution characteristic diagram of the packaging structure according to the three-dimensional distribution data includes: Inputting the three-dimensional distribution data of the fluorescence response signal into an optical system model, wherein the optical system model includes a lateral resolution parameter determined by the numerical aperture of the objective lens and a longitudinal resolution parameter determined by the excitation wavelength; Generate a Gaussian point spread function according to the optical system model, perform a deconvolution operation on the three-dimensional distribution data and the Gaussian point spread function, and perform Wiener filtering on the result of the deconvolution operation to obtain three-dimensional fluorescence distribution data after eliminating the influence of the optical system; Measure the photoelastic coefficient tensor and refractive index anisotropy coefficient of the material corresponding to the packaging structure, and establish the photoelastic effect response equation of the material stress field and refractive index change; Establishing a nonlinear relationship equation between the refractive index and the fluorescence intensity during two-photon excitation, combining the nonlinear relationship equation with the photoelastic effect response equation to obtain a quantitative mapping model between the stress field and the fluorescence intensity; Selecting an initial stress field, calculating a theoretical fluorescence intensity according to the quantitative mapping model, comparing the theoretical fluorescence intensity with the measured fluorescence intensity of the three-dimensional fluorescence distribution data, and calculating a mean square error as an optimization objective function; performing a gradient calculation, calculating a gradient vector of the optimization objective function, determining a conjugate search direction, performing a conjugate search, and iteratively searching for an optimal step size; The stress field distribution is updated using the optimal step size. When the convergence value of the optimization objective function is greater than a preset convergence threshold, the gradient calculation and conjugate search are continued. When the convergence value of the optimization objective function is less than or equal to the convergence threshold, the current stress field distribution is used as a stress distribution characteristic diagram of the packaging structure.

3. The method according to claim 1, characterized in that Marking the area where the stress gradient exceeds a preset threshold in the stress distribution characteristic map, and generating the stress anomaly location map includes: Using the Sobel operator to calculate the gradient of the stress distribution characteristic diagram in the horizontal direction to obtain a first gradient component matrix, and using the Sobel operator to calculate the gradient of the stress distribution characteristic diagram in the vertical direction to obtain a second gradient component matrix; Calculating the square sum of the first gradient component matrix and the second gradient component matrix, and performing a square root operation to obtain a gradient matrix of the stress distribution characteristic diagram; performing a power-of-2 multi-scale decomposition on the gradient matrix, wherein the multi-scale decomposition ranges from 2×2 to 16×16, and obtaining sub-matrices of multiple scales; Calculating the grayscale mean and grayscale standard deviation of each local area in the sub-matrices of the multiple scales respectively, determining the segmentation threshold of the corresponding local area according to the grayscale mean and the grayscale standard deviation, and determining the adaptive segmentation threshold; Based on an adaptive segmentation threshold, the sub-matrices of the multiple scales are binarized respectively to obtain binarized sub-matrices of the multiple scales; the grayscale variance of each binarized sub-matrix of the multiple scales is calculated respectively, and the grayscale variance is normalized to obtain the fusion weight of each scale; the binarized sub-matrices of the multiple scales are multiplied and superimposed with the corresponding fusion weights respectively to obtain a fused segmentation matrix; the fused segmentation matrix is ​​binarized to obtain a mask matrix, and in the mask matrix, the stress abnormality area is determined based on the area with a value of 1; the stress abnormality area is marked in the stress distribution characteristic map to generate a stress abnormality location map.

4. The method according to claim 1, characterized in that: Determining a detection area according to the stress gradient classification map, arranging a quantum nitrogen vacancy probe array in the detection area, and setting radio frequency pulse sequence parameters according to the spatial distribution of the quantum nitrogen vacancy probe array so that the quantum nitrogen vacancy probe array achieves a corresponding detection sensitivity in each detection area includes: Constructing a quantum state evolution measurement system, applying an initialization radio frequency pulse to the quantum state nitrogen vacancy probe array, and initializing the quantum state nitrogen vacancy probe array to a spin ground state; Using a quantum state readout sequence to obtain the spin state information of the quantum state nitrogen vacancy probe array, record the spin state evolution time, measure the quantum coherence time of the quantum state nitrogen vacancy probe array in each detection area, and establish a quantum state response model; Using a quantum state subspace decomposition method, the quantum state response model is characteristically decomposed to obtain a quantum state characteristic vector of each detection area; a quantum decoherence noise estimator is constructed based on the quantum state characteristic vector, the influence of environmental noise on the quantum state nitrogen vacancy probe of each detection area is measured, and the quantum state decoherence rate of the detection area is calculated; According to the quantum state decoherence rate, the minimum detection sensitivity corresponding to each detection area is determined, and a multi-channel quantum state manipulation timing scheme is established to achieve synchronous quantum state manipulation between detection areas; For each detection area, a quantum state characteristic transformation matrix is ​​constructed to extract quantum coherence information of the quantum state nitrogen vacancy probe to obtain the quantum state response characteristics of the detection area; Inputting the quantum state response characteristics into the quantum information fusion framework, and calculating the sensitivity parameters of each detection area based on the quantum state evolution analysis algorithm; Based on the sensitivity parameters, the radio frequency pulse phase accumulation time corresponding to each detection area is calculated, the radio frequency pulse repetition period is determined, and radio frequency pulse sequence parameters that meet the detection sensitivity requirements of each detection area are generated.

5. The method according to claim 1, characterized in that: Collecting the electromagnetic field response signal generated by the quantum nitrogen vacancy probe array under stress, reconstructing the electromagnetic field distribution diagram of the packaging structure, and spatially corresponding it with the stress gradient classification diagram to form a stress-electromagnetic field coupling model includes: The electromagnetic field response signal of the quantum nitrogen vacancy probe array is collected in each detection area to construct an electromagnetic field response time series; the electromagnetic field response time series is discretely sampled, a signal analytical expression is obtained by Hilbert transform, and the instantaneous value of the electromagnetic field intensity in each detection area is extracted; the instantaneous value of the electromagnetic field intensity is constructed as an electromagnetic field response matrix, and the main eigenvectors of the electromagnetic field distribution are extracted by principal component analysis method; Using continuous wavelet transform to perform time-frequency decomposition on the electromagnetic field response matrix to obtain electromagnetic field spectrum components at different scales; performing adaptive threshold filtering on the electromagnetic field spectrum components to eliminate measurement noise and obtain frequency domain characteristics of the electromagnetic field signal; A compressed sensing reconstruction algorithm is used to sparsely reconstruct the frequency domain features, a mathematical model of electromagnetic field spatial distribution is established, an electromagnetic field inversion equation is constructed based on the mathematical model, and a conjugate gradient iteration method is used to solve the electromagnetic field distribution to obtain an electromagnetic field distribution map of the packaging structure; the electromagnetic field distribution map is affinely transformed and registered with the stress gradient classification map to establish a spatial coordinate correspondence relationship; Based on the spatial coordinate correspondence, a kernel function regression method is used to extract the mapping relationship between the electromagnetic field intensity and the stress gradient, and the electromagnetic field-stress response data is obtained; based on the electromagnetic field-stress response data, a nonlinear coupling equation group is established, and the coupling coefficient is solved by the Newton iteration method; The coupling coefficient is substituted into the Maxwell stress-electromagnetic field tensor equation to form a stress-electromagnetic field coupling model of the packaging structure.

6. The method according to claim 1, characterized in that A group of thermomechanical and electromagnetic coupling equations for the packaging structure is established based on the stress-electromagnetic field coupling model; the evolution characteristics of the group of thermomechanical and electromagnetic coupling equations in each detection area are calculated to obtain the interaction relationship between the physical fields; based on the interaction relationship, the failure probability of each detection area under the multi-field coupling is calculated to establish a stability evaluation system including mechanical failure risk and electromagnetic failure risk; Determine the quality grade of the packaging structure by comparing with the stability evaluation system, and generate a quality assessment report of the packaging structure including: The temperature field is introduced into the stress-electromagnetic field coupling model to construct a stress distribution function including a thermal expansion coefficient; a corresponding relationship between temperature and conductivity is established based on the stress distribution function, and an exponential function including a parameter of the ratio of activation energy to the Boltzmann constant is used to determine the change characteristics of conductivity with temperature; The stress distribution function is combined with a predetermined electromagnetic field basic equation group, a trial function is selected to construct a weight function, a discretization equation is established by piecewise linear interpolation, and a thermal electromagnetic coupling equation group of the packaging structure is obtained by minimizing the residual error; The solution time interval is divided into multiple subintervals, and four iterative calculations including midpoint estimation and full interval correction are performed in each subinterval. The calculation step size is dynamically adjusted based on the local error of the subinterval to obtain the numerical solution of the thermo-electromagnetic coupling equations. Performing signal decomposition operation on the numerical solution to extract characteristic mode functions that characterize the inherent change characteristics of the physical field of the packaging structure; calculating the ratio of the joint probability density to the edge probability density of the characteristic mode functions to obtain the interaction strength between the physical fields in the packaging structure; By calculating the predicted force of the known physical field on the target physical field, a quantitative characterization model of the interaction between the physical fields of the packaging structure is established, and based on the quantitative characterization model, a failure judgment function including the physical field intensity threshold and the load response is constructed; Based on the expectation value minimization criterion, an optimized probability density function is constructed, and multiple random sampling calculations are performed according to the failure judgment function to obtain the failure probability of the packaging structure; Establishing a stability evaluation index system based on the failure probability of the packaging structure, and calculating a membership relationship matrix of each stability evaluation index to the quality grade; The weight coefficient vector of the stability evaluation index is calculated by pairwise comparison, and the weight coefficient vector is multiplied by the membership matrix to obtain the stability evaluation result of the packaging structure.

7. The method according to claim 6, characterized in that The solution time interval is divided into multiple subintervals. Four iterative calculations including midpoint estimation and full interval correction are performed in each subinterval. The calculation step size is dynamically adjusted based on the local error of the subinterval. The numerical solutions of the thermo-electromagnetic coupling equations are obtained, including: Divide the solution time interval into a plurality of subintervals, set an initial time step at the initial time of the subinterval; at the initial time of the subinterval, use an explicit Euler format to calculate the initial value of the physical field of the thermo-electromagnetic coupling equation group; Calculate the predicted value of the physical field at the intermediate time, the formula is as follows: ； in, t i Indicates the current moment, t i+1 / 2 represents the intermediate moment of advancing half a step, u * (·) represents the predicted value of the physical field, u (·) represents the physical field value, Δ t represents the time step, f (·) represents the value of the physical field function; Calculate the predicted value of the physical field at the end time, the formula is as follows: ； in, u * ( t i+1 ) represents the predicted value of the physical field at the termination time; Calculate the physical field correction value at the end time, the formula is as follows: ； in, u ' (·) represents the physical field correction value; Calculate the local error estimate as follows: ； in, E represents the local error estimate, ||·|| represents the norm; Compare the local error estimate with a preset error threshold, and when the local error estimate is greater than the error threshold, halve the initial time step and return to execute the calculation step; When the local error estimate is less than one tenth of the error threshold, the initial time step is doubled and the physical field correction value at the termination time is used as the physical field initial value of the next subinterval; When the local error estimate is between one tenth of the error threshold and the error threshold, the initial time step is kept unchanged and the physical field correction value at the termination time is used as the physical field initial value of the next sub-interval; The calculation steps are repeatedly performed on the subintervals until the calculation of the entire solution time interval is completed to obtain the numerical solution of the thermal electromagnetic coupling equations.

8. An integrated circuit packaging quality assessment system based on non-destructive testing technology, used to implement the method described in any one of claims 1 to 7, characterized in that: include: The first unit is used to place the package structure of the integrated circuit to be detected on the sample stage of the two-photon excitation microscope system, use a femtosecond laser to generate a two-photon excitation pulse, and use an acousto-optic modulator to control the energy density of the two-photon excitation pulse in real time, so that the package structure generates a stable fluorescence response signal; control the confocal scanning unit to scan the package structure layer by layer along the Z-axis direction to obtain three-dimensional distribution data of the fluorescence response signal; establish a stress distribution characteristic map of the package structure according to the three-dimensional distribution data, and mark the area where the stress gradient exceeds a preset threshold in the stress distribution characteristic map to generate a stress anomaly location map; Calculate the stress gradient value for each area in the stress anomaly location map and establish a stress gradient classification map; The second unit is used to determine the detection area according to the stress gradient grading map, arrange the quantum nitrogen vacancy probe array in the detection area, and set the radio frequency pulse sequence parameters according to the spatial distribution of the quantum nitrogen vacancy probe array, so that the quantum nitrogen vacancy probe array can achieve the corresponding detection sensitivity in each detection area; collect the electromagnetic field response signal generated by the quantum nitrogen vacancy probe array under the action of stress, reconstruct the electromagnetic field distribution map of the packaging structure, and spatially correspond it with the stress gradient grading map to form a stress-electromagnetic field coupling model; The third unit is used to establish a group of thermal electromagnetic coupling equations for the packaging structure based on the stress-electromagnetic field coupling model; Calculate the evolution characteristics of the thermo-electromagnetic coupling equations in each detection area, obtain the interaction relationship between the physical fields, calculate the failure probability of each detection area under the multi-field coupling effect according to the interaction relationship, and establish a stability evaluation system including mechanical failure risk and electromagnetic failure risk; The quality grade of the packaging structure is determined by comparing the stability evaluation system, and a quality evaluation report of the packaging structure is generated.

9. An electronic device, characterized in that: include: processor; a memory for storing processor-executable instructions; The processor is configured to call the instructions stored in the memory to execute the method described in any one of claims 1 to 7.

10. A computer-readable storage medium having computer program instructions stored thereon, characterized in that: When the computer program instructions are executed by a processor, the method according to any one of claims 1 to 7 is implemented.

Images