Prediction of operating dynamics in floating-zone crystal growth using Gaussian mixture model

ABSTRACT We have applied a Gaussian mixture regression to the prediction of operation dynamics in floating zone crystal growth as an example of a materials process. From only five demonstration trajectories, we successfully predicted the operating dynamics using the Gaussian mixture model with better precision than obtained by using linear regression or neural networks. The current results indicate that the Gaussian mixture regression is suitable for predicting the operation dynamics of materials processes in which it is preferable to avoid large changes from stable operating conditions. Furthermore, precise prediction by the Gaussian mixture regression will lead to the optimization of operation trajectories and automatic control of materials processes. GRAPHICAL ABSTRACT


Introduction
The advances and application of informatics has ushered in a new era in materials processing [1]. Materials processes have been upgraded, automated, and efficiently optimized by applying various kinds of informatics algorithms such as Bayesian optimization, and neural networks [2][3][4][5][6][7][8][9]. Application of Bayesian optimization to materials processes has successfully reduced the number of trials for optimization and suggested better operation conditions for high performance and yields [10][11][12][13][14][15]. Surrogate modelling of materials processes by neural networks can efficiently search for the optimal conditions, [2,[16][17][18]. Although the parameters for most materials processes are independently given, which means that a set of input parameters gives a single output, some materials processes are manually controlled according to the information obtained during operation. For example, in crystal growth by the optical floating zone (FZ) method, an operator monitors the status of the melt in a furnace by a camera and changes the input parameters to maintain suitable conditions for singlecrystal growth (details of crystal growth by the FZ method are described later). In this case, it is necessary to conduct a huge number of trials when conventional Bayesian optimization is applied, and massive numbers of simulation results are required to construct a surrogate model by neural networks [18]. CONTACT  In other research fields such as robotics, chemical engineering and aerospace engineering, trajectory optimization has been studied for a long time for the development of various kinds of methods [19][20][21][22][23][24]. Commonly, trajectory optimization begins with predicting the dynamics of a system from demonstration experiments. In the present study, we predicted the dynamics of crystal growth by the optical FZ method as an example using Gaussian mixture regression (GMR) with consideration of the characteristics of the trajectory of materials processing in manufacturing.

FZ method and the characteristic of its trajectory
The FZ method without using a crucible was developed to grow high-purity silicon single crystals without the molten zone coming in contact with any foreign materials [25,26]. Silicon wafers up to 200 mm diameter are manufactured by the FZ method using RF heating, for the production of semiconductors. Not only silicon crystals but also numerous other crystals such as semiconductors [27], metals [28], alloys [29], intermetallic compounds [30][31][32], oxides [33][34][35][36], borides [37,38], carbides [39] and silicides [40][41][42], especially those with a high melting temperature, have been grown by the optical FZ method [43]. Figure 1 shows a schematic illustration and photograph of a typical optical FZ apparatus. The FZ crystal growth process is usually monitored by a video camera installed in front of a hole in the focusing mirrors. Two polycrystalline rods are mounted so that their tips meet at the focal point of the ellipsoidal mirrors. Halogen or xenon lamps are placed at the other focal points of these mirrors. The optical FZ crystal growth begins with melting the tips of the polycrystalline rods and bringing them together to create a liquid melt called the floating zone between the lower (seed) rod and the upper (feed) rod. After the floating zone is created, the zone begins to move upwards by moving the seed and feed downward (or by moving the mirrors upward), the liquid melt cools, and the crystal finally grows on the seed rod. During crystal growth, the rods rotate either in the same direction or opposite directions at certain speeds. An operator controls the input parameters such as the power of the lamps and the movement speed of the feed so that the liquid melt is not separated or dripped off (the speed at which the seed moves, which is assumed to be crystal growth rate, is usually fixed during crystal growth). If the lamp power is too low, the tips of the seed and feed collapse and the liquid melt can easily separate. On the other hand, if the power is too high, the liquid melt has difficulty maintaining the shape of the floating zone by its surface tension, and it can easily drip off. In this situation, the operator attempts to grow a single crystal by forming a certain shape in which the crystal diameter is first reduced (called 'necking') and then gradually increased to the necessary size.
The reason why the operator controls the input parameters according to the state of the liquid melt is that the same operation trajectories never yield the exact same melt conditions leading to different results due to the slight deviations, such as the condition and configuration of the feed and seed rods, quartz tube, halogen lamp and other components in the FZ furnace. Although the response of the melt conditions to the input parameters is almost the same, the slight differences in melt conditions are accumulated and the results become gradually different in each operation. This is true of other materials processes that are manually controlled by human operation according to the monitored conditions.
In summary, the trajectories of the crystal growth process by the FZ method and other materials processes manually controlled according to monitored conditions have the following characteristics that should be reflected in the algorithm: (1) The number of demonstrations (operation trajectories for model construction) is limited. (2) Each trajectory is similar, and is sufficient to predict the dynamics within a limited parameter space. Since materials processes such as FZ crystal growth experiments generally take a long time for preparation and it is difficult to generate 'big data', a sampleefficient algorithm is necessary. On the other hand, it is not necessary to predict the dynamics in all parameter space because for these materials processes, especially for mass production, it is preferable not to deviate greatly from the stable operation trajectory. Considering these characteristics of the trajectories for these materials processes, prediction of the dynamics using GMR is a good choice since GMR can efficiently construct nonlinear dynamics within a neighbourhood around the space covered by the demonstrations, which was successfully demonstrated in learning the non-linear dynamics of robot motion [23,44].

Emulation of FZ crystal growth trajectory for learning
Since the real dynamics of FZ crystal growth are unknown, prediction of the dynamics was validated by virtual experiments using an emulator of the FZ crystal growth experiment with a given set of dynamics. Note that whether the given dynamics are consistent with the real dynamics is not relevant to validation of the prediction. Thus, we first made an emulator imitating the FZ crystal growth experiment in which the lamp power (P) and movement speeds of the feed (u) and seed (v) can be controlled as input parameters, and the height of the liquid melt (h) and the diameter of the grown crystal (d) were determined according to the following dynamics: where a is a constant and d 0 is the diameter of the feed.
Here we assumed P is independent of time and the volume of a feed, melt and a crystal are represented by the 2-dimensional area size. The left-hand side of Equation (2) represents the time differential of melt area, and the right-hand side represents the speed of increment of melt by changing the feed into the melt and the speed of decrement of melt by changing the melt into the crystal. The value of h, which is independent of time, is determined by P, and the value of d depending on time is determined by Equation (2). The system states of the FZ crystal growth with the fixed values of P and u were described by d and _ d. Figure 2 shows a schematic illustration of the model for FZ crystal growth and its parameters, and a screenshot of the emulator program with a graphical user interface (GUI). The program can emulate the FZ crystal growth experiments and acquire the trajectories of input and output parameters. We manually prepared 10 trajectories for learning (training) (#1~5) and validation (testing) (#6~10) of the prediction of the dynamics with the actual FZ crystal growth experiment in mind. In the present study, we fixed both the movement speed of the feed and the lamp power as 1.0. Figure 3 shows the 5 trajectories for learning in which we aimed to create the following crystal shape: (1) At the initial stage (100 < t), the diameter of the crystal is kept the same as that of the feed rods (d = 1.0) (2) Then (100 < t < 300), the diameter is decreased to 0.1 corresponding to the "necking" process. (3) After the necking process (300 < t < 500), the diameter is increased to 1.0. (4) Finally, the crystal diameter is kept at 1.0 until the crystal growth experiment is over (500 < t < 1000).
Since the trajectories were manually prepared using the emulator program, they were different from each other and did not completely fulfil the above aim.

Prediction of dynamics by GMR
For prediction of the dynamics of the FZ crystal growth, we utilized GMR, which was reported to efficiently learn non-linear dynamics for robot motions . Here we described the prediction of the dynamics by GMR assuming FZ crystal growth. The state of the liquid melt at time (t+1), which is composed of the height and diameter of the liquid melt and described as x t +1 = (h t + 1 , d t + 1 ), is assumed to be determined by the state of the liquid melt (x t ) and inputs composed of the lamp power and moving speeds of the feed and seed (y t = (P t , u t , v t )) at time t as follows: The function f was represented by a locally linear function as follows: F z and F 0 are the coefficient matrix and vector, respectively. In linear regression (LR), these coefficients are determined independently of z t so as to minimize the sum of the squared errors. In GMR, however, these coefficients are estimated by a Gaussian mixture model (GMM), which is a non-linear combination of a finite set of Gaussian kernels constructed from the set of demonstration trajectories. GMM defines a joint probability distribution P(z t , x t +1 ) of the current state, input and next state over a set of demonstration trajectories as a mixture of a finite set of n Gaussians G(z t ,x t +1 | μ k , Σk) (k = 1, 2, . . . , n) with μ k and Σk being the mean value and covariance matrix of the k th Gaussian as follows; where π k is a prior. Gaussian mixtures were optimized to fit the demonstration trajectories by an expectation maximization (EM) algorithm. The initial values of μ k , Σ k and π k were set. In the present study, the mean value and covariance matrix for a demonstration trajectory, and a uniform prior distribution was used as the initial the initial values of μ k , Σ k and π k . The values of μ k , Σ k and π k were optimized to fit the demonstration trajectories by an expectation maximization (EM) algorithm [46] . Taking the posterior mean estimate of P(x t +1 |z t ), one can estimate the function f as a nonlinear sum of linear dynamics as follows [45] : Equation (9) shows that the next state (x t+1 ) is described by the current state (x t ), the inputs (y t ) and the optimized Gaussian mixtures by EM algorithm.
In the present study we set the number of Gaussian mixtures as n = 50 for prediction of the dynamics of FZ crystal growth as a result of parameter search.

Results and discussion
Among the 10 trajectories prepared by the emulator, the GMM was trained from the 5 demonstration trajectories (Train-1~5) and validated by the other remaining 5 trajectories (Test-1 ~ 5). Figure 4 shows the Gaussian mixtures optimized to the 5 demonstration trajectories. The Gaussian mixtures are formed to cover the demonstration trajectories, which implies that the dynamics predicted by the GMM reproduce actual dynamics near the demonstration trajectories. Note that some Gaussians have very narrow width indicating the overfitting to the training trajectories. However, the overfitting is acceptable since we aim to make a model which can precisely predict the dynamics near the stable operation trajectories. Figure 5 shows the trajectory and absolute error of the diameter of the grown crystal (d) predicted by the GMM in comparison with the actual trajectories as well as the trajectories predicted by linear regression and a neural network (NN). The GMM  reproduced the actual trajectories better than the neural network and linear regression. From the viewpoint of crystal growth, the trajectory predicted by using GMM is acceptable except for Test-4, in which the values of d are deviating to the negative values at around t = 400. The negative values of d correspond to the melt is separated. Most of the trajectories predicted by NN tend to deviate from the true value when the value of d decreases at around t = 150 or increases at around t = 400, and the value of d tends to increase rapidly, which make the quality of single crystal worse. The trajectories predicted by linear regression seems to much larger deviation than GMM. Figure 6 shows the mean absolute error (MAE) for the trajectories predicted by GMM. Compared to the other methods, the GMM can more accurately predict the trajectories. Figure 7 shows the relative errors in the value of d in the next time step predicted from the value of d and v in the current time step accommodating the training trajectories. Relative error was calculated as the difference between the value of d in the next time step calculated from Equation (2) and GMM. The prediction is sufficiently precise, and the relative errors are less than 0.01 near the training trajectories. The non-linear local dynamics of the operation trajectories for the FZ crystal growth were accurately predicted by the GMM near the training trajectories. Although it is difficult to strictly define whether the trajectory in operation is near the training trajectories or not, we can notice the accuracy of the prediction of the dynamics during the operation since we can compare the actual state and predicted state from the previous state and the GMM.
In the present study, we successfully demonstrated that the GMR is an efficient method for the prediction of operation dynamics for FZ crystal growth. From only five operation trajectories, the GMM can predict the actual trajectories of FZ crystal growth with relatively high accuracy. Although the current result was demonstrated by operation trajectories made by an emulator program, it is expected that the dynamics of the actual FZ crystal growth can be well-predicted by the GMM. The present demonstration is the first step of the automation of materials process including FZ crystal growth with optimized trajectory.

Conclusion
We have applied the GMR to the prediction of FZ crystal growth, as an example of materials processes in which it is preferable to avoid deviating greatly from the stable operation trajectory. The GMR can predict the dynamics of the FZ crystal growth process from only five operation trajectories better than neural networks and linear regression. Although there is still a gap between the demonstration in the simple emulation program and its application to the actual materials processes which is complex and probabilistic, the present results imply the feasibility of the prediction of the actual dynamics for the FZ crystal growth process, which will lead to the optimization of operation trajectories and automatic control of materials processes.  Author contributions SS conceptualized the basic idea and SH conceived the application to the materials process. RO constructed algorithm and programs for analysis under the guidance of SS, with the assistance of YT and in continuous discussion with all authors. The manuscript was written by SH in discussion with all authors.