Method for predicting flow rate of silicon melt in czochralski silicon single crystal growth process

Abstract

The application discloses a method for predicting the flow velocity of a silicon melt in a Czochralski silicon single crystal growth process, which comprises the following steps: firstly, a two-dimensional axisymmetric rotation model of the silicon melt in the Czochralski silicon single crystal growth process is established, and the boundary conditions of the fluid model are set; then, a space factor physical information neural network is constructed; finally, the two-dimensional axisymmetric rotation model of the silicon melt in the Czochralski silicon single crystal growth process is solved by using the space factor physical information neural network, and finally the prediction result of the flow velocity of the silicon melt is obtained. The application solves the problem that the flow velocity of the silicon melt in the Czochralski silicon single crystal growth process cannot be predicted in real time in the prior art.

Description

Technical Field

[0001] This invention belongs to the field of semiconductor silicon single crystal material preparation technology, specifically relating to a method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process. Background Technology

[0002] Integrated circuit chips play a crucial role in the development of society and science and technology. Silicon, as the basic material for chip manufacturing, is one of the most important semiconductor materials. The Czochralski method, also known as the CZ method, is the main method for preparing silicon single crystals and is widely used in the growth and preparation of large-diameter, high-quality semiconductor silicon crystals. Flow field, temperature field, and magnetic field are the main physical fields in silicon single crystal growth research. Among them, the flow field serves as the starting point for the study of other physical fields and contains various complex convection currents, such as natural convection, forced convection, and Marangoni convection. The interaction between these complex convection currents affects the temperature distribution, concentration distribution, and impurity distribution of the silicon melt during silicon single crystal growth. Predicting the flow rate of the silicon melt is beneficial for analyzing the internal convection and temperature distribution of the melt, enabling the optimizer to adopt appropriate control strategies, which is of great significance for optimizing the crystal growth process and improving crystal quality. However, in actual engineering, due to the complex furnace environment and excessively high silicon melt temperature, it is impossible to directly measure the flow rate of the silicon melt. Computational fluid dynamics (CFD) methods are commonly used to predict the velocity of molten silicon, requiring discretization of the fluid model before solution processing. The resulting discrete solutions limit applications in real-time prediction and model optimization control. This invention proposes a spatial factorial physical information neural network-based method for predicting the velocity of molten silicon, satisfying the strong form of partial differential equations in the fluid model. After network training, only the spatial coordinates of the target point need to be input to quickly predict the velocity of the molten silicon at that location. This method can be applied to silicon single-crystal growth model control, crystal quality optimization, and other applications. Summary of the Invention

[0003] The purpose of this invention is to provide a method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process, which solves the problem in the prior art that the silicon melt flow rate cannot be predicted in real time during the silicon single crystal growth process.

[0004] The technical solution adopted in this invention is a method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process, which is implemented according to the following steps:

[0005] Step 1: Establish a two-dimensional axisymmetric rotational model of silicon melt during the Czochralski silicon single crystal growth process, and set the boundary conditions of the fluid model;

[0006] Step 2: Construct a spatial factor physical information neural network;

[0007] Step 3: The spatial factor physical information neural network used in Step 2 is used to solve the two-dimensional axisymmetric rotation model of silicon melt in the Czochralski silicon single crystal growth process in Step 1, and finally the prediction result of silicon melt flow rate is obtained.

[0008] The invention is further characterized in that,

[0009] Step 1 is as follows:

[0010] Step 1.1: Construct a fluid model of silicon melt that includes the continuity equation and the momentum equation;

[0011] Step 1.2: Establish a silicon melt fluid model in cylindrical coordinates;

[0012] Step 1.3: Construct the boundary conditions for a two-dimensional axisymmetric rotational model of silicon melt.

[0013] Step 1.1 is as follows:

[0014]

[0015]

[0016] In the formula Represents the divergence of a fluid. It is the Hamiltonian operator. It is the velocity vector in the flow field. Let ρ represent the convection term, P be the fluid density, P be the pressure in the flow field, and μ be the viscosity coefficient of the silicon melt. It is the Laplace operator, and g is the acceleration due to gravity;

[0017] Expanding equations (1) and (2) yields:

[0018]

[0019]

[0020] In the formula, u, v, and w represent the velocity components of the x-axis, y-axis, and z-axis, respectively.

[0021] Step 1.2 is as follows:

[0022] Treating the silicon melt as a Newtonian, incompressible, axisymmetric rotating fluid, with the axis of symmetry of the crystal and crucible as the axis of symmetry, based on the above characteristics, equations (3) and (4) are transformed into cylindrical coordinates to establish a silicon melt fluid model in cylindrical coordinates:

[0023]

[0024]

[0025] In the formula vr v z v θ These represent the radial, axial, and azimuth velocity components, respectively, where r, z, and θ are the coordinates of the radial, axial, and azimuth directions.

[0026] Since the model is symmetric about the rotation axis, equations (5) and (6) satisfy the following for any variable η: The two-dimensional axisymmetric rotational model of silicon melt is expressed as follows:

[0027]

[0028]

[0029] Equations (7) and (8) are the final two-dimensional axisymmetric rotational model of the silicon melt in cylindrical coordinates.

[0030] Step 1.3 is as follows:

[0031] This two-dimensional axisymmetric rotational model of the silicon melt contains five boundaries, and the boundary conditions for each boundary are as follows:

[0032] Axisymmetric boundary: The boundary position is r = 0, 0 ≤ z ≤ H, and the corresponding boundary condition is v. r =0, v θ =0、

[0033] Bottom of crucible: Boundary position is 0≤r≤R c z = 0, the corresponding boundary condition is v r =0, v z =0, v θ =ω c ·r;

[0034] Crucible sidewall: boundary position is r = R c 0≤z≤H, the corresponding boundary condition is v r =0, v z =0, v θ =ω c ·r;

[0035] Solid-liquid interface: the boundary position is 0≤r≤R x z = H, and the corresponding boundary condition is v r =0, v z =0, v θ =ω x ·r;

[0036] Free surface: Boundary location is R x ≤r≤R c The boundary conditions for z = H are:

[0037] Where R is x Crystal radius, R c Let H be the radius of the crucible and H be the height of the crucible.

[0038] Step 2 is as follows:

[0039] Step 2.1: Construct a physical information neural network:

[0040] First, the general expression of the partial differential equation is given:

[0041]

[0042] In the formula, F is a nonlinear operator, and u represents a solution that satisfies the partial differential equation. For spatial domain,

[0043] The boundary conditions are expressed as follows:

[0044]

[0045] In the formula The boundary condition g(X) represents the boundary condition that the solution to the partial differential equation must satisfy.

[0046] The input to the physical information neural network is X = x1, x2, ... x d The output of the physical information neural network is In the formula, θ represents the network parameters, including the weight matrix W, the bias vector b, and the physical information neural network output. The dimension of the physical information neural network remains the same as that of the true solution U(X) of the partial differential equation, and it serves as an alternative solution to U(X); the expression for the physical information neural network is as follows:

[0047] H 1 =σ(W 1 X+b 1 (11)

[0048] H k =σ(W k H k-1 +b k (12)

[0049]

[0050] In the formula H k Let σ represent the output of the k-th layer of the network, σ be the activation function, X be the spatial sample points of the input, and W be the output of the k-th layer of the network. k b k This represents the weight matrix and bias vector of the k-th layer. This represents the final output of the network, where L is the total number of layers in the physical information neural network.

[0051] The loss function L of a physical information neural network consists of two parts: boundary point loss L b Configuration point loss L f Boundary points are sampling points at the boundary of the simulation domain and must satisfy boundary conditions; collocation points are sampling points inside the simulation domain and must satisfy the partial differential equations corresponding to the model, as follows:

[0052] L = L b +L f (14)

[0053]

[0054]

[0055] In the formula, N b N represents the number of boundary points. f Indicates the number of configuration points. This indicates that the network is connected to the i-th boundary point. The output value, Represents the i-th boundary point The corresponding actual value, This indicates that the network is configured at the i-th configuration point. The output value, This represents the error of the m-th equation in a system of partial differential equations;

[0056] Step 2.2, Spatial Factor Physical Information Neural Network:

[0057] The original spatial input X is reintroduced into each hidden layer of the network to construct a spatial factor physical information neural network. The hidden layer equation (12) is updated as follows:

[0058] H k =σ(W k H k-1 +b k )·X+b e (17)

[0059] In the formula, b e For the additional bias vector, · denotes vector multiplication;

[0060] Two encoders, M and N, are used instead of X and b. e The hidden layer (17) is updated as follows:

[0061] H k =σ(W k H k-1 +b k )·M+N (18)

[0062] In the formula, M = σ(W) M X+b M ), N = σ(W N X+b N ), W M b M These are the weight matrix and bias vector of encoder M, respectively, and W. M b M These are the weight matrix and bias vector of encoder N, respectively;

[0063] The expression for the spatial factor physical information neural network is:

[0064] H 1 =σ(W 1 X+b 1 (19)

[0065] H k =σ(W k H k-1 +b k )·M+N (20)

[0066]

[0067] The loss function of the spatial factor physical information neural network is composed of the boundary point loss L. b Configuration point loss L f Composed of two parts, during the training process, in order to balance the boundary point loss and the placement point loss, when the boundary point loss reaches the set target value, the optimization of the placement point loss is strengthened. Therefore, the loss function (14) to (16) of the spatial parameter physical information neural network is improved as follows:

[0068]

[0069]

[0070]

[0071] In the formula, tol is the threshold of the boundary loss, and ε is the switching function, satisfying:

[0072]

[0073] Step 3 is as follows:

[0074] Step 3.1: Set the specific parameters of the spatial factor physical information neural network:

[0075] The spatial factor physical information neural network has a 2-neuron input layer for inputting the spatial coordinates radius r and height z, and a 4-neuron output layer for outputting the radial velocity v. raxial velocity v z azimuth angular velocity v θ The training process involves a pressure P, 5 hidden layers with 128 neurons per layer, an initial learning rate of 0.002, and a learning rate that decreases to 0.9 times the initial value every 1000 iterations. The entire training process consists of 30,000 iterations.

[0076] Step 3.2: Construct the training set:

[0077] The training set includes the boundary point set X. b Set of configuration points X f The boundary points in the boundary point set need to be sampled at the boundaries of the simulation region. There are five boundaries in total: axisymmetric boundary, crucible bottom, crucible sidewall, solid-liquid interface, and free surface, which are used to satisfy the set boundary conditions. The sampling method is equal-interval sampling, and the number of sampling points for each boundary is 400. The configuration points in the configuration point set are sampled inside the simulation region to satisfy the partial differential equation. The sampling method is Latin hypercube sampling, and the number of samples is 6000.

[0078] Step 3.3: Construct the loss function for the two-dimensional axisymmetric rotational model of the silicon melt:

[0079] Based on the established two-dimensional axisymmetric rotation model of silicon melt, a loss function L for the spatial factor physical information neural network is constructed, comprising two parts: boundary point loss L. b Configuration point loss L f Boundary point loss function L b as follows:

[0080]

[0081] Location point loss function L f as follows:

[0082]

[0083]

[0084]

[0085]

[0086]

[0087] The total loss L is:

[0088]

[0089] Step 3.4: Substitute the experimental parameters:

[0090] Physical parameters of silicon melt: ρ = 2530 kg / m 3 μ = 0.0008 kg / m·s;

[0091] Assuming the crystal radius is R x =0.15m, crucible radius is R c =0.3m, crucible height is H=0.18m, crucible rotation speed w c =0.003 r / min, crystal rotation speed w x =0.003 r / min;

[0092] Step 3.5: Train the spatial factor physical information neural network and make predictions:

[0093] The experiment was conducted on the Python framework Tensorflow-GPU 1.15.2. The computer configuration was as follows: CPU: AMD Ryzen 55600G, frequency: 3.90GHz, RAM: 16GB; Graphics card: RTX 2060, VRAM: 6GB.

[0094] The ADAM gradient descent method was used to optimize the network parameters with an exponential decay rate of 0.99. Training was stopped when the number of iterations reached 30,000 or the loss value fell below 1e-5. After the network was trained, the radius r and height z of the spatial coordinates were used as inputs to predict the flow velocity in the simulated region.

[0095] The beneficial effect of this invention is that it enables the solution of the silicon melt fluid model during the Czochralski silicon single crystal growth process without discretizing the model. The spatial factor physical information neural network can express any strongly nonlinear relationship and has a fast convergence speed. The trained spatial factor physical information neural network has high accuracy and strong real-time performance; it only requires input spatial coordinates and does not require interpolation to quickly predict the flow velocity of the silicon melt at the corresponding location, and can be used in model control, optimization, and other application scenarios. Attached Figure Description

[0096] Figure 1 This is a flowchart of a method for predicting the flow rate of silicon melt;

[0097] Figure 2 This is a schematic diagram of a single crystal furnace based on the Czochralski method;

[0098] Figure 3 It is the simulated region for the flow of molten silicon;

[0099] Figure 4 This is a schematic diagram of the structure of a spatial factor physical information neural network;

[0100] Figure 5 It is a graph showing the velocity distribution and error predicted by a physical information neural network;

[0101] Figure 6 It is a velocity distribution and error map predicted by a spatial factor physical information neural network;

[0102] Figure 7 It is a velocity streamline diagram predicted by a spatial factor physical information neural network;

[0103] Figure 8 It is the convergence curve of the spatial factor physical information neural network and the physical information neural network;

[0104] Figure 9 These are the convergence curves of the spatial factor physical information neural network and the physical information neural network during the 25,000-30,000 training iterations.

[0105] In the diagram, 1. Crystal, 2. Heater, 3. Heat shield, 4. Silicon melt, 5. Crucible, 6. Crucible shaft. Detailed Implementation

[0106] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments.

[0107] The present invention provides a method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process. The overall process is as follows: Figure 1 As shown, please follow these steps:

[0108] Step 1: Establish a two-dimensional axisymmetric rotational model of silicon melt during the Czochralski silicon single crystal growth process, and set the boundary conditions of the fluid model;

[0109] Step 1 is as follows:

[0110] The overall structure of a single crystal furnace based on the Czochralski method is as follows: Figure 2 As shown, heater 2 and heat shield 3 are located around crucible 5, silicon melt 4 is located inside crucible 5, and crystal 1 is located above silicon melt 4. The simulated flow area of ​​silicon melt 4 is as follows. Figure 3 As shown. The upper part of the column is crystal 1 with radius R. x Rotate clockwise at a speed of Ω x The lower half is crucible 5, with a radius of R. c Rotating counterclockwise at a speed of Ω c The gray cross-section represents the simulated region of the silicon melt.

[0111] Step 1.1: Construct a fluid model of silicon melt that includes continuity and momentum equations;

[0112] Step 1.1 is as follows:

[0113]

[0114]

[0115] In the formula Represents the divergence of a fluid. It is the Hamiltonian operator. It is the velocity vector in the flow field. Let ρ represent the convection term, P be the fluid density, P be the pressure in the flow field, and μ be the viscosity coefficient of the silicon melt. It is the Laplace operator, and g is the acceleration due to gravity;

[0116] Expanding equations (1) and (2) yields:

[0117]

[0118]

[0119] In the formula, u, v, and w represent the velocity components of the x-axis, y-axis, and z-axis, respectively.

[0120] Step 1.2: Establish a silicon melt fluid model in cylindrical coordinates;

[0121] Step 1.2 is as follows:

[0122] Treating the silicon melt as a Newtonian, incompressible, axisymmetric rotating fluid, with the axis of symmetry of the crystal and crucible as the axis of symmetry, based on the above characteristics, equations (3) and (4) are transformed into cylindrical coordinates to establish a silicon melt fluid model in cylindrical coordinates:

[0123]

[0124]

[0125] In the formula v r v z v θ These represent the radial, axial, and azimuth velocity components, respectively, where r, z, and θ are the coordinates of the radial, axial, and azimuth directions.

[0126] Since the model is symmetric about the rotation axis, equations (5) and (6) satisfy the following for any variable η: The two-dimensional axisymmetric rotational model of silicon melt is expressed as follows:

[0127]

[0128]

[0129] Equations (7) and (8) are the final two-dimensional axisymmetric rotational model of the silicon melt in cylindrical coordinates.

[0130] Step 1.3: Construct the boundary conditions for the two-dimensional axisymmetric rotation model of the silicon melt in sequence.

[0131] Step 1.3 is as follows:

[0132] This two-dimensional axisymmetric rotational model of the silicon melt contains five boundaries, and the boundary conditions for each boundary are as follows:

[0133] Axisymmetric boundary: The boundary position is r = 0, 0 ≤ z ≤ H, and the corresponding boundary condition is v. r =0, v θ =0、

[0134] Bottom of crucible: Boundary position is 0≤r≤R c z = 0, the corresponding boundary condition is v r =0, v z =0, v θ =ω c ·r;

[0135] Crucible sidewall: boundary position is r = R c 0≤z≤H, the corresponding boundary condition is v r =0, v z =0, v θ =ω c ·r;

[0136] Solid-liquid interface: Boundary position is 0≤r≤R x z = H, and the corresponding boundary condition is v r =0, v z =0, v θ =ω x ·r;

[0137] Free surface: Boundary location is R x ≤r≤R c The boundary conditions for z = H are:

[0138] Where R is x Crystal radius, R c Let H be the radius of the crucible and H be the height of the crucible.

[0139] The specific boundary conditions for each boundary are shown in Table 1:

[0140] Table 1 shows the boundary conditions for each boundary.

[0141]

[0142] In the formula R x R is the crystal radius. c Let H be the radius of the crucible and H be the height of the crucible.

[0143] Step 2: Construct a spatial factor physical information neural network;

[0144] Step 2 is as follows:

[0145] Step 2.1: Construct a physical information neural network:

[0146] First, the general expression of the partial differential equation is given:

[0147]

[0148] In the formula, F is a nonlinear operator, and u represents a solution that satisfies the partial differential equation. For spatial domain,

[0149] The boundary conditions are expressed as follows:

[0150]

[0151] In the formula The boundary condition g(X) represents the boundary condition that the solution to the partial differential equation must satisfy.

[0152] The input to the physical information neural network is X = x1, x2, ... x d The output of the physical information neural network is In the formula, θ represents the network parameters, including the weight matrix W, the bias vector b, and the physical information neural network output. The dimension of the physical information neural network remains the same as that of the true solution U(X) of the partial differential equation, and it serves as an alternative solution to U(X); the expression for the physical information neural network is as follows:

[0153] H 1 =σ(W 1 X+b 1 (11)

[0154] H k =σ(W k H k-1 +b k (12)

[0155]

[0156] In the formula H k Let σ represent the output of the k-th layer of the network, σ be the activation function, X be the spatial sample points of the input, and W be the output of the k-th layer of the network. k b k This represents the weight matrix and bias vector of the k-th layer. This represents the final output of the network, where L is the total number of layers in the physical information neural network.

[0157] The loss function L of a physical information neural network consists of two parts: boundary point loss Lb Configuration point loss L f Boundary points are sampling points at the boundary of the simulation domain and must satisfy boundary conditions; collocation points are sampling points inside the simulation domain and must satisfy the partial differential equations corresponding to the model, as follows:

[0158] L = L b +L f (14)

[0159]

[0160]

[0161] In the formula, N b N represents the number of boundary points. f Indicates the number of configuration points. This indicates that the network is connected to the i-th boundary point. The output value, Represents the i-th boundary point The corresponding actual value, This indicates that the network is configured at the i-th configuration point. The output value, This represents the error of the m-th equation in a system of partial differential equations;

[0162] Step 2.2, Spatial Factor Physical Information Neural Network:

[0163] Physical information neural networks (PINs) encode partial differential equations into their loss functions, requiring the use of automatic differentiation techniques from deep learning frameworks to obtain derivative information of the input spatial variables and incorporate it into the loss function calculation. This means that the loss function of a PIN is not only determined by the deep features of the network but also influenced by the spatial features of the original input—a significant difference between PINs and traditional deep neural networks. This leads to a problem: as the network depth increases, the derivative information in the partial differential equation portion of the loss becomes weaker, diminishing the representation of the original features in the loss function. This situation is similar to the vanishing gradient problem when optimizing network parameters.

[0164] To address this problem, this invention proposes an optimized physical information neural network architecture, which reintroduces the original spatial input X into each hidden layer of the network to construct a spatial factor physical information neural network, thereby strengthening the correlation between each layer of the network and the original spatial input. The hidden layer equation (12) is updated as follows:

[0165] H k =σ(W k H k-1 +b k )·X+b e (17)

[0166] In the formula, b e For the additional bias vector, · denotes vector multiplication;

[0167] To enhance the additional introduction of X and b e To enhance expressive power, this invention employs two encoders, M and N, instead of X and b. e Another function of the encoder is to change the dimension of the original spatial input X to match the hidden layer. The hidden layer equation (17) is updated as follows:

[0168] H k =σ(W k H k-1 +b k )·M+N (18)

[0169] In the formula, M = σ(W) M X+b M ), N = σ(W N X+b N ), W M b M These are the weight matrix and bias vector of encoder M, respectively, and W. M b M These are the weight matrix and bias vector of encoder N, respectively;

[0170] The structure of the spatial factor physical information neural network is as follows: Figure 3 As shown, the expression is:

[0171] H 1 =σ(W 1 X+b 1 (19)

[0172] H k =σ(W k H k-1 +b k )·M+N (20)

[0173]

[0174] The loss function of the spatial factor physical information neural network is composed of the boundary point loss L. b Configuration point loss L f Composed of two parts, during the training process, in order to balance the boundary point loss and the placement point loss, when the boundary point loss reaches the set target value, the optimization of the placement point loss is strengthened. Therefore, the loss function (14) to (16) of the spatial parameter physical information neural network is improved as follows:

[0175]

[0176]

[0177]

[0178] In the formula, tol is the threshold of the boundary loss, and ε is the switching function, satisfying:

[0179]

[0180] In summary, the spatial factor physical information neural network improves the physical information neural network from two aspects: network structure and loss function.

[0181] The spatial factor physical information neural network proposed in this invention has the following advantages: the spatial factor physical information neural network reintroduces spatial information during the forward propagation process of the network, which strengthens the correlation between the network and the original spatial input and is conducive to recovering the derivative information of the initial spatial input during the backward propagation process; at the same time, the encoder introduced into the network structure enhances the expressive power of the network.

[0182] Step 3: The spatial factor physical information neural network used in Step 2 is used to solve the two-dimensional axisymmetric rotation model of silicon melt in the Czochralski silicon single crystal growth process in Step 1, and finally the prediction result of silicon melt flow rate is obtained.

[0183] Step 3 is as follows:

[0184] Step 3.1: Set the specific parameters of the spatial factor physical information neural network:

[0185] The spatial factor physical information neural network has a 2-neuron input layer for inputting the spatial coordinates radius r and height z, and a 4-neuron output layer for outputting the radial velocity v. r axial velocity v z azimuth angular velocity v θ The training process involves a pressure P, 5 hidden layers with 128 neurons per layer, an initial learning rate of 0.002, and a learning rate that decreases to 0.9 times the initial value every 1000 iterations. The entire training process consists of 30,000 iterations.

[0186] Step 3.2: Construct the training set:

[0187] The training set includes the boundary point set X. b Set of configuration points X fThe boundary points in the boundary point set need to be sampled at the boundaries of the simulation region. There are five boundaries in total: axisymmetric boundary, crucible bottom, crucible sidewall, solid-liquid interface, and free surface. Their corresponding boundary locations are shown in Table 1. To satisfy the set boundary conditions, the sampling method is equal-interval sampling, with 400 sampling points for each boundary. The configuration points in the configuration point set are sampled inside the simulation region to satisfy the partial differential equations. The sampling method is Latin hypercube sampling, with 6000 sampling points.

[0188] Step 3.3: Construct the loss function for the two-dimensional axisymmetric rotational model of the silicon melt:

[0189] Based on the established two-dimensional axisymmetric rotation model of silicon melt, a loss function L for the spatial factor physical information neural network is constructed, comprising two parts: boundary point loss L. b Configuration point loss L f Boundary point loss function L b as follows:

[0190]

[0191] Location point loss function L f as follows:

[0192]

[0193]

[0194]

[0195]

[0196]

[0197] The total loss L is:

[0198]

[0199] Step 3.4: Substitute the experimental parameters:

[0200] Physical parameters of silicon melt: ρ = 2530 kg / m 3 μ = 0.0008 kg / m·s;

[0201] Assuming the crystal radius is R x =0.15m, crucible radius is R c =0.3m, crucible height is H=0.18m, crucible rotation speed w c =0.003 r / min, crystal rotation speed w x =0.003 r / min;

[0202] Step 3.5: Train the spatial factor physical information neural network and make predictions:

[0203] The experiment was conducted on the Python framework Tensorflow-GPU 1.15.2. The computer configuration was as follows: CPU: AMD Ryzen 55600G, frequency: 3.90GHz, RAM: 16GB; Graphics card: RTX 2060, VRAM: 6GB.

[0204] The ADAM gradient descent method was used to optimize the network parameters with an exponential decay rate of 0.99. Training was stopped when the number of iterations reached 30,000 or the loss value fell below 1e-5. After the network was trained, the radius r and height z of the spatial coordinates were used as inputs to predict the flow velocity in the simulated region.

[0205] The flow velocity and absolute error predicted by the physical information neural network and the spatial factor physical information neural network are respectively as follows: Figure 5 , Figure 6 As shown in the formula, the horizontal axis represents the radius r, and the vertical axis represents the height z. The spatial factor physical information neural network predicts streamlines as follows: Figure 7 As shown, two forced convections can be observed, one induced by the crucible and the other by the crystal.

[0206] Figure 8 The graphs show the convergence of the physical information neural network and the spatial factor physical information neural network during 30,000 training iterations. The spatial factor physical information neural network converges smoothly and rapidly, exhibiting significantly stronger convergence than the physical information neural network. Figure 9 The convergence plot is shown for iterations 25000-30000. Table 2 shows the loss values ​​and L2 errors of the physical information neural network and the spatial factor physical information neural network. The L2 error expression is as follows:

[0207]

[0208] Table 2

[0209]

[0210] Compared to the physical information neural network, the spatial factor physical information neural network reduces the loss value and L2 error of velocity by two orders of magnitude and the L2 error of pressure by one order of magnitude, showing significantly higher prediction accuracy than the physical information neural network.

Claims

1. A method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process, characterized in that, The specific steps are as follows: Step 1: Establish a two-dimensional axisymmetric rotational model of silicon melt during the Czochralski silicon single crystal growth process, and set the boundary conditions of the fluid model; Step 1 is described in detail as follows: Step 1.1: Construct a fluid model of silicon melt that includes continuity and momentum equations; Step 1.2: Establish a silicon melt fluid model in cylindrical coordinates; Step 1.3: Construct the boundary conditions for a two-dimensional axisymmetric rotational model of the silicon melt; Step 2: Construct a spatial factor physical information neural network; Step 2 includes the following steps: Step 2.1: Construct a physical information neural network: Loss function of physical information neural network It consists of two parts: boundary point loss Configuration point loss Boundary points are sampling points at the boundary of the simulation domain and must satisfy boundary conditions; collocation points are sampling points inside the simulation domain and must satisfy the partial differential equations corresponding to the model, as follows: （14） （15） （16） In the formula, Indicates the number of boundary points. Indicates the number of configuration points. Indicates the network's response to the first... boundary points The output value, Indicates the first boundary points The corresponding actual value, Indicates the network's response to the first... Configuration points The output value, In a system of partial differential equations, the first... The error of each equation; Step 2.2, Spatial Factor Physical Information Neural Network: The loss function of the spatial factor physical information neural network is composed of boundary point loss. Configuration point loss Composed of two parts, during the training process, in order to balance the boundary point loss and the placement point loss, when the boundary point loss reaches the set target value, the optimization of the placement point loss is strengthened. Therefore, the loss function (14)~(16) of the spatial parameter physical information neural network is improved as follows: （22） （23） （24） In the formula The threshold for boundary loss, For a switching function, satisfying: （25）； Step 3: The spatial factor physical information neural network used in Step 2 is used to solve the two-dimensional axisymmetric rotation model of silicon melt in the Czochralski silicon single crystal growth process in Step 1, and finally the prediction result of silicon melt flow rate is obtained.

2. The method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process according to claim 1, characterized in that, Step 1.1 is as follows: （1） （2） In the formula Represents the divergence of a fluid. It is the Hamiltonian operator. It is the velocity vector in the flow field. Represents the convection term. It is the density of the fluid. It is the pressure in the flow field. It is the viscosity coefficient of the silicon melt. It is the Laplace operator. It is gravitational acceleration; Expanding equations (1) and (2) yields: （3） （4） In the formula , , They represent axis, axis, The velocity component of the shaft.

3. The method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process according to claim 2, characterized in that, Step 1.2 is as follows: Treating the silicon melt as a Newtonian, incompressible, axisymmetric rotating fluid, with the axis of symmetry of the crystal and crucible as the axis of symmetry, based on the above characteristics, equations (3) and (4) are transformed into cylindrical coordinates to establish a silicon melt fluid model in cylindrical coordinates: （5） In the formula , , These represent the radial, axial, and azimuth velocity components, respectively. , , These are the radial, axial, and azimuth coordinates. Since the model is symmetric about the rotation axis, equations (5) and (6) are valid for any variable. satisfy The two-dimensional axisymmetric rotational model of silicon melt is expressed as follows: （7） （8） Equations (7) and (8) are the final two-dimensional axisymmetric rotational model of the silicon melt in cylindrical coordinates.

4. The method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process according to claim 3, characterized in that, Step 1.3 is as follows: This two-dimensional axisymmetric rotational model of the silicon melt contains five boundaries, and the boundary conditions for each boundary are as follows: Axisymmetric boundary: The boundary position is , The corresponding boundary conditions are , , ; Bottom of the crucible: The boundary position is , The corresponding boundary conditions are , , ; Crucible sidewall: boundary position is , The corresponding boundary conditions are , , ; Solid-liquid interface: the boundary location is , The corresponding boundary conditions are , , ; Free surface: boundary location is , The corresponding boundary conditions are , , ; Among them Crystal radius, Where is the radius of the crucible. This represents the height of the crucible.

5. The method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process according to claim 4, characterized in that, Step 2 is described in detail below: Step 2.1: Construct a physical information neural network: First, the general expression of the partial differential equation is given: （9） In the formula It is a nonlinear operator. This represents a solution that satisfies the partial differential equation. For spatial domain, ; The boundary conditions are expressed as follows: （10） In the formula The boundary conditions represent the boundary conditions that the solution to the partial differential equation must satisfy. ; The input to the physical information neural network is The output of the physical information neural network is In the formula Network parameters, including the weight matrix. Bias vector Physical information neural network output With the true solution of the partial differential equation The dimensions remain the same, and are used as The alternative solution; the expression for the physical information neural network is as follows: （11） （12） （13） In the formula Indicates the network number The output of the layer, For activation function, For the input spatial sample points, , Indicates the first The weight matrix and bias vector of the layer, This represents the final output of the network. This represents the total number of layers in the physical information neural network. Step 2.2, Spatial Factor Physical Information Neural Network: Input the original space By reintroducing each hidden layer of the network and constructing a spatial factor physical information neural network, the hidden layer equation (12) is updated as follows: （17） In the formula, For additional bias vectors, Represents vector multiplication; Two encoders are used , replace , The hidden layer (17) is updated as follows: （18） In the formula , , , They are encoders The weight matrix and bias vector, , They are encoders The weight matrix and bias vector; The expression for the spatial factor physical information neural network is: （19） （20） （21）。 6. The method for predicting the silicon melt flow rate during the Czochralski silicon single crystal growth process according to claim 5, characterized in that, Step 3 is described in detail below: Step 3.1: Set the specific parameters of the spatial factor physical information neural network: The input layer of the spatial factor physical information neural network consists of two neurons, used to input the radius of the spatial coordinates. and height The output layer has four neurons, each outputting radial velocity. axial velocity azimuth angular velocity ,pressure The hidden layers consist of 5 layers, each with 128 neurons. The initial learning rate is set to 0.

002. Every 1000 iterations, the learning rate is reduced to 0.9 times the previous value. The entire training process involves 30,000 iterations. Step 3.2: Construct the training set: The training set includes the set of boundary points. Set of configuration points The boundary points in the boundary point set need to be sampled at the boundaries of the simulation region. There are five boundaries in total: axisymmetric boundary, crucible bottom, crucible sidewall, solid-liquid interface, and free surface, which are used to satisfy the set boundary conditions. The sampling method is equal-interval sampling, and the number of sampling points for each boundary is 400. The configuration points in the configuration point set are sampled inside the simulation region to satisfy the partial differential equation. The sampling method is Latin hypercube sampling, and the number of samples is 6000. Step 3.3: Construct the loss function for the two-dimensional axisymmetric rotational model of the silicon melt: Based on the aforementioned two-dimensional axisymmetric rotation model of silicon melt, a loss function for a spatial factor physical information neural network is constructed. It consists of two parts: boundary point loss. Configuration point loss , Boundary point loss function as follows: （26） Location point loss function as follows: （27） （28） （29） （30） （31） Total loss for: （32） Step 3.4: Substitute the experimental parameters: Physical parameters of silicon melt: , ; Assuming the crystal radius is The radius of the crucible is The height of the crucible is Crucible rotation speed Crystal rotation speed ; Step 3.5: Train the spatial factor physical information neural network and make predictions: The experiment was conducted on the Python framework Tensorflow-GPU 1.15.

2. The computer configuration was as follows: CPU was AMD Ryzen 5 5600G with a frequency of 3.90 GHz, RAM was 16GB, and graphics card was RTX 2060 with 6GB of VRAM. The ADAM gradient descent method was used to optimize the network parameters, with an exponential decay rate of 0.

99. Training was stopped when the number of iterations reached 30,000 or the loss value fell below 1e-5. After the network training was completed, the radius of the spatial coordinates was... ,high As input to the network, the flow velocity in the simulated region is predicted.

Images