Method and system for adjusting growth process parameters of heteroepitaxial layer of wide bandgap semiconductor
By combining digital twin models with multi-source data monitoring, the heteroepitaxial growth process of wide bandgap semiconductors is monitored and optimized in real time, solving the difficult problems of strain and defect control in traditional processes, achieving efficient technology application, and improving the intelligence and automation level of the process.
Patent Information
- Application Number
- CN202510749332.8
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-06-06
- Publication Date
- 2025-09-23
- Estimated Expiration
- Not applicable · inactive patent
AI Technical Summary
The existing wide bandgap semiconductor heteroepitaxial growth process lacks real-time monitoring and prediction capabilities, making it difficult to accurately control strain and defect evolution. The process parameter optimization cycle is long, the quality is unstable, and there is a lack of data-driven intelligent and automated methods.
The digital twin model is combined with multi-source data monitoring. Data is acquired through in-situ ellipsometers, reflective high-energy electron diffraction, and acoustic emission sensors. Deep neural networks are used to calculate the strain tensor field and defect density distribution, perform hierarchical and progressive process parameter optimization, and conduct closed-loop feedback control.
Precise control of the heteroepitaxial growth process is achieved, strain and defect density are reduced, the crystallization quality of the epitaxial layer is improved, and the repeatability and stability of the process are enhanced.
Smart Images

Figure CN120690343A_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to semiconductor technology, and in particular to a method and system for adjusting process parameters of a wide bandgap semiconductor heteroepitaxial layer growth. Background Art
[0002] Currently, wide-bandgap semiconductor heteroepitaxial growth processes rely primarily on empirical parameters and trial and error, lacking precise control methods based on the physical mechanisms of the material. Traditionally, process parameter adjustments are primarily based on offline characterization and empirical formulas, making it difficult to capture the strain evolution and defect formation mechanisms during growth in real time. Furthermore, the growth process involves the coupling of multiple physical and chemical processes, such as surface adsorption, atomic diffusion, lattice matching, and interfacial reactions. The complex interactions between these processes make it difficult for traditional control methods to comprehensively account for their impact.
[0003] The main deficiencies of existing technologies include: First, there is a lack of real-time monitoring and prediction capabilities for strain field and defect evolution during heteroepitaxial growth, making it difficult to precisely control the material according to the actual state of growth. Most process parameter adjustments rely on post-characterization analysis, and are unable to intervene in time with problems that arise during the growth process, resulting in long process optimization cycles and unstable material quality. Secondly, the coupling effects of multiple physical fields in the heteroepitaxial process are difficult to accurately characterize and quantify, the mutual influence relationship between process parameters is unclear, and there is a lack of systematic parameter optimization methods, resulting in a narrow process window and poor repeatability. Thirdly, existing technologies lack effective data-driven methods to link multi-source monitoring data with the physical mechanisms of materials, and are unable to fully utilize the large amount of data generated during the growth process for knowledge mining and process optimization, resulting in a low level of intelligence and automation in process control. Summary of the Invention
[0004] The embodiments of the present invention provide a method and system for adjusting process parameters of a wide bandgap semiconductor heteroepitaxial layer growth, which can solve the problems in the prior art.
[0005] A first aspect of an embodiment of the present invention provides a method for adjusting process parameters of a wide bandgap semiconductor heteroepitaxial layer growth process, comprising: Obtain the surface morphology parameters, lattice parameters, and surface defect density of the wide bandgap semiconductor substrate to be processed, build a digital twin model to predict the strain evolution trend in the early stages of heteroepitaxial growth, and determine the optimal initial process parameter set corresponding to the minimum strain and defect density as the starting process point for heteroepitaxial growth; During the heteroepitaxial layer growth process, an in-situ ellipsometer is used to obtain growth rate data and refractive index data, reflection high-energy electron diffraction is used to obtain surface reconstruction data and lattice constant data, and an acoustic emission sensor is used to obtain stress wave signal data. The multi-source data obtained are input into a deep neural network to calculate the strain tensor field and defect density distribution; Based on the strain tensor field and defect density distribution, a hierarchical and progressive process parameter optimization is performed. In the first layer, the stress wave signal data is used to construct a parameter coupling influence matrix. In the second layer, the base surface reconstruction data and lattice constant data are used to establish a dislocation density evolution prediction model for defect evolution. In the third layer, the process parameter adjustment strategy is determined based on the growth rate data and refractive index data. The growth process is adjusted according to the process parameters after layered optimization, the actual process parameters and monitoring data are fed back to the digital twin model, and the parameter coupling influence matrix and dislocation density evolution prediction model are updated.
[0006] In an optional embodiment, A digital twin model was constructed to predict the strain evolution trend in the early stages of heteroepitaxial growth and determine the optimal initial process parameter set, including: The digital twin model includes a strain energy density function module, a dislocation motion dynamics equation module, and a critical layer thickness calculation module. The surface morphology parameters, lattice parameters, and surface defect density are input into the digital twin model as initial conditions and constraints. The digital twin model constructs a coupled simulation environment of the temperature field distribution equation, the airflow field distribution equation, and the stress field distribution equation, and calculates a parameter coupling coefficient matrix of the temperature parameter, the gas flow parameter, the cavity pressure parameter, and the radio frequency power parameter based on the coupled simulation environment; The parameter coupling coefficient matrix is input into a nonlinear response optimization module, and the nonlinear response optimization module calculates the lattice strain value based on the strain energy density function module and calculates the defect density value based on the dislocation motion dynamics equation module; with the critical layer thickness value as a constraint condition, a weighted objective function is constructed with the lattice strain value and the defect density value, and a constrained gradient descent algorithm is used to iteratively optimize the weighted objective function until a minimum value is obtained, and an optimal initial process parameter group is output, wherein the optimal initial process parameter group includes temperature parameters, gas flow parameters, cavity pressure parameters, and radio frequency power parameters.
[0007] In an optional embodiment, The acquired multi-source data is input into the deep neural network to calculate the strain tensor field and defect density distribution including: Constructing a multi-scale feature extraction network, the multi-scale feature extraction network including a time-scale attention module and a space-scale convolution module, and inputting the multi-source data into the multi-scale feature extraction network to obtain a fused feature vector; Constructing a long short-term memory unit network, inputting the fused feature vector into the long short-term memory unit network, and the long short-term memory unit network calculating a temporal dependency feature based on historical state information and current input information; Constructing a strain prediction network and a defect density prediction network, inputting the timing-dependent features into the strain prediction network to obtain a strain tensor field, wherein the strain tensor field satisfies a continuity equation constraint; inputting the timing-dependent features and the strain tensor field into the defect density prediction network to obtain a defect density distribution, wherein the defect density distribution satisfies an energy conservation constraint; A multi-task joint loss function is constructed based on the strain tensor field and defect density distribution. The multi-task joint loss function includes a strain prediction loss term, a defect density prediction loss term and a regularization loss term. An adaptive learning rate is used to optimize the training of the multi-scale feature extraction network, the long short-term memory unit network, the strain prediction network and the defect density prediction network.
[0008] In an optional embodiment, Based on the strain tensor field and defect density distribution, a hierarchical and progressive process parameter optimization is performed. In the first layer, the stress wave signal data is used to construct a parameter coupling influence matrix. In the second layer, the base surface reconstruction data and lattice constant data are used to establish a dislocation density evolution prediction model for defect evolution. In the third layer, the process parameter adjustment strategy is determined based on the growth rate data and refractive index data. The following are included: By calculating the time-frequency spectrum and energy transfer function of the stress wave signal data, the sensitivity coefficients between the process parameters are obtained, and a parameter coupling influence matrix is established; Establishing a reconstruction response function using surface reconstruction data, calculating a surface reconstruction phase transition characteristic time constant based on the reconstruction response function, calculating a critical thickness value based on lattice constant data, and constructing a dislocation density evolution prediction model based on the surface reconstruction phase transition characteristic time constant and the critical thickness value; The current growth stage is identified based on the growth rate data, the component uniformity is evaluated based on the refractive index data, and the optimization target is determined based on the growth stage and the component uniformity; the adjustment direction of the process parameters is determined based on the optimization target and the parameter coupling influence matrix, and the adjustment range of the process parameters is determined based on the dislocation density evolution prediction model to obtain the adjustment amount of each process parameter.
[0009] In an optional embodiment, Constructing a dislocation density evolution prediction model includes: The surface reconstruction data is Fourier transformed to obtain the reconstruction period characteristics, and the reconstruction response function is established in combination with the reconstruction intensity change trend; the phase change inflection point is determined based on the first-order derivative and the second-order derivative of the reconstruction response function, and the time interval between adjacent inflection points is calculated to obtain the surface reconstruction phase change characteristic time constant; Substituting the lattice constant data into the elastic strain energy equation and the plastic strain energy equation, respectively, the total strain energy distribution of the heterogeneous interface is calculated, and a dislocation nucleation criterion function is established; solving the strain energy balance equation based on the dislocation nucleation criterion function to obtain the critical layer thickness, and considering the influence of the interface step effect and dislocation interaction on the critical layer thickness to obtain the critical thickness value; A dislocation nucleation-extension evolution equation is constructed based on the surface reconstruction phase transition characteristic time constant and the critical thickness value, and a dislocation density evolution prediction model is established.
[0010] In an optional embodiment, Identifying the current growth stage based on the growth rate data, evaluating the composition uniformity based on the refractive index data, and determining the optimization target based on the growth stage and composition uniformity include: Calculating the first-order derivative of the growth rate data to obtain the rate change rate, calculating the mean square error of the growth rate data to obtain the rate fluctuation; calculating the gradient of the refractive index data in the depth direction to obtain the component distribution function, and calculating the second-order derivative of the refractive index data in the plane to obtain the component uniformity; constructing a nucleation period discriminant function based on the rate change rate and the rate fluctuation, constructing a stable period discriminant function based on the rate change rate, the rate fluctuation and the component uniformity, and identifying the current growth stage based on the nucleation period discriminant function and the stable period discriminant function; An interface quality control function is constructed based on the square integral of the gradient of the component distribution function and the rate fluctuation during the nucleation period, and a process stability control function is constructed based on the component uniformity and the second-order derivative of the rate change rate during the stabilization period; According to the identified growth stage, the weight coefficients of the interface quality control function and the process stability control function are dynamically adjusted using an exponential decay function, and the weighted control functions are combined to construct an optimization objective function.
[0011] In an optional embodiment, Feeding back actual process parameters and monitoring data to the digital twin model and updating the parameter coupling influence matrix and dislocation density evolution prediction model include: Establishing a deviation evaluation function between actual process parameters and theoretically predicted parameters, and updating the strain energy density function module and dislocation motion dynamics equation module in the digital twin model based on the output result of the deviation evaluation function; Calculating a cross-correlation function between the actually monitored stress wave signal and the predicted stress wave signal, and updating a sensitivity coefficient in the parameter coupling influence matrix based on peak and delay characteristics of the cross-correlation function; The deviation between the actually observed surface reconstruction period and the predicted period is substituted into the correction equation of the reconstruction response function, and the deviation between the actually measured critical thickness and the predicted critical thickness is substituted into the correction equation of the dislocation nucleation criterion function. The evolution equation parameters in the dislocation density evolution prediction model are updated based on the corrected function.
[0012] A second aspect of an embodiment of the present invention provides a system for adjusting process parameters for growing a wide bandgap semiconductor heteroepitaxial layer, comprising: The first unit is used to obtain the surface morphology parameters, lattice parameters, and surface defect density of the wide bandgap semiconductor substrate to be processed, build a digital twin model to predict the strain evolution trend in the early stage of heteroepitaxial growth, and determine the optimal initial process parameter set corresponding to the minimum strain and defect density as the starting process point for heteroepitaxial growth; The second unit is used to obtain growth rate data and refractive index data using in-situ ellipsometers during the growth of heteroepitaxial layers, obtain surface reconstruction data and lattice constant data using reflection high-energy electron diffraction, and obtain stress wave signal data using acoustic emission sensors. The obtained multi-source data are input into a deep neural network to calculate the strain tensor field and defect density distribution; A third unit is configured to perform hierarchical and progressive process parameter optimization based on the strain tensor field and defect density distribution, constructing a parameter coupling influence matrix using stress wave signal data in the first layer, establishing a dislocation density evolution prediction model for defect evolution using base surface reconstruction data and lattice constant data in the second layer, and determining a process parameter adjustment strategy based on growth rate data and refractive index data in the third layer; The fourth unit is used to adjust the growth process according to the process parameters after layered optimization, feed back the actual process parameters and monitoring data to the digital twin model, and update the parameter coupling influence matrix and dislocation density evolution prediction model.
[0013] According to a third aspect of an embodiment of the present invention, an electronic device is provided, including: processor; a memory for storing processor-executable instructions; The processor is configured to call the instructions stored in the memory to execute the aforementioned method.
[0014] According to a fourth aspect of an embodiment of the present invention, 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 method described above is implemented.
[0015] The wide bandgap semiconductor heteroepitaxial layer growth process parameter adjustment method provided by the present invention combines digital twin technology with in-situ multi-source data monitoring to achieve precise control of the heteroepitaxial growth process, significantly reduce the strain and defect density in the wide bandgap semiconductor heterostructure, and improve the crystallization quality of the epitaxial layer.
[0016] The present invention establishes a multi-source data fusion analysis framework based on deep neural networks, which monitors and accurately calculates the strain tensor field and defect density distribution in real time. Through a hierarchical and progressive process parameter optimization strategy, it realizes dynamic regulation of the growth process, effectively solving the regulation difficulty caused by parameter coupling in traditional processes.
[0017] The process parameter adjustment method of the present invention realizes closed-loop feedback control, continuously updates the actual monitoring data to the digital twin model, continuously improves the model prediction accuracy, enhances the repeatability and stability of the process, and provides technical support for the large-scale preparation of high-quality wide-bandgap semiconductor heterostructures. BRIEF DESCRIPTION OF THE DRAWINGS
[0018] Figure 1 The figure is a flow chart of a method for adjusting process parameters of a wide bandgap semiconductor heteroepitaxial layer growth according to an embodiment of the present invention. DETAILED DESCRIPTION
[0019] To make the objectives, technical solutions, and advantages of the embodiments of the present invention more clear, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the accompanying 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 making creative efforts shall fall within the scope of protection of the present invention.
[0020] The following specific embodiments are used to describe the technical solution of the present invention in detail. 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.
[0021] Figure 1 FIG. 1 is a flow chart of a method for adjusting process parameters of a wide bandgap semiconductor heteroepitaxial layer growth according to an embodiment of the present invention. Figure 1 As shown, the method includes: Obtain the surface morphology parameters, lattice parameters, and surface defect density of the wide bandgap semiconductor substrate to be processed, build a digital twin model to predict the strain evolution trend in the early stages of heteroepitaxial growth, and determine the optimal initial process parameter set corresponding to the minimum strain and defect density as the starting process point for heteroepitaxial growth; During the heteroepitaxial layer growth process, an in-situ ellipsometer is used to obtain growth rate data and refractive index data, reflection high-energy electron diffraction is used to obtain surface reconstruction data and lattice constant data, and an acoustic emission sensor is used to obtain stress wave signal data. The multi-source data obtained are input into a deep neural network to calculate the strain tensor field and defect density distribution; Based on the strain tensor field and defect density distribution, a hierarchical and progressive process parameter optimization is performed. In the first layer, the stress wave signal data is used to construct a parameter coupling influence matrix. In the second layer, the base surface reconstruction data and lattice constant data are used to establish a dislocation density evolution prediction model for defect evolution. In the third layer, the process parameter adjustment strategy is determined based on the growth rate data and refractive index data. The growth process is adjusted according to the process parameters after layered optimization, the actual process parameters and monitoring data are fed back to the digital twin model, and the parameter coupling influence matrix and dislocation density evolution prediction model are updated.
[0022] In an optional embodiment, a digital twin model is constructed to predict the strain evolution trend in the early stage of heteroepitaxial growth and determine the optimal initial process parameter set, including: The digital twin model includes a strain energy density function module, a dislocation motion dynamics equation module, and a critical layer thickness calculation module. The surface morphology parameters, lattice parameters, and surface defect density are input into the digital twin model as initial conditions and constraints. The digital twin model constructs a coupled simulation environment of the temperature field distribution equation, the airflow field distribution equation, and the stress field distribution equation, and calculates a parameter coupling coefficient matrix of the temperature parameter, the gas flow parameter, the cavity pressure parameter, and the radio frequency power parameter based on the coupled simulation environment; The parameter coupling coefficient matrix is input into a nonlinear response optimization module, and the nonlinear response optimization module calculates the lattice strain value based on the strain energy density function module and calculates the defect density value based on the dislocation motion dynamics equation module; with the critical layer thickness value as a constraint condition, a weighted objective function is constructed with the lattice strain value and the defect density value, and a constrained gradient descent algorithm is used to iteratively optimize the weighted objective function until a minimum value is obtained, and an optimal initial process parameter group is output, wherein the optimal initial process parameter group includes temperature parameters, gas flow parameters, cavity pressure parameters, and radio frequency power parameters.
[0023] For example, the construction process of the digital twin model first needs to establish a strain energy density function module, a dislocation motion dynamics equation module and a critical layer thickness calculation module. The strain energy density function module uses the elastic theory framework to calculate the strain energy density caused by lattice mismatch during heteroepitaxial growth. For given lattice parameters a1 and a2, the strain energy density is related to the lattice mismatch f=(a1-a2) / a2. In actual calculations, when the lattice mismatch is 2%, the strain energy density is about 23.5 mJ / m 2 . The dislocation motion dynamics equation module describes the dislocation nucleation and slip process, taking into account parameters such as shear stress, Poisson's ratio and critical shear stress. For example, when the external stress reaches 42 MPa, the dislocation begins to move, and the speed is nonlinearly related to the stress. The critical layer thickness calculation module is based on the energy balance principle. When the epitaxial layer thickness exceeds the critical value, the strain will be released through dislocation formation. For germanium-silicon heterostructures, when the germanium content is 25%, the critical layer thickness is about 9.7 nm.
[0024] The temperature field distribution equation uses the finite element method to solve the heat conduction equation, taking into account the thermal conductivity, density, and specific heat capacity of the material. In the heteroepitaxial reaction chamber, when the center temperature is set to 800°C, the temperature at the edge of the chamber is approximately 780°C, resulting in a temperature gradient of 20°C. The gas flow field distribution equation is based on the Navier-Stokes equations, taking into account parameters such as gas viscosity, density, and pressure. In actual gas flow simulations, when the gas flow rate is 200 sccm and the chamber pressure is 50 Pa, the gas flow rate at the center of the reaction zone is approximately 0.75 m / s. The stress field distribution equation considers the combined effects of thermal stress and lattice mismatch stress and is solved using the continuum mechanics equations. In actual calculations, when the temperature is 850°C and the lattice mismatch is 1.5%, the average stress on the epitaxial layer surface is approximately 275 MPa.
[0025] The parameter coupling coefficient matrix is constructed by designing multiple sets of experiments and recording the effects of different input parameters on the output results. For example, a 10°C change in temperature affects strain by 0.12%, and a 10 sccm change in gas flow affects strain by 0.05%. The complete parameter coupling coefficient matrix is a 4×2 matrix, containing the coefficients for the effects of temperature, gas flow, chamber pressure, and RF power on lattice strain and defect density.
[0026] The nonlinear response optimization module first calculates the lattice strain value and defect density value under different process parameter combinations. The lattice strain is calculated by the strain energy density function. The current parameters are substituted into the strain energy density, which is then converted into the lattice strain value. For example, under the conditions of a temperature of 830°C, a gas flow rate of 180 sccm, a chamber pressure of 45 Pa, and an RF power of 120 W, the calculated lattice strain value is 0.98%. The defect density is calculated by the dislocation motion dynamics equation. The current process parameters are substituted into the dislocation nucleation rate and expansion rate formula to calculate the number of defects per unit area. Under the above conditions, the calculated defect density is approximately 5.6×10 8 cm -2 .
[0027] During the construction of the weighted objective function, the lattice strain and defect density need to be standardized to make their dimensions consistent. After standardization, the weight coefficients w1 and w2 are set for the lattice strain and defect density, respectively, usually w1+w2=1. In this embodiment, w1=0.6, w2=0.4, and more emphasis is placed on the control of lattice strain. The weighted objective function F=0.6×(S / S0)+0.4×(D / D0), where S is the current lattice strain, S0 is the target lattice strain, D is the current defect density, and D0 is the allowed defect density threshold.
[0028] During the optimization process of the constrained gradient descent algorithm, the critical layer thickness is set as a constraint. If the layer thickness caused by the current process parameters exceeds the critical value, a penalty term is applied to increase the objective function value. The initial step size of the algorithm is set to 0.05, and the convergence threshold is 0.001. In each iteration, the partial derivatives of the objective function with respect to each parameter are calculated, and the parameters are updated in the opposite direction of the gradient. If the new parameters violate the constraints, a boundary penalty is applied. In the actual optimization process, after 32 iterations of convergence, the optimal initial process parameter group was finally obtained: the temperature parameter was 825±2℃, the gas flow parameter was 175±5 sccm, the chamber pressure parameter was 48±1 Pa, and the RF power parameter was 125±3 W.
[0029] Verification experiments show that when the optimal initial process parameter set is used for heteroepitaxial growth, the lattice strain of the epitaxial layer is 0.92% and the defect density is 3.2×10 8 cm -2 , and the critical layer thickness is 12.5 nm, which are better than conventional parameters. In the comparative experiment, under conventional process parameters (temperature 850℃, gas flow rate 200 sccm, chamber pressure 40 Pa, RF power 150 W), the lattice strain is 1.35%, and the defect density is 8.7×10 8 cm -2, with a critical layer thickness of 8.9 nm. This digital twin model can predict the impact of process parameters on strain evolution in advance, optimize initial process parameters, reduce the number of actual experiments, and lower R&D costs. The model is applicable to a variety of semiconductor heterostructures, such as SiGe / Si and GaAs / Si heteroepitaxial systems.
[0030] The present invention realizes nonlinear optimization of process parameters by constructing a digital twin model that includes strain energy density, dislocation motion and critical layer thickness calculation, and combines the coupled simulation of temperature field, airflow field and stress field. It can accurately predict the strain evolution trend in the early stage of heteroepitaxial growth and provide a theoretical basis for determining the optimal initial process parameters.
[0031] In an optional embodiment, inputting the acquired multi-source data into a deep neural network to calculate the strain tensor field and defect density distribution includes: Constructing a multi-scale feature extraction network, the multi-scale feature extraction network including a time-scale attention module and a space-scale convolution module, and inputting the multi-source data into the multi-scale feature extraction network to obtain a fused feature vector; Constructing a long short-term memory unit network, inputting the fused feature vector into the long short-term memory unit network, and the long short-term memory unit network calculating a temporal dependency feature based on historical state information and current input information; Constructing a strain prediction network and a defect density prediction network, inputting the timing-dependent features into the strain prediction network to obtain a strain tensor field, wherein the strain tensor field satisfies a continuity equation constraint; inputting the timing-dependent features and the strain tensor field into the defect density prediction network to obtain a defect density distribution, wherein the defect density distribution satisfies an energy conservation constraint; A multi-task joint loss function is constructed based on the strain tensor field and defect density distribution. The multi-task joint loss function includes a strain prediction loss term, a defect density prediction loss term and a regularization loss term. An adaptive learning rate is used to optimize the training of the multi-scale feature extraction network, the long short-term memory unit network, the strain prediction network and the defect density prediction network.
[0032] Exemplarily, the multi-source data preprocessing steps include normalization, noise filtering, data alignment, etc.
[0033] The multiscale feature extraction network consists of a time-scale attention module and a spatial-scale convolution module. The time-scale attention module consists of three fully connected layers, with 512 neurons in the input layer, 256 neurons in the middle layer, and 512 neurons in the output layer. Each layer is followed by a Reluctant Relu (ReLU) activation function. This module calculates a weight matrix W_t of dimension 512×512, which is used to weight the temporal features. The spatial-scale convolution module consists of four convolutional layers with kernel sizes of 7×7, 5×5, 3×3, and 3×3, respectively, with identical padding and a stride of 1. The number of channels per layer is 64, 128, 256, and 512, respectively. Each layer is followed by batch normalization and a LeakyReLU activation function. After processing by these two modules, a fused feature vector V_f of dimension 512 is obtained. In practical applications, for fatigue damage monitoring of steel materials, high-resolution surface topography images, stress-strain curve data, and ambient temperature data can be input into this network to extract a fused feature vector containing multiscale information about the material state.
[0034] A long short-term memory (LSTM) network is used to capture the temporal dependencies of material evolution. This network consists of two layers of LSTM units, each with 128 hidden units. The input is the fused feature vector V_f obtained in the previous step, and the output is the temporal dependency feature V_t, with a dimension of 128. Each LSTM unit contains an input gate, a forget gate, an output gate, and a cell state. Through these gating mechanisms, the network is able to learn long-term dependencies. When addressing fatigue crack growth in aluminum alloys, this network effectively memorizes crack growth history and predicts the next stage of growth behavior. Experimental results show that the prediction accuracy reaches 92%.
[0035] The strain prediction network is responsible for calculating the strain tensor field based on temporal dependency features. This network employs a fully convolutional architecture consisting of five convolutional layers with kernel sizes of 3×3 and channel numbers of 128, 64, 32, 16, and 6, respectively. The last layer outputs six channels, corresponding to the six independent components of the strain tensor. To ensure that the strain tensor field satisfies the continuity constraints of the continuity equation, a physical constraint layer is added at the end of the network. By explicitly calculating the spatial derivatives of the strain components and applying a penalty term, the prediction results are made consistent with physical laws. In strain prediction tests on carbon fiber composites, this network achieved a relative error of less than 3.5%.
[0036] The defect density prediction network predicts the internal defect distribution of materials based on time-dependent features and the strain tensor field. The network consists of six convolutional layers. The first three layers have a 3×3 kernel size and 256, 128, and 64 channels, respectively. The last three layers use transposed convolutions for upsampling, with a 4×4 kernel size, a stride of 2, and 32, 16, and 1 channels, respectively. The input is the channel-wise concatenation of the time-dependent features V_t and the strain tensor field, and the output is a single-channel defect density field. To ensure energy conservation, an energy calculation module is incorporated into the network to calculate the total system energy of the predicted defect distribution and compare it with the theoretical value to ensure that the deviation is within an acceptable range.
[0037] The construction of a multi-task joint loss function is crucial for network training. The strain prediction loss term uses mean squared error (MSE) to calculate the difference between the predicted and true values, with a weight of 0.4. The defect density prediction loss term combines the MSE and structural similarity metrics, with a weight of 0.4. The regularization loss term includes an L2 norm penalty and a physical constraint term, with a weight of 0.2. Adaptive learning rate is implemented using the Adam optimizer, with an initial learning rate of 0.001, which is decayed by 0.8 every 50 epochs, for a total of 300 training epochs. Mini-batch training is used with a batch size of 32.
[0038] The present invention adopts a multi-scale feature extraction network to process multi-source data, combines it with a long short-term memory unit network to capture timing characteristics, and realizes the precise calculation of the strain tensor field and defect density distribution through the joint optimization of the strain prediction network and the defect density prediction network, thereby improving the accuracy and robustness of the prediction model.
[0039] In an optional embodiment, a hierarchical and progressive process parameter optimization is performed based on the strain tensor field and defect density distribution. In the first layer, stress wave signal data is used to construct a parameter coupling influence matrix. In the second layer, a dislocation density evolution prediction model for defect evolution is established based on base surface reconstruction data and lattice constant data. In the third layer, a process parameter adjustment strategy is determined based on growth rate data and refractive index data. The following steps are included: By calculating the time-frequency spectrum and energy transfer function of the stress wave signal data, the sensitivity coefficients between the process parameters are obtained, and a parameter coupling influence matrix is established; Establishing a reconstruction response function using surface reconstruction data, calculating a surface reconstruction phase transition characteristic time constant based on the reconstruction response function, calculating a critical thickness value based on lattice constant data, and constructing a dislocation density evolution prediction model based on the surface reconstruction phase transition characteristic time constant and the critical thickness value; The current growth stage is identified based on the growth rate data, the component uniformity is evaluated based on the refractive index data, and the optimization target is determined based on the growth stage and the component uniformity; the adjustment direction of the process parameters is determined based on the optimization target and the parameter coupling influence matrix, and the adjustment range of the process parameters is determined based on the dislocation density evolution prediction model to obtain the adjustment amount of each process parameter.
[0040] For example, in the first-level optimization, stress wave signal data is used to construct a parameter coupling influence matrix. Stress wave signal data is collected using an acoustic emission sensor with a sampling frequency of 10 kHz and a signal sampling length of 1024 points. A short-time Fourier transform is performed on the collected stress wave signal data to obtain time-frequency spectrum characteristics. A Hamming window is selected as the time window, with a window length of 128 points and a 50% overlap ratio. The energy transfer function is calculated based on the time-frequency spectrum by integrating the energy distribution across different frequency bands. In actual implementation, the energy is primarily concentrated in the 200 Hz to 2 kHz frequency band.
[0041] Subsequently, by varying process parameters such as temperature (±10°C), gas flow rate (±5 sccm), chamber pressure (±0.1 Torr), and RF power (±10 W), the rate of change of the stress wave energy transfer function (ESTF) was observed and the sensitivity coefficient was calculated. For example, if the rate of change of the ETF was 3.5% when the temperature increased from 800°C to 810°C, the corresponding sensitivity coefficient for the temperature parameter was 0.35% / °C. A similar method was used to obtain the sensitivity coefficients between the various process parameters and construct a parameter coupling influence matrix. In the GaN / AlN heteroepitaxial growth process, this matrix showed a coupling coefficient of 0.28 between temperature and gas flow rate, indicating a strong mutual influence between the two parameters.
[0042] In the second-level optimization, a dislocation density evolution prediction model for defect evolution was established based on surface reconstruction data and lattice constant data. A reconstruction response function was established using the surface reconstruction data, and the time-varying curve of the reconstruction intensity was fitted to a sigmoid function. Based on this reconstruction response function, its first and second derivatives were calculated to determine the inflection points of the response curve. In practice, two inflection points typically appear during the transition from 1×1 to 2×2 reconstruction on the GaN surface: one approximately 3.2 seconds after the start of reconstruction, and the other approximately 7.8 seconds. The time interval between these two adjacent inflection points was calculated to obtain the characteristic time constant of the surface reconstruction phase transition, which in this example was 4.6 seconds. The lattice constant data were substituted into the elastic strain energy equation and the plastic strain energy equation, respectively, to calculate the total strain energy distribution at the heterojunction interface. At the GaN / AlN interface, strain energy is primarily concentrated within a 5nm radius near the interface. A dislocation nucleation criterion function was established based on the strain energy distribution. Dislocations begin to form when the strain energy density at the interface exceeds a critical value (approximately 0.85 J / m²).
[0043] Based on the dislocation nucleation criterion function, the strain energy balance equation is solved to obtain the critical layer thickness. In the GaN / AlN system, the initial calculated critical layer thickness is about 2.3nm. Further considering the interface step effect (step height is about 0.25nm, density is about 1.2×10 6 / cm) and dislocation interaction (interaction distance is about 5b, b is the Burgers vector) on the critical layer thickness, and finally the corrected critical thickness value is 2.8nm.
[0044] Based on the characteristic time constant and critical thickness of the surface reconstruction phase transition, a dislocation nucleation-propagation evolution equation was constructed, and a dislocation density evolution prediction model was established. This model can predict the variation trend of dislocation density with growth time and thickness under given growth conditions. In the implementation verification, for a GaN / AlN heterostructure with a growth rate of 0.8μm / h, the model predicted that the dislocation density near the interface would increase from an initial 3.2×10 9 / cm 2 Gradually decreases to 8.5×10 8 / cm 2 , and the deviation from the experimental measurement results is less than 15%.
[0045] In the third level of optimization, the process parameter adjustment strategy is determined based on the growth rate and refractive index data. First, the growth rate and refractive index data are acquired using an in-situ ellipsometer with an operating wavelength range of 400-800nm, a measurement angle of 70 degrees, and a sampling frequency of 1Hz.
[0046] The current growth stage is identified based on growth rate data. The specific method is to calculate the first-order derivative of the growth rate data to obtain the rate change rate, and the mean square error of the growth rate data to obtain the rate fluctuation. During the initial growth stage of GaN, the rate change rate is typically greater than 0.05μm / h·min, and the rate fluctuation is greater than 0.08μm / h.
[0047] At the same time, the compositional uniformity is evaluated based on the refractive index data. The compositional distribution function is calculated by calculating the gradient of the refractive index data in the depth direction, and the compositional uniformity is calculated by calculating the second-order derivative of the refractive index data in the plane. For high-quality GaN films, the coefficient of variation of the compositional uniformity in the plane should be less than 0.03.
[0048] A discriminant function for the nucleation phase was constructed by combining the rate change rate and rate fluctuation. When the function value was greater than 0.7, the nucleation phase was determined. A discriminant function for the stability phase was constructed by combining the rate change rate, rate fluctuation, and component uniformity. When the function value was greater than 0.8, the stability phase was determined. The current growth stage was identified based on the discriminant results.
[0049] Different optimization objectives are determined for different growth stages. During the nucleation phase, an interface quality control function is constructed based on the squared integral of the gradient of the component distribution function and the rate fluctuation, focusing on controlling interface flatness and component uniformity. During the stabilization phase, a process stability control function is constructed based on the second-order derivative of the component uniformity and rate change rate, focusing on controlling the stability of the growth process and the consistency of film quality.
[0050] Based on the determined optimization objectives and the previously constructed parameter coupling influence matrix, the process parameter adjustments are determined. For example, if uneven distribution of interface components is observed during the nucleation phase, the parameter coupling influence matrix indicates that increasing the temperature by 2°C and reducing the gas flow by 1.5 sccm can improve this problem. Furthermore, the dislocation density evolution prediction model is used to determine the extent of process parameter adjustments to avoid over-adjustments that could lead to new defects. Ultimately, the adjustment amounts for each process parameter are determined, enabling fine-tuning of the process parameters.
[0051] This method organically combines physical models with data-driven approaches, fully leveraging multi-source monitoring data to accurately understand the coupling relationships between process parameters and the evolution of defects, as well as fine-grained control based on the growth stage. Compared with traditional empirical control methods, it can significantly improve the quality and uniformity of heteroepitaxial layers, reduce defect density, and enhance device performance. It also shortens process optimization cycles and reduces R&D costs, which is of great significance for the performance improvement and industrialization of wide-bandgap semiconductor devices.
[0052] In an optional embodiment, constructing a dislocation density evolution prediction model includes: The surface reconstruction data is Fourier transformed to obtain the reconstruction period characteristics, and the reconstruction response function is established in combination with the reconstruction intensity change trend; the phase change inflection point is determined based on the first-order derivative and the second-order derivative of the reconstruction response function, and the time interval between adjacent inflection points is calculated to obtain the surface reconstruction phase change characteristic time constant; Substituting the lattice constant data into the elastic strain energy equation and the plastic strain energy equation, respectively, the total strain energy distribution of the heterogeneous interface is calculated, and a dislocation nucleation criterion function is established; solving the strain energy balance equation based on the dislocation nucleation criterion function to obtain the critical layer thickness, and considering the influence of the interface step effect and dislocation interaction on the critical layer thickness to obtain the critical thickness value; A dislocation nucleation-extension evolution equation is constructed based on the surface reconstruction phase transition characteristic time constant and the critical thickness value, and a dislocation density evolution prediction model is established.
[0053] For example, surface reconstruction data is processed to obtain reconstruction periodic characteristics. A reconstruction response function is established based on the temporal trend of the reconstruction intensity. The RHEED fringe intensity variation data is extracted and normalized to fit an S-shaped curve function. During the GaN surface reconstruction process, this response function exhibits a typical sigmoid morphology, with slow initial growth, rapid growth in the middle phase, and finally leveling off. Based on this reconstruction response function, its first-order and second-order derivatives are calculated to determine the location of the phase transition inflection point. The extreme value of the first-order derivative indicates the moment when the response function changes at the highest rate, corresponding to the time point when the phase transition is most active; the zero point of the second-order derivative corresponds to the inflection point where the rate of change changes from increasing to decreasing. During the GaN reconstruction process from 1×1 to 2×2, two significant inflection points are typically observed, located approximately 2.8 seconds and 7.3 seconds after the start of reconstruction. The time interval between these two adjacent inflection points is calculated to obtain the characteristic time constant of the surface reconstruction phase transition, which is 4.5 seconds in this example. This time constant reflects the kinetics of surface atomic rearrangement and is an important parameter for constructing dislocation evolution models.
[0054] The lattice constant data is processed to calculate the strain energy distribution. The in-plane lattice constant of the epitaxial layer can be accurately measured by the change in the spacing of the RHEED diffraction fringes. Combined with the out-of-plane lattice constants measured by X-ray diffraction (XRD), complete lattice parameter information is obtained. In the GaN / AlN heterostructure, the lattice constants a of GaN are 3.189Å and c are 5.185Å, while the lattice constants a of AlN are 3.112Å and c are 4.982Å. Substitute these lattice constant data into the elastic strain energy equation to calculate the elastic strain energy distribution caused by lattice mismatch. At the same time, substitute the plastic strain energy equation to calculate the plastic strain energy released after dislocation formation. The difference between the two is the total strain energy distribution. At the GaN / AlN interface, the calculated initial elastic strain energy density is approximately 1.2J / m 2 .
[0055] Based on the above strain energy distribution, a dislocation nucleation criterion function is established. This criterion function expresses the critical condition for dislocation formation in the material, that is, when the strain energy density at the interface exceeds a specific threshold of the material, dislocations begin to nucleate. For the GaN / AlN system, this threshold is approximately 0.9 J / m 2 As the epitaxial layer thickness increases, the accumulated strain energy also increases until it exceeds the threshold. Solving the strain energy balance equation based on this dislocation nucleation criterion function yields the critical layer thickness. For GaN / AlN heterostructures, a preliminary calculation indicates a critical layer thickness of approximately 2.5 nm.
[0056] Further considering the influencing factors in actual material growth, the critical layer thickness is modified. First, the interface step effect is considered. On the actual grown GaN / AlN interface, there are atomic-level steps with a step height of about 0.25nm and a step density of about 5×106 / cm. These steps will reduce the energy barrier for dislocation nucleation, making the actual critical layer thickness smaller than the theoretical value. After correction, the critical layer thickness is reduced by about 0.3nm. Secondly, considering the interaction of dislocations, when multiple dislocations are formed, the stress fields between them will affect each other and form a dislocation network. In GaN, the effective interaction distance between dislocations is about 4 times the Burgers vector (about 1nm). After considering this effect, the critical layer thickness increases by about 0.5nm. Combining these two effects, the final corrected critical thickness value is 2.7nm.
[0057] According to the characteristic time constant of surface reconstruction phase transition and the critical thickness value, the dislocation nucleation-extension evolution equation is constructed, and a dislocation density evolution prediction model is established. This model links the surface dynamics process (surface reconstruction) with the bulk phase phase transition relaxation process (dislocation formation). In the specific implementation, the relationship function between the dislocation nucleation rate and the surface reconstruction characteristic time constant is first established. The nucleation rate is inversely proportional to the characteristic time constant. The shorter the characteristic time constant, the higher the surface atomic activity and the greater the dislocation nucleation rate. In the GaN / AlN system, when the characteristic time constant is 4.5 seconds, the corresponding dislocation nucleation rate is about 3×10 8 / m 2 ·s.
[0058] A function was established to determine the relationship between the dislocation propagation rate and the critical thickness. The critical thickness reflects the material's ability to withstand strain. A smaller critical thickness indicates a greater tendency for dislocations to form and a higher dislocation propagation rate. At a thickness 1.5 times the critical thickness (approximately 4 nm), the dislocation propagation rate in the GaN / AlN heterostructure is approximately 0.5 μm / s.
[0059] By combining the dislocation nucleation rate and propagation rate, a complete dislocation density evolution prediction model is constructed. This model can predict the evolution trend of dislocation density with growth time and thickness under given growth conditions. In actual verification, for a GaN / AlN heterostructure grown at a temperature of 780°C and a growth rate of 0.7μm / h, the model predicts that the dislocation density near the interface is approximately 5×10 9 / cm 2 , gradually decreases with increasing thickness. When the thickness reaches 500 nm, the dislocation density drops to about 1×10 8 / cm 2 , and the deviation from the transmission electron microscopy (TEM) and cathodoluminescence (CL) measurements is less than 12%.
[0060] The present invention's method for constructing a dislocation density evolution prediction model establishes a correlation between surface dynamics and bulk defect evolution by fusing surface reconstruction data with lattice constant data. This method fully utilizes in-situ monitoring techniques, accurately capturing the dynamic evolution characteristics of the material during growth and accounting for practical factors such as interface step effects and dislocation interactions. Compared to traditional methods that rely on post-characterization, this method enables real-time prediction of dislocation density evolution, providing a scientific basis for optimizing process parameters during heteroepitaxial growth and effectively reducing defect density.
[0061] In an optional embodiment, identifying the current growth stage based on the growth rate data, evaluating the composition uniformity based on the refractive index data, and determining the optimization target based on the growth stage and composition uniformity includes: Calculating the first-order derivative of the growth rate data to obtain the rate change rate, calculating the mean square error of the growth rate data to obtain the rate fluctuation; calculating the gradient of the refractive index data in the depth direction to obtain the component distribution function, and calculating the second-order derivative of the refractive index data in the plane to obtain the component uniformity; constructing a nucleation period discriminant function based on the rate change rate and the rate fluctuation, constructing a stable period discriminant function based on the rate change rate, the rate fluctuation and the component uniformity, and identifying the current growth stage based on the nucleation period discriminant function and the stable period discriminant function; An interface quality control function is constructed based on the square integral of the gradient of the component distribution function and the rate fluctuation during the nucleation period, and a process stability control function is constructed based on the component uniformity and the second-order derivative of the rate change rate during the stabilization period; According to the identified growth stage, the weight coefficients of the interface quality control function and the process stability control function are dynamically adjusted using an exponential decay function, and the weighted control functions are combined to construct an optimization objective function.
[0062] For example, the acquired growth rate data is processed and its first-order derivative is calculated to obtain the rate change rate. Using the central difference method, the rate change rate at time point t is equal to the growth rate at time point t+Δt minus the growth rate at time point t-Δt, divided by 2Δt, where Δt is the sampling interval, typically set to 0.5 seconds. During the initial nucleation stage of GaN, the rate change rate is typically large, reaching 0.1μm / h·min; during the stable growth stage, the rate change rate drops to below 0.01μm / h·min.
[0063] The rate fluctuation is calculated by calculating the mean square error of the growth rate data within a certain time window. The time window length is usually set to 30 seconds. For each time point, the mean square error of the growth rate within the previous 30 seconds is calculated. Specifically, the average of all growth rate data within the time window is calculated, and then the square of the difference between each data point and the average is calculated. Finally, these square values are averaged and the square root is taken to obtain the rate fluctuation. During the nucleation stage of GaN, the rate fluctuation is relatively large, about 0.05-0.15μm / h; during the stable growth stage, the rate fluctuation drops to below 0.02μm / h.
[0064] For the refractive index data, its gradient in the depth direction is calculated to obtain the component distribution function. The distribution of refractive index with depth can be obtained by combining the measurement results of the ellipsometer at different incident angles with the inversion of the optical model. For the refractive index data in the depth direction, a data point is taken every 10nm, and the refractive index difference between adjacent points is calculated and divided by the distance to obtain the local gradient. These local gradients are connected to form the component distribution function. At the high-quality GaN / AlN interface, the gradient value is usually between 2×10 -3 / nm; in samples with poor interface quality, the gradient value can reach 5×10 -3 / nm or above.
[0065] The second-order derivative of the refractive index data in the plane (xy plane) is calculated to obtain the component uniformity. The in-plane refractive index distribution is obtained by scanning the sample at different positions with an ellipsometer. The second-order derivative of the refractive index is calculated for the x-direction and the y-direction respectively, and then the mean square value is taken as the component uniformity index. In the specific implementation, a 5×5 measurement grid is set on the sample with a grid spacing of 5mm, and the finite difference method is applied to the refractive index value of each measurement point to calculate the second-order derivative. For highly uniform GaN films, the component uniformity value is usually less than 5×10 -7 / mm²; for samples with poor uniformity, this value can reach 2×10 -6 / mm 2 above.
[0066] A nucleation phase discriminant function is constructed based on the rate change rate and rate fluctuation. This function takes the form of a weighted sum: the sum of the rate change rate multiplied by a weight coefficient of 0.7 and the rate fluctuation multiplied by a weight coefficient of 0.3. When the discriminant function value is greater than the preset threshold of 0.08, the current nucleation phase is determined. For example, when the rate change rate is 0.09 μm / h·min and the rate fluctuation is 0.11 μm / h, the nucleation phase discriminant function value is 0.096, which is greater than the threshold of 0.08 and is therefore determined to be the nucleation phase.
[0067] At the same time, a stable period discriminant function is constructed based on the rate change rate, rate fluctuation and component uniformity. This function also adopts the weighted sum form, with the weight of the rate change rate being 0.2 (after negation), the weight of the rate fluctuation being 0.3 (after negation), and the weight of the component uniformity being 0.5 (after negation and normalization). When the value of the stable period discriminant function is greater than the preset threshold of 0.85, it is determined that the current period is stable. For example, when the rate change rate is 0.008 μm / h·min, the rate fluctuation is 0.015 μm / h, and the component uniformity is 3×10 -7 / mm 2 When , the discriminant function value of the stable period is 0.91, which is greater than the threshold value of 0.85, and it is determined to be a stable period.
[0068] Based on the calculation results of the nucleation phase discriminant function and the stability phase discriminant function, the current growth stage is identified. If the nucleation phase discriminant function value is greater than its threshold and the stability phase discriminant function value is less than its threshold, the nucleation phase is determined; if the stability phase discriminant function value is greater than its threshold and the nucleation phase discriminant function value is less than its threshold, the stability phase is determined; if both are less than their respective thresholds, the transition phase is determined; if both are greater than their respective thresholds (which rarely occurs), the nucleation phase is prioritized.
[0069] Different control functions are constructed for different growth stages. During the nucleation period, an interface quality control function is constructed based on the square integral of the gradient of the component distribution function and the rate fluctuation. In specific implementation, the square integral of the gradient of the component distribution function in the interface area (usually 20nm thickness) is calculated, that is, the square of the gradient value of each point is summed within the interface thickness range. The interface quality control function is constructed by multiplying the integral value by the weight coefficient 0.65 and the rate fluctuation by the weight coefficient 0.35. The smaller the function value, the higher the interface quality. For high-quality GaN / AlN interfaces, the function value is usually less than 0.04.
[0070] During the stable period, a process stability control function is constructed based on the composition uniformity and the second-order derivative of the rate change rate. The second-order derivative of the rate change rate reflects the acceleration of the growth rate change and is obtained by recalculating the derivative of the rate change rate data. The process stability control function is formed by multiplying the composition uniformity by a weighting factor of 0.6 and the absolute value of the second-order derivative of the rate change rate by a weighting factor of 0.4. The smaller the value of this function, the more stable the growth process. For stably grown GaN films, this function value is typically less than 0.02.
[0071] Based on the identified growth stage, an exponential decay function is used to dynamically adjust the weight coefficients of the interface quality control function and the process stability control function. The exponential decay function is expressed as a weight coefficient equal to the initial weight multiplied by an exponential decay term. The exponential decay term is equal to e-t divided by the power of τ, where t is the growth time and τ is the characteristic decay time, typically set to a parameter related to the epitaxial layer thickness. For example, for the GaN / AlN system, τ is set to 50 seconds multiplied by the growth rate (μm / h).
[0072] In the initial nucleation phase, the weight coefficient of the interface quality control function is set to 0.9, and the weight coefficient of the process stability control function is set to 0.1. As growth progresses, the weight coefficient of the interface quality control function gradually decreases, while the weight coefficient of the process stability control function gradually increases. During the stable period, the weight coefficient of the interface quality control function drops below 0.2, while the weight coefficient of the process stability control function rises to above 0.8.
[0073] The weighted interface quality control function and process stability control function are combined to construct an optimization objective function. This optimization objective function guides the direction and magnitude of process parameter adjustments, achieving precise control over the heteroepitaxial growth process. For example, in the GaN / AlN heteroepitaxial process, when the nucleation phase is detected and the interface quality control function value is 0.06 (high), the optimization objective function will instruct the system to reduce the chamber pressure by 5% and increase the substrate temperature by 10°C to improve the interface quality.
[0074] The present invention's method for identifying growth stages and determining optimization targets based on growth rate and refractive index data enables intelligent monitoring and precise control of wide-bandgap semiconductor heteroepitaxial growth processes. By quantitatively analyzing the changing characteristics of growth rate and refractive index data, the method accurately identifies different growth stages, such as the nucleation phase and the stabilization phase, and constructs corresponding control functions for each phase, enabling dynamic adjustment of the optimization target.
[0075] In an optional embodiment, feeding back actual process parameters and monitoring data to the digital twin model, and updating the parameter coupling influence matrix and dislocation density evolution prediction model includes: Establishing a deviation evaluation function between actual process parameters and theoretically predicted parameters, and updating the strain energy density function module and dislocation motion dynamics equation module in the digital twin model based on the output result of the deviation evaluation function; Calculating a cross-correlation function between the actually monitored stress wave signal and the predicted stress wave signal, and updating a sensitivity coefficient in the parameter coupling influence matrix based on peak and delay characteristics of the cross-correlation function; The deviation between the actually observed surface reconstruction period and the predicted period is substituted into the correction equation of the reconstruction response function, and the deviation between the actually measured critical thickness and the predicted critical thickness is substituted into the correction equation of the dislocation nucleation criterion function. The evolution equation parameters in the dislocation density evolution prediction model are updated based on the corrected function.
[0076] For example, a deviation evaluation function was established to compare actual process parameters with theoretically predicted parameters. The actual process parameters controlled in the MOCVD equipment, such as temperature, gas flow rate, chamber pressure, and RF power, were compared with the optimal parameters predicted by the digital twin model. For temperature, the acceptable deviation range was set to ±5°C; for gas flow, the acceptable deviation range was ±3%; for chamber pressure, the acceptable deviation range was ±0.05 Torr; and for RF power, the acceptable deviation range was ±2%.
[0077] The relative deviation of each parameter is calculated, that is, the actual parameter value minus the predicted parameter value divided by the predicted parameter value. For example, when the predicted temperature is 850°C and the actual temperature is 856°C, the relative deviation is +0.7%. The deviation values of each parameter are then weighted and summed to form a deviation evaluation function. The weight coefficient of each parameter is determined according to its importance to material growth: the temperature parameter weight is 0.4, the gas flow parameter weight is 0.3, the cavity pressure parameter weight is 0.2, and the RF power parameter weight is 0.1. During the GaN / AlN heteroepitaxial growth process, if the calculated deviation evaluation function value is greater than 0.05, the model update mechanism is triggered.
[0078] Based on the output of the deviation evaluation function, the strain energy density function module and the dislocation motion dynamics equation module in the digital twin model are updated. For the strain energy density function module, the elastic constant and thermal expansion coefficient are adjusted according to the deviation of the temperature parameters and the cavity pressure parameters. For example, when the actual temperature value is higher than the predicted value, the elastic constant decreases by about 0.02% / °C and the thermal expansion coefficient increases by about 0.01% / °C. For GaN, the typical elastic constant C11 is 390GPa, C12 is 145GPa, C13 is 106GPa, and C33 is 398GPa; the thermal expansion coefficient in the a direction is about 5.6×10 -6 / ℃, and about 3.2×10 -6 / ℃.
[0079] For the dislocation motion dynamics equation module, the dislocation slip activation energy and dislocation multiplication coefficient are adjusted according to the deviation of temperature parameters and RF power parameters. When the actual temperature value is 10°C higher than the predicted value, the dislocation slip activation energy decreases by about 2% and the dislocation multiplication coefficient increases by about 5%. In GaN, the typical dislocation slip activation energy is 1.5-2.5eV and the dislocation multiplication coefficient is about 5×10 4 / m. Through the dynamic adjustment of these parameters, the digital twin model can more accurately reflect the actual growth environment.
[0080] The cross-correlation function between the actual monitored stress wave signal and the predicted stress wave signal is calculated to update the sensitivity coefficients in the parameter coupling influence matrix. The actual monitored stress wave signal is collected using an acoustic emission sensor with a sampling frequency of 20 kHz and a signal sampling length of 2048 points. The predicted stress wave signal is calculated based on the stress field distribution equation in the digital twin model.
[0081] When calculating the cross-correlation function, the two signals are first normalized to eliminate the influence of amplitude differences. The correlation between the two signals at different time delays is then calculated to produce a cross-correlation curve. Two key features are extracted from this curve: the peak value and the delay feature. The peak value represents the maximum correlation between the two signals, reflecting the similarity between the predicted signal and the actual signal. The delay feature indicates the time delay when the maximum correlation is achieved, reflecting the time response accuracy of the prediction model.
[0082] During GaN / AlN heteroepitaxial growth, if the peak value of the cross-correlation function is less than 0.75 or the delay characteristic exceeds 50ms, it indicates that the parameter coupling influence matrix needs to be updated. This update involves adjusting the overall sensitivity coefficient ratio based on the peak value, with the overall sensitivity coefficient increasing by 5% for every 0.1 decrease in the peak value. The off-diagonal elements of the sensitivity coefficient matrix are also adjusted based on the delay characteristic, with the off-diagonal elements increasing by 3% for every 10ms increase in the delay characteristic to enhance the coupling strength between the parameters. For example, if the coupling coefficient between temperature and gas flow parameters is originally 0.35, but the peak value is 0.65 and the delay characteristic is 70ms, the updated coupling coefficient is 0.35 × 1.05 × 1.06 = 0.39.
[0083] Finally, the dislocation density evolution prediction model is updated by comparing the deviation between the observed and predicted values. First, the deviation between the observed and predicted surface reconstruction period is substituted into the correction equation of the reconstruction response function. The relative deviation is calculated by comparing the surface reconstruction period observed by the RHEED system with the model prediction. For example, for a GaN surface reconstruction from 1×1 to 2×2, if the predicted period is 6.5Å and the observed period is 6.3Å, the relative deviation is -3.1%.
[0084] Based on this deviation, key parameters in the reconstructed response function are modified. Specifically, if the relative deviation is negative (actual value is less than predicted), the surface diffusion barrier parameter in the reconstructed response function is increased; if the relative deviation is positive (actual value is greater than predicted), the parameter is decreased. On the GaN surface, the typical surface diffusion barrier is approximately 0.7-1.2 eV, and every 1% periodic deviation results in an approximately 0.5% adjustment to the barrier value.
[0085] Similarly, the deviation between the actual measured critical thickness and the predicted critical thickness is substituted into the correction equation for the dislocation nucleation criterion function. The critical thickness can be measured using transmission electron microscopy (TEM) or X-ray diffraction (XRD). For example, if the predicted critical thickness is 2.7 nm and the actual measured value is 3.0 nm, the relative deviation is +11.1%.
[0086] Based on this deviation, the strain energy balance parameter in the dislocation nucleation criterion function is modified. Specifically, if the relative deviation is positive (actual value is greater than predicted), the dislocation nucleation energy parameter is increased; if the relative deviation is negative (actual value is less than predicted), the parameter is decreased. For GaN / AlN heterostructures, the typical dislocation nucleation energy is approximately 2.5-3.5 eV / nm, and every 5% deviation in critical thickness results in a roughly 3% adjustment in the nucleation energy value.
[0087] Based on the revised reconstruction response function and dislocation nucleation criterion function, the evolution equation parameters in the dislocation density evolution prediction model are updated. Specific adjustments include the dislocation nucleation rate coefficient and the dislocation expansion rate coefficient. When the surface diffusion barrier of the reconstruction response function increases, the dislocation nucleation rate coefficient decreases accordingly; when the nucleation energy of the dislocation nucleation criterion function increases, the dislocation expansion rate coefficient decreases accordingly. In the GaN / AlN system, the typical dislocation nucleation rate coefficient is about 4×10 8 / cm 2 ·s, and the dislocation growth rate coefficient is about 0.6μm / s.
[0088] By dynamically updating these parameters, the dislocation density evolution prediction model can more accurately reflect the defect evolution patterns during actual material growth. In actual verification, the error between the dislocation density predicted by the updated model and the direct measurement results was reduced to within 10%, significantly improving prediction accuracy.
[0089] This method establishes a mechanism for assessing the deviation between actual process parameters and theoretically predicted parameters, combining stress wave signal cross-correlation analysis with surface reconstruction period and critical thickness measurements to dynamically update the digital twin model, parameter coupling influence matrix, and dislocation density evolution prediction model. This method feeds real-time monitoring data from the material growth process back into the theoretical model, forming a closed-loop control system that continuously approximates the model parameters to the actual physical process, significantly improving the accuracy and adaptability of heteroepitaxial layer growth process parameter adjustment.
[0090] A second aspect of an embodiment of the present invention provides a system for adjusting process parameters for growing a wide bandgap semiconductor heteroepitaxial layer, comprising: The first unit is used to obtain the surface morphology parameters, lattice parameters, and surface defect density of the wide bandgap semiconductor substrate to be processed, build a digital twin model to predict the strain evolution trend in the early stage of heteroepitaxial growth, and determine the optimal initial process parameter set corresponding to the minimum strain and defect density as the starting process point for heteroepitaxial growth; The second unit is used to obtain growth rate data and refractive index data using in-situ ellipsometers during the growth of heteroepitaxial layers, obtain surface reconstruction data and lattice constant data using reflection high-energy electron diffraction, and obtain stress wave signal data using acoustic emission sensors. The obtained multi-source data are input into a deep neural network to calculate the strain tensor field and defect density distribution; A third unit is configured to perform hierarchical and progressive process parameter optimization based on the strain tensor field and defect density distribution, constructing a parameter coupling influence matrix using stress wave signal data in the first layer, establishing a dislocation density evolution prediction model for defect evolution using base surface reconstruction data and lattice constant data in the second layer, and determining a process parameter adjustment strategy based on growth rate data and refractive index data in the third layer; The fourth unit is used to adjust the growth process according to the process parameters after layered optimization, feed back the actual process parameters and monitoring data to the digital twin model, and update the parameter coupling influence matrix and dislocation density evolution prediction model.
[0091] According to a third aspect of an embodiment of the present invention, an electronic device is provided, including: processor; a memory for storing processor-executable instructions; The processor is configured to call the instructions stored in the memory to execute the aforementioned method.
[0092] According to a fourth aspect of an embodiment of the present invention, 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 method described above is implemented.
[0093] 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.
[0094] 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 above embodiments, those skilled in the art should understand that they can still modify the technical solutions described in the above 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. A method for adjusting process parameters for growing a wide bandgap semiconductor heteroepitaxial layer, characterized in that: include: Obtain the surface morphology parameters, lattice parameters, and surface defect density of the wide bandgap semiconductor substrate to be processed, build a digital twin model to predict the strain evolution trend in the early stages of heteroepitaxial growth, and determine the optimal initial process parameter set corresponding to the minimum strain and defect density as the starting process point for heteroepitaxial growth; During the heteroepitaxial layer growth process, an in-situ ellipsometer is used to obtain growth rate data and refractive index data, reflection high-energy electron diffraction is used to obtain surface reconstruction data and lattice constant data, and an acoustic emission sensor is used to obtain stress wave signal data. The multi-source data obtained are input into a deep neural network to calculate the strain tensor field and defect density distribution; Based on the strain tensor field and defect density distribution, a hierarchical and progressive process parameter optimization is performed. In the first layer, the stress wave signal data is used to construct a parameter coupling influence matrix. In the second layer, the base surface reconstruction data and lattice constant data are used to establish a dislocation density evolution prediction model for defect evolution. In the third layer, the process parameter adjustment strategy is determined based on the growth rate data and refractive index data. The growth process is adjusted according to the process parameters after layered optimization, the actual process parameters and monitoring data are fed back to the digital twin model, and the parameter coupling influence matrix and dislocation density evolution prediction model are updated.
2. The method according to claim 1, characterized in that A digital twin model was constructed to predict the strain evolution trend in the early stages of heteroepitaxial growth and determine the optimal initial process parameter set, including: The digital twin model includes a strain energy density function module, a dislocation motion dynamics equation module, and a critical layer thickness calculation module. The surface morphology parameters, lattice parameters, and surface defect density are input into the digital twin model as initial conditions and constraints. The digital twin model constructs a coupled simulation environment of the temperature field distribution equation, the airflow field distribution equation, and the stress field distribution equation, and calculates a parameter coupling coefficient matrix of the temperature parameter, the gas flow parameter, the cavity pressure parameter, and the radio frequency power parameter based on the coupled simulation environment; The parameter coupling coefficient matrix is input into a nonlinear response optimization module, and the nonlinear response optimization module calculates the lattice strain value based on the strain energy density function module and calculates the defect density value based on the dislocation motion dynamics equation module; with the critical layer thickness value as a constraint condition, a weighted objective function is constructed with the lattice strain value and the defect density value, and a constrained gradient descent algorithm is used to iteratively optimize the weighted objective function until a minimum value is obtained, and an optimal initial process parameter group is output, wherein the optimal initial process parameter group includes temperature parameters, gas flow parameters, cavity pressure parameters, and radio frequency power parameters.
3. The method according to claim 1, characterized in that The acquired multi-source data is input into the deep neural network to calculate the strain tensor field and defect density distribution including: Constructing a multi-scale feature extraction network, the multi-scale feature extraction network including a time-scale attention module and a space-scale convolution module, and inputting the multi-source data into the multi-scale feature extraction network to obtain a fused feature vector; Constructing a long short-term memory unit network, inputting the fused feature vector into the long short-term memory unit network, and the long short-term memory unit network calculating a temporal dependency feature based on historical state information and current input information; Constructing a strain prediction network and a defect density prediction network, inputting the timing-dependent features into the strain prediction network to obtain a strain tensor field, wherein the strain tensor field satisfies a continuity equation constraint; inputting the timing-dependent features and the strain tensor field into the defect density prediction network to obtain a defect density distribution, wherein the defect density distribution satisfies an energy conservation constraint; A multi-task joint loss function is constructed based on the strain tensor field and defect density distribution. The multi-task joint loss function includes a strain prediction loss term, a defect density prediction loss term and a regularization loss term. An adaptive learning rate is used to optimize the training of the multi-scale feature extraction network, the long short-term memory unit network, the strain prediction network and the defect density prediction network.
4. The method according to claim 1, wherein Based on the strain tensor field and defect density distribution, a hierarchical and progressive process parameter optimization is performed. In the first layer, the stress wave signal data is used to construct a parameter coupling influence matrix. In the second layer, the base surface reconstruction data and lattice constant data are used to establish a dislocation density evolution prediction model for defect evolution. In the third layer, the process parameter adjustment strategy is determined based on the growth rate data and refractive index data. The following are included: By calculating the time-frequency spectrum and energy transfer function of the stress wave signal data, the sensitivity coefficients between the process parameters are obtained, and a parameter coupling influence matrix is established; Establishing a reconstruction response function using surface reconstruction data, calculating a surface reconstruction phase transition characteristic time constant based on the reconstruction response function, calculating a critical thickness value based on lattice constant data, and constructing a dislocation density evolution prediction model based on the surface reconstruction phase transition characteristic time constant and the critical thickness value; The current growth stage is identified based on the growth rate data, the component uniformity is evaluated based on the refractive index data, and the optimization target is determined based on the growth stage and the component uniformity; the adjustment direction of the process parameters is determined based on the optimization target and the parameter coupling influence matrix, and the adjustment range of the process parameters is determined based on the dislocation density evolution prediction model to obtain the adjustment amount of each process parameter.
5. The method according to claim 4, characterized in that Constructing a dislocation density evolution prediction model includes: The surface reconstruction data is Fourier transformed to obtain the reconstruction period characteristics, and the reconstruction response function is established in combination with the reconstruction intensity change trend; the phase change inflection point is determined based on the first-order derivative and the second-order derivative of the reconstruction response function, and the time interval between adjacent inflection points is calculated to obtain the surface reconstruction phase change characteristic time constant; Substituting the lattice constant data into the elastic strain energy equation and the plastic strain energy equation, respectively, the total strain energy distribution of the heterogeneous interface is calculated, and a dislocation nucleation criterion function is established; solving the strain energy balance equation based on the dislocation nucleation criterion function to obtain the critical layer thickness, and considering the influence of the interface step effect and dislocation interaction on the critical layer thickness to obtain the critical thickness value; A dislocation nucleation-extension evolution equation is constructed based on the surface reconstruction phase transition characteristic time constant and the critical thickness value, and a dislocation density evolution prediction model is established.
6. The method according to claim 4, characterized in that Identifying the current growth stage based on the growth rate data, evaluating the composition uniformity based on the refractive index data, and determining the optimization target based on the growth stage and composition uniformity include: Calculating the first-order derivative of the growth rate data to obtain the rate change rate, calculating the mean square error of the growth rate data to obtain the rate fluctuation; calculating the gradient of the refractive index data in the depth direction to obtain the component distribution function, and calculating the second-order derivative of the refractive index data in the plane to obtain the component uniformity; constructing a nucleation period discriminant function based on the rate change rate and the rate fluctuation, constructing a stable period discriminant function based on the rate change rate, the rate fluctuation and the component uniformity, and identifying the current growth stage based on the nucleation period discriminant function and the stable period discriminant function; An interface quality control function is constructed based on the square integral of the gradient of the component distribution function and the rate fluctuation during the nucleation period, and a process stability control function is constructed based on the component uniformity and the second-order derivative of the rate change rate during the stabilization period; According to the identified growth stage, the weight coefficients of the interface quality control function and the process stability control function are dynamically adjusted using an exponential decay function, and the weighted control functions are combined to construct an optimization objective function.
7. The method according to claim 5, characterized in that Feeding back actual process parameters and monitoring data to the digital twin model and updating the parameter coupling influence matrix and dislocation density evolution prediction model include: Establishing a deviation evaluation function between actual process parameters and theoretically predicted parameters, and updating the strain energy density function module and dislocation motion dynamics equation module in the digital twin model based on the output result of the deviation evaluation function; Calculating a cross-correlation function between the actually monitored stress wave signal and the predicted stress wave signal, and updating a sensitivity coefficient in the parameter coupling influence matrix based on peak and delay characteristics of the cross-correlation function; The deviation between the actually observed surface reconstruction period and the predicted period is substituted into the correction equation of the reconstruction response function, and the deviation between the actually measured critical thickness and the predicted critical thickness is substituted into the correction equation of the dislocation nucleation criterion function. The evolution equation parameters in the dislocation density evolution prediction model are updated based on the corrected function.
8. A wide bandgap semiconductor heteroepitaxial layer growth process parameter adjustment system, used to implement the method according to any one of claims 1 to 7, characterized in that: include: The first unit is used to obtain the surface morphology parameters, lattice parameters, and surface defect density of the wide bandgap semiconductor substrate to be processed, build a digital twin model to predict the strain evolution trend in the early stage of heteroepitaxial growth, and determine the optimal initial process parameter set corresponding to the minimum strain and defect density as the starting process point for heteroepitaxial growth; The second unit is used to obtain growth rate data and refractive index data using in-situ ellipsometers during the growth of heteroepitaxial layers, obtain surface reconstruction data and lattice constant data using reflection high-energy electron diffraction, and obtain stress wave signal data using acoustic emission sensors. The obtained multi-source data are input into a deep neural network to calculate the strain tensor field and defect density distribution; A third unit is configured to perform hierarchical and progressive process parameter optimization based on the strain tensor field and defect density distribution, constructing a parameter coupling influence matrix using stress wave signal data in the first layer, establishing a dislocation density evolution prediction model for defect evolution using base surface reconstruction data and lattice constant data in the second layer, and determining a process parameter adjustment strategy based on growth rate data and refractive index data in the third layer; The fourth unit is used to adjust the growth process according to the process parameters after layered optimization, feed back the actual process parameters and monitoring data to the digital twin model, and update the parameter coupling influence matrix and dislocation density evolution prediction model.
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 according to 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.
Citation Information
Cited By
Accurate control method of VCSEL chip oxidation pore channel
CN121332291A
Single crystal furnace crystal pulling ending parameter control method and system based on self-adaptive regulation and control
CN121381156A
Method and system for producing high-resistance epitaxial wafer with light boron-doped substrate
CN121451286A