A method and system for early warning of a laser crystal manufacturing process

Abstract

A kind of early warning method and system of laser crystal preparation process, including obtaining multimodal sensor data, time alignment and preprocessing are carried out, the multi-dimensional features of each modal data are extracted respectively, and the features of each modal in the same time window are fused, to generate fusion feature sequence;The fusion feature sequence is input to the light-weight time series prediction network, multi-scale time series features are extracted, and attention mask matrix is generated based on the preset time delay prior knowledge, the multi-scale time series features are weighted and fused, the time series features after weighting and fusion are globally averaged to obtain context vector, respectively input to prediction decoder, abnormal risk score binary classification network and fault type multi-classification network, respectively generate the prediction variable sequence of future period, the abnormal risk score of current time and the probability of multiple class faults;When abnormal risk score is greater than the first preset threshold, generate early warning signal;When abnormal risk score is greater than the second preset threshold, generate emergency shutdown signal.

Description

Technical Field

[0001] This invention relates to the field of laser crystal fabrication technology, and in particular to an early warning method and system for the laser crystal fabrication process. Background Technology

[0002] The Czochralski growth process of laser crystals (such as Nd:YAG, Ti:Sapphire, etc.) is extremely sensitive to fluctuations in process parameters such as temperature field distribution, pulling rate, rotation speed, and protective atmosphere. Furthermore, the morphology of the solid-liquid interface front cannot be directly observed within the high-temperature, sealed cavity, and its abnormal evolution exhibits strong nonlinearity and multi-physics coupling characteristics. Currently, monitoring methods for the laser crystal growth process mainly rely on alarms for single or limited parameters, such as thermocouple temperature detection, diameter estimation by weighing, and motor current monitoring. These methods have the following shortcomings: First, monitoring only a few parameters makes it difficult to fully reflect the multi-physics coupling state during the growth process. Second, even if some systems collect data from multiple sensors simultaneously, the data is mostly stored independently or analyzed afterward, lacking a mechanism for time alignment and real-time feature fusion of multimodal data. The causal and time delay relationships between different sensors are not effectively utilized, resulting in low model prediction accuracy and a high false alarm rate. Third, there is a lack of parallel probabilistic diagnostic capabilities for various typical crystal growth faults, and fault location relies on manual experience, resulting in low response efficiency.

[0003] In summary, existing technologies suffer from the lack of multimodal data fusion mechanisms and the inability to achieve parallel diagnosis of multiple faults. Summary of the Invention

[0004] To address the aforementioned technical problems, this invention provides an early warning method for the laser crystal fabrication process, comprising the following steps: Acquire multimodal sensor data, including infrared thermal images, lifting shaft vibration signals, crucible weight signals, atmosphere composition data, and process parameters; the process parameters include lifting rate, rotation rate, and heating power; The multimodal sensor data is time-aligned and preprocessed, multidimensional features of each modal data are extracted, and the modal features within the same time window are fused to generate a fused feature sequence. The fused feature sequence is input into a pre-trained lightweight temporal prediction network; Multi-scale temporal features are extracted using the lightweight temporal prediction network, and an attention mask matrix is ​​generated based on preset time delay prior knowledge. The multi-scale temporal features are then weighted and fused, and the weighted and fused temporal features are obtained by global average pooling to obtain a context vector. The context vector is then input into the prediction decoder, the anomaly risk scoring binary classification network, and the fault type multi-class classification network to generate the prediction variable sequence for future time periods, the anomaly risk score for the current time period, and the probabilities of multiple fault types. The multiple fault types include two or more of the following: seed crystal fracture, interface instability, component undercooling, pulling mechanism failure, temperature runaway, and abnormal atmosphere. When the abnormal risk score is greater than a first preset threshold, an early warning signal is generated; when the abnormal risk score is greater than a second preset threshold, an emergency shutdown signal is generated; the second preset threshold is greater than the first preset threshold.

[0005] Optionally, multidimensional features can be extracted for each modality of data, including: The infrared thermal image is subjected to non-uniformity correction and bad pixel compensation, and the average surface temperature of the melt, surface temperature distribution variance, temperature gradient, isotherm profile radius of curvature and symmetry index are calculated. The vibration signal of the lifting shaft is subjected to bandpass filtering, Fourier transform and envelope demodulation to obtain the effective amplitude, vibration kurtosis and characteristic frequency amplitude of the vibration. The crucible weight signal is processed by Kalman filtering, and the crystal diameter and diameter fluctuation rate are calculated based on the filtered crucible weight change rate. Atmospheric parameters are extracted based on atmospheric composition data, including at least oxygen content and water content.

[0006] Optionally, the infrared thermal image is subjected to non-uniformity correction and bad pixel compensation, and the average surface temperature of the melt, surface temperature distribution variance, temperature gradient, isotherm profile radius of curvature and symmetry index are calculated, including: Non-uniformity correction is performed on the infrared thermal image, pixels that deviate from the neighborhood mean by more than a preset multiple of the standard deviation are detected, and the neighborhood median is used for replacement. Noise reduction is performed by Gaussian filtering. Based on the image segmentation algorithm, the region of interest in the infrared thermal image is extracted, which includes the melt surface, the crystal side surface and the solid-liquid interface. Calculate the arithmetic mean of the temperature values ​​of all pixels in the region of interest as the average temperature of the melt surface, and calculate the temperature gradient; The surface temperature distribution variance is calculated based on the average temperature of the melt surface, which is used to characterize the radial temperature uniformity. Extract all pixels in the infrared thermal image whose pixel temperature is equal to a preset temperature threshold, and form isotherm contours; For each pixel grid of a preset size, the state corresponding to the grid vertex is determined according to the preset temperature threshold, and the crossing edge of the isotherm is determined according to the combination of the vertex states, and the position of the isotherm point is calculated. Connect all isothermal points in sequence to form a closed profile curve, and calculate the radius of curvature, which is equal to the reciprocal of the curvature of the closed profile curve. The symmetry index is calculated based on the left and right boundary coordinates and reference width of the isotherm profile. The symmetry index is used to characterize the degree of symmetry of the isotherm profile relative to the central axis.

[0007] Optionally, the vibration signal of the lifting shaft is subjected to bandpass filtering, Fourier transform, and envelope demodulation to obtain the effective amplitude, kurtosis, and characteristic frequency amplitude of the vibration, including: Based on the vibration signal of the lifting shaft, the least squares method is used to fit the polynomial trend term and the polynomial trend is subtracted from the original vibration signal; The DC component and high-frequency electronic noise are filtered out by a Butterworth bandpass filter, which is designed based on the mechanical characteristic frequency of the crystal growth pulling mechanism. A Hanning window is applied to the filtered lifting shaft vibration signal, followed by a fast Fourier transform to calculate the power spectral density and extract the energy distribution of the preset frequency band. The envelope spectrum of the vibration signal of the lifting shaft is extracted using Hilbert transform to identify the modulation frequency components; Based on the filtered time-domain signal, the square root of the mean square of the vibration acceleration at the sampling points and the vibration kurtosis are calculated; the square root of the power spectrum amplitude is taken as the characteristic frequency amplitude.

[0008] Optionally, the crucible weight signal is subjected to Kalman filtering, and the crystal diameter and diameter fluctuation rate are calculated based on the filtered crucible weight change rate, including: Kalman filtering is used to estimate the state of the crucible weight signal and establish a state space model. The state variables of the state space model include crucible weight and crucible weight change rate. The process noise covariance matrix is ​​dynamically adjusted according to the lifting speed and rotation speed to estimate the optimal weight value and the corresponding crucible weight change rate in real time. Using the synchronously acquired vibration signal of the lifting shaft as a reference, the periodic interference component introduced by the vibration of the lifting shaft is removed from the crucible weight signal. The long-term trend of the crucible weight signal after removing periodic interference components was extracted using a low-pass filter; The crystal diameter is calculated in real time based on the physical relationship between the crucible weight change rate and the crystal growth rate. The standard deviation of the crystal diameter sequence is calculated using a sliding window to identify diameter fluctuation anomalies and to calculate the diameter fluctuation rate.

[0009] Optionally, an attention mask matrix is ​​generated based on preset time delay prior knowledge, and the multi-scale temporal features are weighted and fused. The weighted and fused temporal features are then subjected to global average pooling to obtain a context vector, including: Obtain preset time delay prior knowledge, which includes at least a first time delay interval for temperature change to affect crystal diameter change, a second time delay interval for pulling speed change to affect crystal diameter change, a third time delay interval for rotation speed change to affect crystal diameter change, and a fourth time delay interval for weight change rate to affect crystal diameter change. For the time difference between the current moment and each historical moment, determine whether the time difference falls into the time delay interval corresponding to the feature type of the historical moment, and generate a binary attention mask matrix based on the judgment result. The standard attention weight matrix is ​​calculated based on the multi-scale temporal features. The standard attention weight matrix is ​​then multiplied element-wise with the binary attention mask matrix, and the multiplication result is normalized to obtain the constrained attention weights. The multi-scale temporal features are weighted and summed using the constrained attention weights to generate weighted fused temporal features; The weighted and fused temporal features are then subjected to global average pooling along the time dimension to obtain the context vector. After inputting the context vector into the prediction decoder, the sequence of predictor variables for future time periods is generated through the following steps: The context vector is enhanced by a fully connected layer and then initialized as the hidden state of the gated loop unit in the decoder; The last time-step feature of the fused feature sequence is used as the initial input of the gated recurrent unit, which then generates a sequence of predictive variables for a preset number of future steps through autoregression.

[0010] Optionally, the lightweight temporal prediction network is trained using a multi-task joint loss function, which includes: The prediction loss term is used to constrain the mean squared error between the predicted variable and the true value at multiple future time points. Anomaly detection loss term, used to constrain the binary cross-entropy between anomaly risk score and true anomaly label; The classification loss term is used to constrain the multi-class cross-entropy between the multi-class fault probabilities and the true fault labels; The physical constraint loss term is used to penalize deviations of the predicted variables from the values ​​calculated by the physical model of crystal growth kinetics, which includes the time delay relationship between diameter and temperature.

[0011] Optionally, the method further includes a knowledge management step: in response to the completion of handling a warning event or an emergency shutdown event, obtaining the associated data packet corresponding to the warning event or emergency shutdown event and storing it as an abnormal case; The associated data packet includes at least: multimodal sensor data, fused feature sequence, predictive variable sequence, abnormal risk score, fault type probability distribution, and handling result label within a preset time window before and after the occurrence of the warning signal or emergency shutdown signal.

[0012] Optionally, the method further includes an edge-cloud collaboration step performed by the edge node: The multimodal sensor data and fused feature sequences are stored, and stored in directories by date with a preset retention period. At preset time intervals, the multimodal sensor data and statistical summary data of the fused feature sequences are sent to the cloud. In response to the abnormal risk score being greater than the first preset threshold or the second preset threshold, multimodal sensor data, fused feature sequences, and predictive variable sequences within a preset time window before and after the occurrence of the warning signal or emergency shutdown signal are extracted, compressed, and sent to the cloud. The model parameters are received from the cloud and the lightweight time series prediction network is updated using a hot update method.

[0013] Corresponding to the early warning method for the laser crystal fabrication process, the present invention provides an early warning system for the laser crystal fabrication process, comprising: The multimodal sensor data acquisition module is used to acquire multimodal sensor data, including infrared thermal images, lifting shaft vibration signals, crucible weight signals, atmosphere composition data, and process parameters; the process parameters include lifting rate, rotation rate, and heating power. The feature fusion module is used to perform time alignment and preprocessing on the multimodal sensor data, extract multidimensional features of each modal data, fuse the modal features within the same time window to generate a fused feature sequence, and input the fused feature sequence into a pre-trained lightweight temporal prediction network. A lightweight temporal prediction network is used to extract multi-scale temporal features based on the fused feature sequence, and generate an attention mask matrix based on preset time delay prior knowledge. The multi-scale temporal features are weighted and fused, and the weighted and fused temporal features are obtained by global average pooling to obtain a context vector. The context vector is input to the prediction decoder, the anomaly risk scoring binary classification network, and the fault type multi-class classification network, respectively, to generate the prediction variable sequence for future time periods, the anomaly risk score for the current time, and the probability of multiple fault types. The multiple fault types include two or more of the following: seed crystal fracture, interface instability, component undercooling, pulling mechanism failure, temperature runaway, and atmosphere anomaly. The early warning module is used to generate an early warning signal when the abnormal risk score is greater than a first preset threshold; and to generate an emergency shutdown signal when the abnormal risk score is greater than a second preset threshold; wherein the second preset threshold is greater than the first preset threshold.

[0014] Compared with the prior art, the present invention has the following beneficial effects: (1) By simultaneously acquiring infrared thermal images, lifting shaft vibration signals, crucible weight signals, atmosphere composition data, and process parameters, comprehensive monitoring of the multi-physics coupling state during laser crystal preparation was achieved; time alignment and preprocessing were performed on the multi-modal sensor data, and after extracting the multi-dimensional features of each modal data, the modal features within the same time window were fused to generate a fused feature sequence, which promoted the formation of a fused representation with complementary information under a unified temporal framework of different sensor data; the fused feature sequence was input into a pre-trained lightweight temporal prediction network to extract multi-scale temporal features, and an attention mask matrix was generated based on the preset time delay prior knowledge, and the multi-scale temporal features were weighted and fused, and the weighted fused time sequence was then used to generate a fused feature sequence. The sequence features are global average pooled to obtain a context vector, which is then input into the predictive decoder to generate a sequence of predictive variables for future time periods. This allows the model to focus on the causal time-delay relationship between temperature, pulling speed, rotation speed, and weight changes in crystal growth state, thereby improving prediction accuracy and reducing false alarm rate. The context vector is input into a fault type multi-classification network to obtain the probabilities of multiple fault types at the current moment, enabling parallel probabilistic diagnosis of various typical crystal growth faults. The context vector is then input into an anomaly risk scoring binary classification network to obtain an anomaly risk score, which is compared with a first preset threshold and a second preset threshold to generate warning signals and emergency shutdown signals, respectively, enabling rapid graded response to anomalies. The synergistic effect of multimodal data fusion, time-delay prior constraint-based time-series prediction, and parallel probabilistic diagnosis of multiple faults significantly improves the comprehensiveness, prediction accuracy, and response efficiency of anomaly monitoring in the laser crystal fabrication process.

[0015] (2) Non-uniformity correction and bad pixel compensation were performed on the infrared thermal image, and the average surface temperature, surface temperature distribution variance, temperature gradient, isotherm profile radius of curvature and symmetry index of the melt were calculated, realizing the quantitative characterization of the temperature field distribution morphology characteristics, and providing image-level feature support for judging the stability of the solid-liquid interface; Bandpass filtering, Fourier transform and envelope demodulation were performed on the lifting shaft vibration signal to obtain the effective amplitude, vibration kurtosis and characteristic frequency amplitude of the vibration, realizing the effective separation and quantification of the mechanical operation state of the lifting mechanism; Kalman filtering was performed on the crucible weight signal, and the crystal diameter and diameter fluctuation rate were calculated, suppressing the noise in the weight signal and improving the accuracy of diameter estimation; Oxygen content and water content were extracted from the atmosphere composition data, realizing the effective monitoring of atmospheric environment anomalies.

[0016] (3) Non-uniformity correction, bad pixel compensation and Gaussian filtering are performed on infrared thermal images to suppress the inherent noise and bad pixel interference of infrared images and improve the accuracy of melt surface temperature measurement. Image segmentation algorithm is used to extract regions of interest such as melt surface, crystal side and solid-liquid interface, so that temperature field analysis focuses on key areas. The radius of curvature and symmetry index of isotherm contour are calculated to transform the morphological characteristics of solid-liquid interface into quantifiable geometric indicators, so that anomalies such as interface instability and temperature field distortion can be identified in time through the symmetry and curvature changes of isotherm contour, thus providing accurate temperature field morphological characteristics for the diagnosis of various faults.

[0017] (4) The least squares method was used to fit and remove the polynomial trend term, thus eliminating the slow drift in the vibration signal of the lifting shaft; the Butterworth bandpass filter was used to filter out the DC component and high-frequency electronic noise, and the passband was designed according to the mechanical characteristic frequency of the crystal growth lifting mechanism, effectively preserving the characteristic frequency components related to the lifting mechanism fault; the spectrum and power spectral density of the windowed signal and the energy distribution of the preset frequency band were calculated by fast Fourier transform, and the envelope spectrum was extracted by Hilbert transform to identify the modulation frequency components, effectively separating mechanical fault characteristics such as bearing wear, shaft misalignment, and periodic impact; the square root of the square mean of the vibration acceleration at the sampling point and the vibration kurtosis were calculated; the square root of the power spectral amplitude was used as the characteristic frequency amplitude, realizing the multi-dimensional quantitative characterization of the mechanical operating state of the lifting shaft.

[0018] (5) A state-space model containing crucible weight and crucible weight change rate was established by Kalman filtering, and the process noise covariance matrix was dynamically adjusted according to the lifting speed and rotation speed. This achieved adaptive optimization of the filtering effect of the weight signal under dynamic working conditions, making the real-time estimated optimal weight value and change rate more accurate. The periodic interference component was eliminated by using the synchronously acquired lifting shaft vibration signal as a reference, thus eliminating the coupling interference of mechanical vibration on weight measurement. The long-term trend was extracted by low-pass filter, which enabled accurate tracking of the slow weight change of crystal growth. The crystal diameter was calculated in real time according to the physical relationship between the weight change rate and the crystal growth rate, and the standard deviation of the crystal diameter sequence was calculated by sliding window to identify abnormal diameter fluctuations.

[0019] (6) By acquiring the time delay intervals of temperature change, pulling speed change, rotation speed change and weight change rate affecting crystal diameter change as prior knowledge, the model has a structured understanding of the causal time delay relationship of multi-physics coupling during crystal growth; by judging whether the time difference between the current moment and each historical moment falls into the corresponding time delay interval, a binary attention mask matrix is ​​generated, and physical priors are embedded into the attention mechanism in a hard constraint manner, so that the model only focuses on historical moments with physical causal significance when calculating attention weights, and shields the interference of irrelevant time steps; by multiplying the standard attention weight matrix and the binary attention mask matrix element by element and normalizing, the constrained attention weights retain both the data-driven adaptability and conform to the crystal growth dynamics; by global average pooling, fully connected layer enhancement and gated recurrent unit autoregressive decoding, the weighted fusion time sequence features are efficiently transformed into a sequence of predictive variables for the future preset number of steps.

[0020] (7) By constraining the mean square error between the predicted variables and the true values ​​at multiple future time points through the prediction loss term, the accuracy of the time series prediction task is ensured; by constraining the binary cross-entropy between the anomaly risk score and the true anomaly label through the anomaly detection loss term, the accuracy of anomaly discrimination is improved; by constraining the multi-class cross-entropy between the multi-class fault probability and the true fault label through the classification loss term, the reliability of multi-fault parallel diagnosis is ensured; by penalizing the deviation of the predicted variables from the calculated values ​​of the crystal growth kinetic physical model through the physical constraint loss term, especially the time delay relationship between diameter and temperature, the network prediction results not only fit the historical data, but also conform to the physical laws of crystal growth.

[0021] (8) By acquiring the associated data packets corresponding to the early warning event or emergency shutdown event and storing them as abnormal cases, a structured abnormal case library is constructed. This not only preserves the complete multimodal evidence chain when the abnormality occurs, but also associates the handling results, so that subsequent network parameters can be continuously optimized, the latency prior knowledge and fault diagnosis rules can be improved through case retrieval, similarity matching or incremental learning, so as to achieve continuous accumulation of knowledge and iterative improvement of model performance.

[0022] (9) By storing multimodal sensor data and fused feature sequences locally on the edge nodes and retaining them in a date-based directory for a preset duration, the needs of localized data management and rapid backtracking are met; by sending statistical summary data to the cloud at preset time intervals, the network bandwidth is efficiently utilized while enabling the cloud to grasp the overall operating status; by responding to an abnormal risk score exceeding the threshold, detailed data within a preset time window before and after the occurrence of a warning signal or emergency shutdown signal is compressed and sent to the cloud, achieving accurate cloud uploading of key data of abnormal events and avoiding bandwidth pressure caused by continuous uploading of full data; by receiving the updated model parameters from the cloud and updating the lightweight time series prediction network in a hot update manner, the edge nodes can continuously obtain optimized model capabilities without interrupting the growth process. Attached Figure Description

[0023] The accompanying drawings, which are included to provide a further understanding of the invention and form part of this invention, illustrate exemplary embodiments of the invention and are used to explain the invention, but do not constitute an undue limitation of the invention. In the drawings: Figure 1 This is a simplified flowchart of an embodiment of the early warning method for the laser crystal preparation process of the present invention; Figure 2 This is a framework diagram of an embodiment of the early warning system for the laser crystal preparation process of the present invention. Detailed Implementation

[0024] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, not all embodiments. It should be understood that the specific embodiments described herein are only used to explain the present invention and are not intended to limit the present invention. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0025] like Figure 1 As shown, the present invention provides an early warning method for the laser crystal fabrication process, which includes the following steps: Acquire multimodal sensor data, including infrared thermal images, lifting shaft vibration signals, crucible weight signals, atmosphere composition data, and process parameters; the process parameters include lifting rate, rotation rate, and heating power; Multimodal sensor data is time-aligned and preprocessed, multidimensional features of each modality are extracted, and the modal features within the same time window are fused to generate a fused feature sequence. The fused feature sequence is input into a pre-trained lightweight temporal prediction network, which includes a shared encoder, a prediction decoder (GRU), an anomaly risk scoring binary classification network (MLP+Sigmoid), and a fault type multi-classification network (Softmax). Multi-scale temporal features are extracted using a lightweight temporal prediction network, and an attention mask matrix is ​​generated based on pre-defined time delay prior knowledge. The multi-scale temporal features are then weighted and fused, and the weighted and fused temporal features are obtained by global average pooling to obtain a context vector. The context vector is then input into a prediction decoder, an anomaly risk scoring binary classification network, and a fault type multi-class classification network, respectively, to generate a sequence of predictive variables for future time periods, an anomaly risk score for the current time period, and probabilities of multiple fault types. The multiple fault types include two or more of the following: seed crystal fracture, interface instability, component undercooling, pulling mechanism failure, temperature runaway, and abnormal atmosphere. When the abnormal risk score exceeds a first preset threshold, an early warning signal is generated; when the abnormal risk score exceeds a second preset threshold, an emergency shutdown signal is generated; the second preset threshold is greater than the first preset threshold. Preferably, the second preset threshold is 0.9 and the first preset threshold is 0.7.

[0026] This invention achieves comprehensive monitoring of the laser crystal fabrication process by simultaneously acquiring infrared thermal images, lifting shaft vibration signals, crucible weight signals, atmosphere composition data, and process parameters. It performs time alignment and preprocessing on multimodal sensor data, extracting multidimensional features from each modality, and then fuses the features within the same time window to generate a fused feature sequence. This promotes the formation of a fused representation with complementary information from different sensor data within a unified temporal framework. The fused feature sequence is input into a pre-trained lightweight temporal prediction network to extract multi-scale temporal features. Based on preset time delay prior knowledge, an attention mask matrix is ​​generated to weightedly fuse the multi-scale temporal features. The weighted fused temporal features are then processed by a full-processing network. Local average pooling yields a context vector, which is input into the predictive decoder to generate a sequence of predictive variables for future time periods. This allows the model to focus on the causal time-delay relationship between temperature, pulling speed, rotation speed, and weight changes in crystal growth state, thereby improving prediction accuracy and reducing false alarm rate. The context vector is then input into a fault type multi-classification network to obtain the probabilities of various fault types at the current moment, enabling parallel probabilistic diagnosis of multiple typical crystal growth faults. Finally, the context vector is input into an anomaly risk scoring binary classification network to obtain an anomaly risk score, which is then compared with a first preset threshold and a second preset threshold to generate warning signals and emergency shutdown signals, respectively, achieving rapid graded response to anomalies. The synergistic effect of multimodal data fusion, time-delay prior constraint-based temporal prediction, and parallel probabilistic diagnosis of multiple faults significantly improves the comprehensiveness, prediction accuracy, and response efficiency of anomaly monitoring in the laser crystal fabrication process.

[0027] In this embodiment, edge nodes synchronize with a global synchronization clock server via NTP (Network Time Protocol) or PTP (Precision Time Protocol) to ensure that all sensor data have a unified time reference system with a synchronization accuracy better than 1ms. Corresponding driver acquisition programs are written for different sensor interface types (GigE Vision, RS485, IEPE, EtherCAT, analog signals, etc.) to uniformly convert the data into an internal standard format. Each data frame or sampling point is immediately stamped with a global synchronization timestamp (accuracy 0.1ms) at the time of acquisition and stored in a local circular buffer. For image data such as infrared thermal imaging, the frame arrival time is used as the timestamp; for continuously sampled vibration signals, the start time of the sampling sequence is used as the reference.

[0028] Preferably, the preprocessing includes outlier detection and removal, missing value processing, and incremental preprocessing of heating power (the control output generated by the closed-loop calculation of the temperature control system uses an incremental PID controller, which directly calculates the increment of power output).

[0029] In this embodiment, outlier detection and removal can be achieved in the following ways: Statistical thresholding: For slowly changing parameters such as weighing and temperature, calculate the mean and standard deviation within a sliding window (window length 60 seconds) and remove outliers that exceed the mean ± 3 times the standard deviation.

[0030] Physical constraint method: Set hard thresholds based on sensor range and physical possibilities. For example, data is directly judged as invalid when the crucible temperature exceeds 2100℃ or is below room temperature; vibration acceleration amplitude is marked as abnormal when it exceeds the sensor range ±10g.

[0031] Temporal continuity detection: Check if there is a long-term interruption in the data stream (more than twice the preset sampling interval). Data during the interruption period is marked as missing values ​​and interpolated later.

[0032] Spike removal: Median filtering (window length 5-11 points) is used to identify and remove single-point spike noise, and the original value is replaced with the median value within the window.

[0033] In this embodiment, missing value handling can be implemented in the following way: Missing type identification: distinguish between random missing (single point of loss) and continuous missing (packet loss caused by communication interruption).

[0034] Linear interpolation: For cases with random missing values ​​and a missing duration of less than 5 sampling periods, linear interpolation of the effective values ​​before and after the missing values ​​is used to fill the missing value.

[0035] Spline interpolation: For parameters with strong trends, such as temperature and weight, cubic spline interpolation is used to ensure the smoothness of the filled values ​​and the continuity of the first derivative.

[0036] Model prediction imputation: For cases where more than 10 consecutive sampling periods are missing, a short-term prediction imputation is performed using an ARIMA model trained on historical data. At the same time, the data segment with low confidence is marked for reference in subsequent analysis.

[0037] Discarding strategy: For data segments with a missing rate exceeding 30%, discard them directly and log the results to avoid introducing excessive errors.

[0038] In this embodiment, time alignment can be achieved in the following way: Interpolation Alignment: Using the highest sampling rate (vibration 10kHz) as the baseline, cubic spline interpolation is used to upsample the low-frequency data (weighing 10Hz, temperature 1Hz) to the same frequency as the vibration. Interpolation function formula: ;in, Indicates the time point where interpolation is needed. Indicates the first Sampling timestamp, For a moment The measured data values, For a moment rate of change at that point For a moment acceleration at that point The rate of change of curvature within the interval. .

[0039] Sliding window alignment: For data windows that need to be analyzed synchronously, data segments from each sensor within the window are extracted with the timestamp as the center to form an aligned multimodal data matrix.

[0040] Event-triggered alignment: When an abnormal event (such as a sudden change in diameter) is detected, all sensor data from 15 seconds before and after the event are automatically extracted to ensure the synchronization of multimodal data in event analysis.

[0041] In this embodiment, multidimensional features of each modality data are extracted, including: Non-uniformity correction and bad pixel compensation are performed on infrared thermal images, and the average surface temperature of the melt, surface temperature distribution variance, temperature gradient, isotherm profile radius of curvature and symmetry index are calculated. The vibration signal of the lifting shaft is subjected to bandpass filtering, Fourier transform and envelope demodulation to obtain the effective amplitude, vibration kurtosis and characteristic frequency amplitude of the vibration. The crucible weight signal is processed by Kalman filtering, and the crystal diameter and diameter fluctuation rate are calculated based on the filtered crucible weight change rate. Atmospheric parameters are extracted based on atmospheric composition data, including at least oxygen content and water content.

[0042] This invention performs non-uniformity correction and bad pixel compensation on infrared thermal images, and calculates the average surface temperature of the melt, surface temperature distribution variance, temperature gradient, isotherm profile radius of curvature, and symmetry index, realizing the quantitative characterization of the temperature field distribution morphology and providing image-level feature support for judging the stability of the solid-liquid interface. Bandpass filtering, Fourier transform, and envelope demodulation are performed on the lifting shaft vibration signal to obtain the effective amplitude, kurtosis, and characteristic frequency amplitude of the vibration, achieving effective separation and quantification of the mechanical operating state of the lifting mechanism. Kalman filtering is applied to the crucible weight signal, and the crystal diameter and diameter fluctuation rate are calculated, suppressing noise in the weight signal and improving the accuracy of diameter estimation. Oxygen and water content are extracted from the atmosphere composition data, enabling effective monitoring of abnormal atmospheric environments.

[0043] In this embodiment, non-uniformity correction and bad pixel compensation are performed on the infrared thermal image, and the average surface temperature of the melt, surface temperature distribution variance, temperature gradient, isotherm profile radius of curvature, and symmetry index are calculated. This includes: performing non-uniformity correction on the infrared thermal image, as shown in the following formula: ;in, Represents the center pixel coordinates after non-uniformity correction. The corresponding temperature value, This represents the x-coordinate of the center pixel of the focal plane array in an infrared thermal imager. This represents the ordinate of the center pixel of the focal plane array in an infrared thermal imager. This is the original temperature matrix, with a size of [number] pixels. , This is the gain correction factor. , Indicates the high-temperature calibration value. Indicates the low-temperature calibration value; , This indicates a high-temperature dark-field bias, where K represents the total number of frames in the acquired dark-field images, and k represents the k-th frame. This represents the high-temperature dark field bias of the k-th frame. ; , Low-temperature dark field bias is used to eliminate stationary noise. This indicates the low-temperature dark field bias of the k-th frame; Pixels that deviate from the neighborhood mean by more than a preset multiple (preferably 3 times) of the standard deviation are detected and replaced with the neighborhood median. Noise reduction is then performed using Gaussian filtering. The specific formula is as follows: ;in, This represents the temperature value corresponding to the center pixel after Gaussian filtering, where u represents the temperature relative to the center pixel. The horizontal offset, v, represents the horizontal offset relative to the center pixel. Vertical offset, =1.5, k=3; Based on image segmentation algorithms, regions of interest are extracted from infrared thermal images, including melt surfaces, crystal sides, and solid-liquid interfaces. Calculate the arithmetic mean of the temperature values ​​of all pixels within the region of interest as the average temperature of the melt surface, and then calculate the temperature gradient. The specific formula is as follows: Where T represents a point on the infrared thermal image; The corresponding temperature value, Represents a point on an infrared thermal image The corresponding horizontal coordinates This represents the vertical coordinates of a point on an infrared thermal image. The variance of surface temperature distribution is calculated based on the average temperature of the melt surface, which is used to characterize radial temperature uniformity. Extract all pixels in the infrared thermal image whose pixel temperature is equal to a preset temperature threshold, and form an isotherm outline; the preset temperature threshold is the difference between the preset temperature threshold and the supercooling. For each preset size (preferred) The pixel grid is determined according to a preset temperature threshold, which determines the state of the four grid vertices. And based on the combinations of vertex states (a total of 16 cases), the edges crossing the isotherm are determined, and the positions of the isotherms are calculated; where: ; Indicates the preset temperature threshold. represents the column index of the pixel grid vertex, and j represents the row index of the pixel grid vertex; P is the location of the isotherm. This is the starting temperature value. Final temperature value; Connect all isothermal points in sequence to form a closed profile curve, and calculate the radius of curvature, which is equal to the reciprocal of the curvature of the closed profile curve. Based on the coordinates of the left and right boundaries of the isotherm profile and the reference width, the symmetry index is calculated. The symmetry index characterizes the degree of symmetry of the isotherm profile relative to the central axis; the specific formula is as follows: ; in, The symmetry index Isotherms at altitude The coordinates of the left boundary at that location, Isotherms at altitude The right boundary coordinates at that location H represents the reference width (e.g., the nominal radius of the crystal at that height), and H represents the height range for analysis. This indicates the highest row index of the isotherm in the infrared thermal image. This indicates the lowest row index of the isotherm in the infrared thermal image. The x-coordinate represents the axis of symmetry and is the vertical center of the image.

[0044] This invention performs non-uniformity correction, dead pixel compensation, and Gaussian filtering on infrared thermal images to suppress inherent noise and dead pixel interference, thereby improving the accuracy of melt surface temperature measurement. It employs image segmentation algorithms to extract regions of interest such as the melt surface, crystal side surfaces, and solid-liquid interfaces, allowing temperature field analysis to focus on key areas. Furthermore, it calculates the radius of curvature and symmetry index of isotherm profiles, transforming the morphological characteristics of the solid-liquid interface into quantifiable geometric indicators. This enables timely identification of anomalies such as interface instability and temperature field distortion through changes in the symmetry and curvature of isotherm profiles, thus providing precise temperature field morphological characteristics for the diagnosis of various faults.

[0045] In this embodiment, the vibration signal of the lifting shaft (preferably, the vibration acceleration signal at the upper end of the lifting shaft is collected, retaining the vibration components in the range of 0.5Hz-5kHz) is subjected to bandpass filtering, Fourier transform, and envelope demodulation to obtain the effective amplitude, kurtosis, and characteristic frequency amplitude of the vibration, including: Based on the vibration signal of the lifting shaft, the least squares method is used to fit the polynomial trend term (order 2-4) and the polynomial trend is subtracted from the original vibration signal to eliminate the influence of sensor temperature drift and low frequency noise. The DC component and high-frequency electronic noise are filtered out by Butterworth bandpass filters (order 4, passband 0.5Hz-5kHz). The Butterworth bandpass filters are designed according to the mechanical characteristic frequencies of the crystal growth pulling mechanism (spindle rotation frequency 5-50Hz, bearing characteristic frequency 50-500Hz). A Hanning window is applied to the filtered lifting shaft vibration signal to suppress spectral leakage. Then, a fast Fourier transform is performed to calculate the spectrum and power spectral density of the windowed signal and extract the energy distribution of the preset frequency band (preferably 50-500Hz). To identify early bearing failure characteristics, Hilbert transform is used to extract the envelope spectrum of the lifting shaft vibration signal and identify the modulation frequency components. Based on the filtered time-domain signal, the square root of the mean square of the vibration acceleration at the sampling points (to reflect the effective amplitude of the vibration) and the vibration kurtosis are calculated; the square root of the power spectrum amplitude is taken as the characteristic frequency amplitude.

[0046] This invention employs least squares fitting and removes polynomial trend terms to eliminate slow drift in the vibration signal of the lifting shaft; it utilizes a Butterworth bandpass filter to remove DC components and high-frequency electronic noise, and designs a passband based on the mechanical characteristic frequency of the crystal growth lifting mechanism, effectively preserving characteristic frequency components related to lifting mechanism faults; it calculates the spectrum and power spectral density of the windowed signal and the energy distribution of the preset frequency band through fast Fourier transform, and extracts the envelope spectrum to identify modulation frequency components through Hilbert transform, effectively separating mechanical fault characteristics such as bearing wear, shaft misalignment, and periodic impact; it calculates the square root of the square mean of the vibration acceleration at the sampling points and the vibration kurtosis; and it uses the square root of the power spectral amplitude as the characteristic frequency amplitude, realizing a multi-dimensional quantitative characterization of the mechanical operating state of the lifting shaft.

[0047] In this embodiment, the crucible weight signal is subjected to Kalman filtering, and the crystal diameter and diameter fluctuation rate are calculated based on the filtered crucible weight change rate, including: Kalman filtering is used to estimate the state of the crucible weight signal, establishing a state-space model. The state variables of the state-space model include the crucible weight and the rate of change of crucible weight. The process noise covariance matrix is ​​dynamically adjusted based on the lifting speed and rotation speed to estimate the optimal weight value and the corresponding rate of change of crucible weight in real time. The filtering formula is as follows: Prediction step: ; Update steps: ;in, for The actual measurement value of the weight sensor at any given time (crucible weight signal). Indicates the discrete time step number. express Control input vectors at any given time (e.g., lifting speed, rotation speed, temperature). yes The Kalman gain at time t, A is the state transition matrix, and B is the control input matrix. For measurement matrix; Using the synchronously acquired lifting shaft vibration signal as a reference, the periodic interference component introduced by the lifting shaft vibration is removed from the crucible weight signal through adaptive noise cancellation technology. To address the weight fluctuations caused by natural convection of the melt (frequency 0.1-1Hz), a low-pass filter is used to extract the long-term trend of the crucible weight signal after removing periodic interference components. The crystal diameter is calculated in real time based on the physical relationship between the crucible weight change rate and the crystal growth rate; the specific formula is as follows: ;in, ρ is the rate of change of crucible weight, ρ is the crystal density, and ν is the pulling rate. The diameter of the crystal; The standard deviation of the crystal diameter sequence is calculated using a sliding window (window length 30 seconds) to identify diameter fluctuation anomalies and calculate the diameter fluctuation rate.

[0048] This invention employs Kalman filtering to establish a state-space model containing the crucible weight and the rate of change of crucible weight. It dynamically adjusts the process noise covariance matrix based on the lifting speed and rotation speed, achieving adaptive optimization of the filtering effect of the weight signal under dynamic conditions, making the real-time estimated optimal weight value and rate of change more accurate. It uses synchronously acquired lifting shaft vibration signals as a reference to eliminate periodic interference components, thus eliminating the coupling interference of mechanical vibration on weight measurement. A low-pass filter is used to extract long-term trends, enabling accurate tracking of slow weight changes during crystal growth. The crystal diameter is calculated in real-time based on the physical relationship between the rate of change of weight and the crystal growth rate, and a sliding window standard deviation is calculated on the crystal diameter sequence to identify abnormal diameter fluctuations.

[0049] In this embodiment, the temporal convolutional network uses a combination of causal convolution and dilated convolution to ensure the unidirectional flow of temporal information and prevent future information leakage.

[0050] In this embodiment, an attention mask matrix is ​​generated based on preset latency prior knowledge. Multi-scale temporal features are then weighted and fused. The weighted and fused temporal features are then subjected to global average pooling to obtain a context vector, including: Acquire the preset time delay prior knowledge, which includes at least the first time delay interval (30-180 seconds) in which temperature changes affect crystal diameter changes, the second time delay interval (5-20 seconds) in which pulling speed changes affect crystal diameter changes, the third time delay interval (10-30 seconds) in which rotation speed changes affect crystal diameter changes, and the fourth time delay interval (0-5 seconds) in which weight change rate affects crystal diameter changes. Given the time difference between the current moment and each historical moment, determine whether the time difference falls within the time delay interval corresponding to the feature type of that historical moment, and generate a binary attention mask matrix based on the determination result; the specific formula is as follows: ; This represents the binary attention mask matrix. This indicates an indicator function that returns 1 when the condition within the parentheses is true and 0 when it is false; t represents the current time. Representing historical moments, This represents the time interval between the current moment and a historical moment, i.e., the lag time of the temperature effect diameter; This indicates the minimum physical response time (preferably 30s, since it takes at least 30 seconds for a temperature change to affect the diameter). This indicates the maximum effective memory time (preferably 180 seconds; temperature information exceeding 180 seconds has negligible impact on the current diameter). The standard attention weight matrix is ​​calculated based on multi-scale temporal features, using the following formula: , ; This represents the degree of attention the model pays to historical moments at the current moment, with a value ranging from [0,1]. , This is the raw attention score, representing the strength of the correlation between the current moment and previous moments. The feature vector output at the current time, The feature vectors output at historical moments. The dimension of the feature vector. For all historical moments, temporary variables; The standard attention weight matrix is ​​multiplied element-wise with the binary attention mask matrix, and the result of the multiplication is normalized to obtain the constrained attention weights. We use constrained attention weights to perform weighted summation on multi-scale temporal features to generate weighted fused temporal features; The weighted and fused temporal features are then subjected to global average pooling along the time dimension to obtain the context vector; After inputting the context vector into the prediction decoder, the sequence of predictor variables for future time periods is generated through the following steps: After the context vector is enhanced by a fully connected layer, it is initialized as the hidden state of the gated recurrent unit in the decoder; The last time-step feature of the fused feature sequence is used as the initial input of the gated recurrent unit, which then generates a sequence of predictive variables for a preset number of future steps through autoregression.

[0051] This invention acquires prior knowledge of the time delay intervals influencing crystal diameter changes due to temperature variations, pulling speed variations, rotation speed variations, and weight change rates, enabling the model to possess a structured understanding of the causal time delay relationships of multi-physics coupling during crystal growth. A binary attention mask matrix is ​​generated by determining whether the time difference between the current moment and each historical moment falls within the corresponding time delay interval. This hard-constraint embeds physical priors into the attention mechanism, ensuring that the model only focuses on historical moments with physical causal significance when calculating attention weights, shielding it from interference from irrelevant time steps. By multiplying the standard attention weight matrix element-wise with the binary attention mask matrix and normalizing, the constrained attention weights retain both data-driven adaptability and conform to crystal growth dynamics. Global average pooling, fully connected layer enhancement, and gated recurrent unit autoregressive decoding efficiently transform the weighted fused temporal features into a sequence of predictive variables for a predetermined number of future steps.

[0052] In this embodiment, after the context vector is input into the anomaly risk scoring binary classification network, the anomaly risk score is generated through the following steps: (1) Fully connected Linear: , For intermediate features of linear output, The weight matrix corresponding to the context vector. For context vectors, It is a 128-dimensional bias vector. Since the context vector contains a variety of information such as temperature field morphology, mechanical vibration state, diameter growth trend, and atmosphere deviation, the first step can establish a special coordinate system for anomaly detection in the feature space, and transform the general growth state vector into "anomaly discrimination intermediate features".

[0053] (2) Nonlinear activation ReLU: , To find the maximum value element by element, These are the features after the first nonlinear activation. Crystal growth anomalies are often caused by multi-factor nonlinear coupling. For example, a simple temperature fluctuation of 0.3°C is not considered an anomaly, but if it is accompanied by a sudden increase in vibration kurtosis and an increase in diameter fluctuation rate, it is highly suspicious. Therefore, in this embodiment, nonlinearity is introduced through ReLU, enabling the network to learn this complex anomaly pattern of "multi-parameter cross-thresholds," rather than simply performing linear weighting.

[0054] (3) Dropout: , , For traversal The index numbers of all vectors, These are the features output after the first Dropout. In laser crystal fabrication, there are slight differences in raw material purity and ambient temperature and humidity between different furnace batches. Dropout forces the network to learn general anomalies that are robust to process fluctuations, rather than overfitting to the specific noise of a particular furnace. Dropout is turned off during the inference phase.

[0055] (4) Fully connected Linear: , The characteristics of the second linear output, for The corresponding weight matrix, This is a 64-dimensional bias vector. This step further compresses the feature dimensions, which is equivalent to focusing on the most suspicious anomaly clues.

[0056] (5) Nonlinear activation ReLU: , To find the maximum value element by element, These are the features after the first nonlinear activation. In this embodiment, a second nonlinear screening is performed in a lower-dimensional space through the fifth step. For example, the first-layer ReLU may have already identified the combination of "abnormal temperature gradient + abnormal vibration spectrum"; the second-layer ReLU further determines whether this combination has reached a severity level requiring a warning.

[0057] (6) Dropout: , , For traversal The index numbers of all vectors, This is the output feature after the second Dropout. Since this output is about to enter the final judgment, Dropout here ensures that the network does not overly rely on a single sensor channel (such as only considering vibration and ignoring temperature), maintaining the ability to make multimodal balanced decisions.

[0058] (7) Fully connected Linear: , The characteristics of the third linear output, for The corresponding weight matrix, This is a 1-dimensional bias vector. This step summarizes all the refined anomalous features into a scalar anomaly index. At this point, if all the previous modes are normal, this scalar value will tend to be negative; if early signs such as seed crystal microcracks or misalignment of the pulling axis exist, the scalar value will increase significantly.

[0059] (8) Abnormal risk score Sigmoid: , The abnormal risk score is defined as s, where e is a natural constant. This step maps the scalar index to a probabilistic abnormal risk score of 0-1. When s is less than 0.7, normal growth continues and monitoring continues; when 0.7 (first preset threshold) ≤ s ≤ 0.9, an early warning signal is triggered to alert the operator; when s ≥ 0.9 (second preset threshold), an emergency shutdown signal is triggered to prevent crystal failure or equipment damage. In this embodiment, the probabilities of multiple fault types are calculated using the following formula: ;in, This is a probability distribution of fault types, used to predict the probability of various faults occurring; Activation functions are used to transform any vector into a probability distribution; For the fault classification matrix, learnable parameters that map context vectors to the fault category space, with dimension [missing information]. ; Indicates the type of fault, preferred. (These are seed crystal fracture, interface instability, compositional undercooling, pulling mechanism failure, temperature runaway, and abnormal atmosphere, respectively). The fault classification bias vector is a learnable bias parameter that allows for translation of the classification boundary.

[0060] In this embodiment, the input data format of the lightweight time series prediction network is as follows: (1) Time window: Multimodal feature sequence of the past T time points, T=120 (corresponding to a 60-second time window, sampling frequency 2Hz) (2) Feature dimension: Each time step contains D-dimensional features, D=20 (including multimodal features such as temperature, vibration, weighing, and images). (3) Input tensor shape: ;in, This is the batch size (64 during training, 1 during inference).

[0061] In this embodiment, the output of the lightweight time-series prediction network includes: a sequence of predictor variables for a future time period (preferably the next 15 seconds), an anomaly risk score for the current time, and probabilities of multiple types of faults. The sequence of predictor variables includes the predicted crystal diameter for the next time period, the crystal growth rate, and the degree of deviation of the temperature field from the ideal curve.

[0062] Preferably, the ideal curve is formed by weighted fusion of the following three parts: ;in, Represents an ideal curve; This represents a statistical template curve fitted based on historical normal batch data. This indicates its corresponding weight, preferably 0.5; This represents the theoretical temperature curve calculated based on a thermodynamic model. This indicates its corresponding weight, preferably 0.3; These represent the adaptive curves that are dynamically adjusted based on the preceding data of the current furnace. This indicates its corresponding weight, which is preferably 0.2.

[0063] In this embodiment, Where N is the number of normal furnace runs, For the first The temperature value of the furnace at the current time t.

[0064] In this embodiment, L represents a sliding window with a length of 60 time steps. for The corresponding weighting coefficient is preferably 0.3. It is a temperature trend predicted based on historical data; Temperature deviation is defined as the current actual temperature. Compared to the ideal curve The difference, .

[0065] In this embodiment, the lightweight temporal prediction network is trained using a multi-task joint loss function. include: The prediction loss term is used to constrain the mean squared error between the predicted variable and the true value at multiple future time points. Anomaly detection loss term, used to constrain the binary cross-entropy between anomaly risk score and true anomaly label; The classification loss term is used to constrain the multi-class cross-entropy between the multi-class fault probabilities and the true fault labels. The physical constraint loss term is used to penalize deviations of the predicted variables from the values ​​calculated by the physical model of crystal growth kinetics, which includes the time delay relationship between diameter and temperature.

[0066] ; Predicting losses ; It refers to the number of prediction steps in CrystalNet. =30 (predicting the next 15 seconds) Indicates the first The true value of the step. Indicates the first The predicted value of the step, This is the corresponding loss weight, preferably 1.0.

[0067] Anomaly detection loss We use binary cross-entropy loss to measure the accuracy of anomaly detection. Indicates whether an anomaly truly exists at the current moment. Indicates an anomaly. This indicates that it is normal. This is the anomaly probability predicted by the model, with a value range of [0, 1]. This is the corresponding loss weight, preferably 0.5.

[0068] Classification loss Using multi-class cross-entropy loss, This represents the true probability distribution, with 1 corresponding to the true fault category and 0 for the rest. To predict the probability distribution, preferably, it is analyzed based on six types of faults (seed crystal fracture, interface instability, compositional undercooling, pulling mechanism failure, temperature runaway, and abnormal atmosphere). ;in, It is the encoding of a certain type of fault's true label. This is the corresponding loss weight, preferably 0.3.

[0069] Physical constraint loss , Indicates in The actual observed value at that moment, It is a physical delay, representing the time delay of a causal relationship (such as the time delay between temperature and diameter). A physical model representing the crystal growth kinetics (diameter-temperature relationship). This is the corresponding loss weight, preferably 0.2.

[0070] This invention ensures the accuracy of time-series prediction tasks by constraining the mean square error between predicted variables and true values ​​at multiple future time points through the prediction loss term; improves the accuracy of anomaly identification by constraining the binary cross-entropy between anomaly risk scores and true anomaly labels through the anomaly detection loss term; ensures the reliability of parallel diagnosis of multiple faults by constraining the multi-class cross-entropy between multi-class fault probabilities and true fault labels through the classification loss term; and penalizes the deviation of predicted variables from the calculated values ​​of the crystal growth kinetic physical model through the physical constraint loss term, especially the time delay relationship between diameter and temperature, so that the network prediction results not only fit historical data but also conform to the physical laws of crystal growth.

[0071] In this embodiment, the optimizer used is AdamW. =0.9, =0.999; Initial learning rate: Learning rate scheduling: cosine annealing, 50 epochs; weight decay: .

[0072] Preferably, the validation set prediction loss is monitored during training, and training is stopped if there is no improvement after 20 consecutive epochs.

[0073] In this embodiment, the lightweight temporal prediction network employs structured pruning (based on batch normalization scaling factor γ to prune redundant channels), quantization-aware training (simulating INT8 quantization), and operator fusion techniques, achieving a model compression rate of over 50% and improving inference speed by 3-4 times.

[0074] Model pruning: Based on batch normalization scaling factor γ, pruning rate 30%, accuracy loss after fine-tuning <0.5%.

[0075] Quantization-aware training: Simulate INT8 quantization, compressing the model size from 20MB to 5MB.

[0076] Inference latency: Less than 200ms per inference run on the Jetson Orin NX platform.

[0077] In this embodiment, the method further includes a knowledge management step: in response to the completion of handling a warning event or an emergency shutdown event, the associated data packet corresponding to the warning event or emergency shutdown event is obtained and stored as an abnormal case to realize root cause analysis assistance and process parameter recommendation; The associated data package includes at least: multimodal sensor data, fused feature sequence, predictive variable sequence, abnormal risk score, fault type probability distribution, and handling result label within a preset time window before and after the occurrence of the warning signal or emergency shutdown signal.

[0078] This invention constructs a structured abnormal case library by acquiring the associated data packets corresponding to early warning events or emergency shutdown events and storing them as abnormal cases. This not only preserves the complete multimodal evidence chain when the abnormality occurs, but also associates the handling results, enabling subsequent optimization of prediction network parameters, improvement of latency prior knowledge and fault diagnosis rules through case retrieval, similarity matching or incremental learning, thereby achieving continuous knowledge accumulation and iterative improvement of model performance.

[0079] In this embodiment, the method further includes an edge-cloud collaboration step performed by the edge node: Multimodal sensor data and fused feature sequences are stored, and stored in directories by date and retained for a preset duration (preferably 7 days). At preset time intervals (preferably 5 minutes), the statistical summary data of the multimodal sensor data and the fused feature sequence are sent to the cloud. In response to an abnormal risk score exceeding a first or second preset threshold, multimodal sensor data, fused feature sequences, and predictive variable sequences within a preset time window before and after the occurrence of an early warning signal or emergency shutdown signal are extracted, compressed, and sent to the cloud. It receives updated model parameters from the cloud and updates the lightweight time series prediction network using a hot update method.

[0080] This invention meets the needs of localized data management and rapid backtracking by storing multimodal sensor data and fused feature sequences locally at edge nodes and retaining them in date-based directories for a preset duration. It also ensures efficient network bandwidth utilization by sending statistical summary data to the cloud at preset time intervals, enabling the cloud to monitor the overall operational status. Furthermore, in response to anomaly risk scores exceeding thresholds, it compresses and sends detailed data within a preset time window before and after warning signals or emergency shutdown signals to the cloud, achieving accurate uploading of key data for abnormal events and avoiding bandwidth pressure from continuous full-data uploads. Finally, by receiving updated model parameters from the cloud and updating the lightweight time-series prediction network via hot updates, edge nodes can continuously acquire optimized model capabilities without interrupting the growth process.

[0081] The technical effects of the present invention are verified through the following experimental comparison.

[0082] I. Experimental Environment Table 1. Examples of Crystal Growth Equipment Parameters

[0083] Table 2. Examples of Raw Material Parameters

[0084] Table 3. Examples of Sensor Configurations

[0085] II. Predicted baseline comparison Table 4 Examples of Prediction Comparison Methods

[0086] III. Details of Experimental Results Table 5. Examples of comparison of crystal diameter prediction accuracy

[0087] As shown in Table 5, CrystalNet (lightweight time series prediction network) has a MAE of 0.32mm, which is 42% lower than LSTM (0.55mm) and 33% lower than Transformer.

[0088] Table 6. Examples of Comparison of Prediction Accuracy for Crystal Growth Rate

[0089] As shown in Table 6, the MAE or MAPE index of CrystalNet is reduced by 38% compared to LSTM.

[0090] Table 7 Comparison of Early Warning Lead Time

[0091] As shown in Table 7, CrystalNet achieves an average early warning lead time of 12-15 minutes, which is about twice that of LSTM (5-8 minutes), and represents a fundamental shift from passive response to proactive early warning compared to threshold alarms.

[0092] Table 8 Examples of inference delay comparisons

[0093] As shown in Table 8, after INT8 quantization, the CrystalNet inference latency is reduced to 62ms, which meets the real-time early warning requirements (<200ms), and the model size is compressed to 5.2MB, making it suitable for edge deployment.

[0094] like Figure 2 As shown, the present invention also provides an early warning system for the laser crystal fabrication process, comprising: The multimodal sensor data acquisition module 10 is used to acquire multimodal sensor data, including infrared thermal images, lifting shaft vibration signals, crucible weight signals, atmosphere composition data, and process parameters; the process parameters include lifting rate, rotation rate, and heating power. The feature fusion module 20 is used to perform time alignment and preprocessing on multimodal sensor data, extract multidimensional features of each modality data, fuse the modal features within the same time window to generate a fused feature sequence, and input the fused feature sequence into a pre-trained lightweight temporal prediction network 30. A lightweight temporal prediction network 30 is used to extract multi-scale temporal features based on fused feature sequences, and to generate an attention mask matrix based on preset time delay prior knowledge. The multi-scale temporal features are weighted and fused, and the weighted and fused temporal features are obtained by global average pooling to obtain a context vector. The context vector is input to the prediction decoder, the anomaly risk scoring binary classification network, and the fault type multi-class classification network, respectively, to generate the prediction variable sequence for future time periods, the anomaly risk score for the current time, and the probability of multiple fault types. The multiple fault types include two or more of the following: seed crystal fracture, interface instability, component undercooling, pulling mechanism failure, temperature runaway, and atmosphere anomaly. The early warning module 40 is used to generate an early warning signal when the abnormal risk score is greater than a first preset threshold; and to generate an emergency shutdown signal when the abnormal risk score is greater than a second preset threshold; the second preset threshold is greater than the first preset threshold.

[0095] It should be noted that the various embodiments in this specification are described in a progressive manner, with each embodiment focusing on the differences from other embodiments. Similar or identical parts between embodiments can be referred to interchangeably. For system embodiments, since they are basically similar to method embodiments, the description is relatively simple; relevant parts can be referred to the descriptions in the method embodiments.

[0096] Furthermore, in this document, the terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Without further limitation, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes said element.

[0097] The foregoing description illustrates and describes preferred embodiments of the present invention. It should be understood that the present invention is not limited to the forms disclosed herein and should not be construed as excluding other embodiments. It can be used in various other combinations, modifications, and environments, and can be altered within the scope of the inventive concept by means of the foregoing teachings or techniques or knowledge in related fields. Any modifications and variations made by those skilled in the art that do not depart from the spirit and scope of the present invention should be within the protection scope of the appended claims.

Claims

1. A method for early warning during laser crystal fabrication, characterized in that, Includes the following steps: Acquire multimodal sensor data, including infrared thermal images, lifting shaft vibration signals, crucible weight signals, atmosphere composition data, and process parameters; the process parameters include lifting rate, rotation rate, and heating power; The multimodal sensor data is time-aligned and preprocessed, multidimensional features of each modal data are extracted, and the modal features within the same time window are fused to generate a fused feature sequence. The fused feature sequence is input into a pre-trained lightweight temporal prediction network; Multi-scale temporal features are extracted using the lightweight temporal prediction network, and an attention mask matrix is ​​generated based on preset time delay prior knowledge. The multi-scale temporal features are then weighted and fused. The weighted and fused temporal features are then subjected to global average pooling to obtain a context vector. This process includes: acquiring preset time delay prior knowledge, which at least includes a first time delay interval where temperature changes affect crystal diameter, a second time delay interval where pulling speed changes affect crystal diameter, a third time delay interval where rotational speed changes affect crystal diameter, and a fourth time delay interval where weight change rate affects crystal diameter; for the time difference between the current moment and each historical moment, determining whether the time difference falls within the time delay interval corresponding to the feature type of that historical moment, and generating a binary attention mask matrix based on the determination result; and then, based on the multi-scale temporal features… Calculate the standard attention weight matrix, multiply the standard attention weight matrix element-wise with the binary attention mask matrix, and normalize the multiplication result to obtain constrained attention weights; use the constrained attention weights to perform weighted summation on the multi-scale temporal features to generate weighted fused temporal features; perform global average pooling along the time dimension on the weighted fused temporal features to obtain a context vector; input the context vector into the prediction decoder, and generate a sequence of predictive variables for future time periods through the following steps: enhance the context vector through a fully connected layer and initialize it as the hidden state of the gated recurrent unit in the decoder; use the feature of the last moment of the fused feature sequence as the initial input of the gated recurrent unit, and generate a sequence of predictive variables for a preset number of future steps through autoregression at each moment; The context vector is input into the prediction decoder, the anomaly risk score binary classification network, and the fault type multi-classification network, respectively, to generate the prediction variable sequence for future time periods, the anomaly risk score for the current time, and the probability of multiple fault types; the multiple fault types include two or more of the following: seed crystal fracture, interface instability, component undercooling, pulling mechanism failure, temperature runaway, and atmosphere anomaly. When the abnormal risk score is greater than a first preset threshold, an early warning signal is generated; when the abnormal risk score is greater than a second preset threshold, an emergency shutdown signal is generated; the second preset threshold is greater than the first preset threshold. Extract the multidimensional features of each modality of data, including: The infrared thermal image undergoes non-uniformity correction and bad pixel compensation. The average surface temperature of the melt, surface temperature distribution variance, temperature gradient, isotherm contour radius of curvature, and symmetry index are calculated. This includes: performing non-uniformity correction on the infrared thermal image; detecting pixels that deviate from the neighborhood mean by more than a preset multiple of the standard deviation and replacing them with the neighborhood median; and performing noise reduction processing using Gaussian filtering. Based on an image segmentation algorithm, the region of interest (ROI) in the infrared thermal image is extracted, including the melt surface, crystal sidewalls, and solid-liquid interface. The arithmetic mean of the temperature values ​​of all pixels within the ROI is calculated as the average surface temperature of the melt, and the temperature gradient is calculated. Based on the average surface temperature of the melt... The surface temperature distribution variance is calculated to characterize radial temperature uniformity. All pixels in the infrared thermal image with a pixel temperature equal to a preset temperature threshold are extracted to form isothermal contours. For each pixel grid of a preset size, the state corresponding to each grid vertex is determined based on the preset temperature threshold, and the crossing edge of the isotherm is determined based on the combination of vertex states. The isothermal point positions are calculated. All isothermal points are connected sequentially to form a closed contour curve, and the radius of curvature is calculated, where the radius of curvature is equal to the reciprocal of the curvature of the closed contour curve. Based on the left and right boundary coordinates of the isothermal contour and the reference width, a symmetry index is calculated, which characterizes the degree of symmetry of the isothermal contour relative to the central axis. The vibration signal of the lifting shaft is subjected to bandpass filtering, Fourier transform, and envelope demodulation to obtain the effective amplitude, kurtosis, and characteristic frequency amplitude of the vibration. This includes: fitting a polynomial trend term using the least squares method and subtracting the polynomial trend from the original vibration signal based on the lifting shaft vibration signal; filtering out DC components and high-frequency electronic noise using a Butterworth bandpass filter, designed according to the mechanical characteristic frequency of the crystal growth lifting mechanism; applying a Hanning window to the filtered lifting shaft vibration signal, followed by a fast Fourier transform to calculate the spectrum and power spectral density of the windowed signal, and extracting the energy distribution of a preset frequency band; extracting the envelope spectrum of the lifting shaft vibration signal using a Hilbert transform to identify the modulation frequency components; calculating the square root of the mean square of the vibration acceleration at the sampling points and the kurtosis based on the filtered time-domain signal; and using the square root of the power spectral amplitude as the characteristic frequency amplitude. The crucible weight signal is processed by Kalman filtering, and the crystal diameter and diameter fluctuation rate are calculated based on the filtered crucible weight change rate. This includes: using Kalman filtering to estimate the state of the crucible weight signal and establishing a state-space model; the state variables of the state-space model include the crucible weight and the crucible weight change rate; dynamically adjusting the process noise covariance matrix based on the pulling speed and rotation speed; estimating the optimal weight value and the corresponding crucible weight change rate in real time; using the synchronously acquired pulling shaft vibration signal as a reference, removing the periodic interference component introduced by the pulling shaft vibration from the crucible weight signal; using a low-pass filter to extract the long-term trend of the crucible weight signal after removing the periodic interference component; calculating the crystal diameter in real time based on the physical relationship between the crucible weight change rate and the crystal growth rate; calculating the standard deviation of the crystal diameter sequence using a sliding window to identify abnormal diameter fluctuations and calculate the diameter fluctuation rate. Atmospheric parameters are extracted based on atmospheric composition data, including at least oxygen content and water content.

2. The early warning method for the laser crystal fabrication process according to claim 1, characterized in that, The lightweight temporal prediction network is trained using a multi-task joint loss function, which includes: The prediction loss term is used to constrain the mean squared error between the predicted variable and the true value at multiple future time points. Anomaly detection loss term, used to constrain the binary cross-entropy between anomaly risk score and true anomaly label; The classification loss term is used to constrain the multi-class cross-entropy between the multi-class fault probabilities and the true fault labels; The physical constraint loss term is used to penalize deviations of the predicted variables from the values ​​calculated by the physical model of crystal growth kinetics, which includes the time delay relationship between diameter and temperature.

3. The early warning method for the laser crystal fabrication process according to claim 1, characterized in that, It also includes knowledge management steps: in response to the completion of the handling of a warning event or emergency shutdown event, obtain the associated data packet corresponding to the warning event or emergency shutdown event, and store it as an exception case; The associated data packet includes at least: multimodal sensor data, fused feature sequence, predictive variable sequence, abnormal risk score, fault type probability distribution, and handling result label within a preset time window before and after the occurrence of the warning signal or emergency shutdown signal.

4. The early warning method for the laser crystal fabrication process according to claim 1, characterized in that, It also includes edge-cloud collaboration steps performed by edge nodes: The multimodal sensor data and fused feature sequences are stored, and stored in directories by date with a preset retention period. At preset time intervals, the multimodal sensor data and statistical summary data of the fused feature sequences are sent to the cloud. In response to the abnormal risk score being greater than the first preset threshold or the second preset threshold, multimodal sensor data, fused feature sequences, and predictive variable sequences within a preset time window before and after the occurrence of the warning signal or emergency shutdown signal are extracted, compressed, and sent to the cloud. The model parameters are received from the cloud and the lightweight time series prediction network is updated using a hot update method.

5. An early warning system for laser crystal fabrication process, characterized in that, include: The multimodal sensor data acquisition module is used to acquire multimodal sensor data, including infrared thermal images, lifting shaft vibration signals, crucible weight signals, atmosphere composition data, and process parameters; the process parameters include lifting rate, rotation rate, and heating power. The feature fusion module is used to perform time alignment and preprocessing on the multimodal sensor data, extract multidimensional features of each modal data, and fuse the modal features within the same time window to generate a fused feature sequence. Furthermore, the fused feature sequence is input into a pre-trained lightweight temporal prediction network; Extract the multidimensional features of each modality of data, including: The infrared thermal image undergoes non-uniformity correction and bad pixel compensation. The average surface temperature of the melt, surface temperature distribution variance, temperature gradient, isotherm contour radius of curvature, and symmetry index are calculated. This includes: performing non-uniformity correction on the infrared thermal image; detecting pixels that deviate from the neighborhood mean by more than a preset multiple of the standard deviation and replacing them with the neighborhood median; and performing noise reduction processing using Gaussian filtering. Based on an image segmentation algorithm, the region of interest (ROI) in the infrared thermal image is extracted, including the melt surface, crystal sidewalls, and solid-liquid interface. The arithmetic mean of the temperature values ​​of all pixels within the ROI is calculated as the average surface temperature of the melt, and the temperature gradient is calculated. Based on the average surface temperature of the melt... The surface temperature distribution variance is calculated to characterize radial temperature uniformity. All pixels in the infrared thermal image with a pixel temperature equal to a preset temperature threshold are extracted to form isothermal contours. For each pixel grid of a preset size, the state corresponding to each grid vertex is determined based on the preset temperature threshold, and the crossing edge of the isotherm is determined based on the combination of vertex states. The isothermal point positions are calculated. All isothermal points are connected sequentially to form a closed contour curve, and the radius of curvature is calculated, where the radius of curvature is equal to the reciprocal of the curvature of the closed contour curve. Based on the left and right boundary coordinates of the isothermal contour and the reference width, a symmetry index is calculated, which characterizes the degree of symmetry of the isothermal contour relative to the central axis. The vibration signal of the lifting shaft is subjected to bandpass filtering, Fourier transform, and envelope demodulation to obtain the effective amplitude, kurtosis, and characteristic frequency amplitude of the vibration. This includes: fitting a polynomial trend term using the least squares method and subtracting the polynomial trend from the original vibration signal based on the lifting shaft vibration signal; filtering out DC components and high-frequency electronic noise using a Butterworth bandpass filter, designed according to the mechanical characteristic frequency of the crystal growth lifting mechanism; applying a Hanning window to the filtered lifting shaft vibration signal, followed by a fast Fourier transform to calculate the spectrum and power spectral density of the windowed signal, and extracting the energy distribution of a preset frequency band; extracting the envelope spectrum of the lifting shaft vibration signal using a Hilbert transform to identify the modulation frequency components; calculating the square root of the mean square of the vibration acceleration at the sampling points and the kurtosis based on the filtered time-domain signal; and using the square root of the power spectral amplitude as the characteristic frequency amplitude. The crucible weight signal is processed by Kalman filtering, and the crystal diameter and diameter fluctuation rate are calculated based on the filtered crucible weight change rate. This includes: using Kalman filtering to estimate the state of the crucible weight signal and establishing a state-space model; the state variables of the state-space model include the crucible weight and the crucible weight change rate; dynamically adjusting the process noise covariance matrix based on the pulling speed and rotation speed; estimating the optimal weight value and the corresponding crucible weight change rate in real time; using the synchronously acquired pulling shaft vibration signal as a reference, removing the periodic interference component introduced by the pulling shaft vibration from the crucible weight signal; using a low-pass filter to extract the long-term trend of the crucible weight signal after removing the periodic interference component; calculating the crystal diameter in real time based on the physical relationship between the crucible weight change rate and the crystal growth rate; calculating the standard deviation of the crystal diameter sequence using a sliding window to identify abnormal diameter fluctuations and calculate the diameter fluctuation rate. Atmospheric parameters are extracted based on atmospheric composition data, including at least oxygen content and water content; A lightweight temporal prediction network is used to extract multi-scale temporal features based on the fused feature sequence, generate an attention mask matrix based on preset time delay prior knowledge, perform weighted fusion of the multi-scale temporal features, and obtain a context vector by global average pooling of the weighted fused temporal features. The network includes: acquiring preset time delay prior knowledge, which at least includes a first time delay interval for temperature changes affecting crystal diameter, a second time delay interval for pulling speed changes affecting crystal diameter, a third time delay interval for rotation speed changes affecting crystal diameter, and a fourth time delay interval for weight change rate affecting crystal diameter; determining whether the time difference between the current moment and each historical moment falls within the time delay interval corresponding to the feature type of the historical moment, and generating a binary attention mask matrix based on the determination result; calculating a standard attention weight matrix based on the multi-scale temporal features, multiplying the standard attention weight matrix element-wise with the binary attention mask matrix, and normalizing the multiplication result to obtain a constrained attention mask matrix. Attention weights are used to weight and sum the multi-scale temporal features to generate weighted fused temporal features. The weighted fused temporal features are then globally averaged along the time dimension to obtain a context vector. After inputting the context vector into the prediction decoder, a sequence of predictive variables for future time periods is generated through the following steps: the context vector is enhanced by a fully connected layer and initialized as the hidden state of a gated recurrent unit in the decoder; the last moment feature of the fused feature sequence is used as the initial input of the gated recurrent unit, which then generates a sequence of predictive variables for a predetermined number of future steps through autoregression at each moment; the context vector is input into the prediction decoder, the anomaly risk scoring binary classification network, and the fault type multi-class classification network to generate a sequence of predictive variables for future time periods, an anomaly risk score for the current moment, and multiple fault probabilities, respectively; the multiple faults include two or more of the following: seed crystal fracture, interface instability, component overcooling, lifting mechanism failure, temperature runaway, and abnormal atmosphere. The early warning module is used to generate an early warning signal when the abnormal risk score is greater than a first preset threshold; and to generate an emergency shutdown signal when the abnormal risk score is greater than a second preset threshold; wherein the second preset threshold is greater than the first preset threshold.

Images