Kinetic Monte Carlo simulation of formation of microstructures in liquid droplets

We study the deposition of indium droplets on a glass surface and the subsequent formation of silicon microcrystals inside these droplets. Kinetic Monte Carlo methods are used to analyse the influence of growth temperature, flux of incoming particles, surface coverage, and in particular an energy parameter simulating the surface tension, upon the morphology of growth. According to the experimental conditions of crystallization, a temperature gradient and diffusion in spherical droplets are included. The simulations explain the formation of silicon crystal structures in good agreement with the experiment. The dependence of their shape and the conditions of formation on the growth parameters are investigated in detail.


Introduction
The formation of liquid droplets of indium on a molybdenum surface from the vapour phase, and the subsequent growth of silicon crystallites inside these droplets by a vapour-liquidsolid mechanism (figure 1) has recently been demonstrated [1] (figure 2). This novel growth method of semiconductor microstructures could be of considerable interest, e.g., for photovoltaic applications. Due to the µm scale of the objects and the different materials involved, the temperature distribution of the substrate surface heated by a lamp is experimentally not accessible. Therefore, the modelling of the process is an essential tool to improve the understanding of the influence of certain process parameters on droplet size and distribution and on the growth of crystallites inside the droplets. This should lead to a more precise control of the growth.
Kinetic Monte Carlo (KMC) simulation is an efficient and established tool to obtain information about the statistical behaviour of growth on an atomic scale. Preferentially, KMC algorithms [2,3] are used for understanding growth kinetics in epitaxial systems, e.g. [4]- [9]. However, to the best of our knowledge, such KMC methods have not been applied to investigate the selective growth of silicon microcrystallites inside liquid droplets deposited on an amorphous substrate. In this paper we consider this situation, and develop novel KMC schemes appropriate for that purpose. We present both a model for the formation of liquid droplets, including an energy parameter simulating the surface tension, and a model for crystallite growth inside these droplets, taking into account the solute diffusion.
These models will be implemented in a KMC scheme in order to investigate the influence of the parameters upon growth.

Thermodynamic aspects
To obtain an understanding of the thermodynamics of liquid droplet formation on a substrate, we discuss the behaviour of free enthalpy (Gibbs free energy), which describes the nucleation energy. We assume that the droplets have the shape of a sphere intersected by the substrate plane. The free enthalpy is given by where O is the surface area of the nucleus, r is the radius of curvature, V is the volume, σ is the surface tension, n is the number of atoms, T is the temperature and S is the oversaturation [10]. In the homogeneous case, resulting in the formation of a sphere, the free enthalpy has a maximum at the critical radius r c . Denoting the atomic volume by v, we obtain for this case  and the critical radius is given by In the case of heterogeneous nucleation on a planar substrate, the free enthalpy depends additionally on the contact angle φ between the substrate and the liquid droplet. Then the maximum free enthalpy condition, as used above, results in where σ In and σ Mb are the surface tensions of the indium droplet and the molybdenum substrate, respectively, By eliminating kT ln S appearing in both equations, we can compare the results of homogeneous and heterogeneous nucleation without knowledge of the oversaturation which is difficult to measure. For the experimentally determined contact angle φ = 5π/12; we get a ratio of approximately 4.3 for the critical radius (figure 3).
We simulate the formation of indium droplets and the crystallization of silicon inside by means of a continuous time KMC method. The model provides information about the behaviour of growth depending on the main parameters, i.e., temperature, flux, and time. In our simulations we employ a cubic lattice, 3 on which the atoms move. These are the elementary processes in our system. To model the deposition and the diffusion on a surface, we use ballistic deposition in the framework of a solid-on-solid simulation [11], neglecting desorption and bulk diffusion. The probability of a diffusion process is described by an Arrhenius equation, where ν 0 is the attempt frequency, k B is Boltzmann's constant, and E D is the diffusion barrier which differs for indium droplet formation and silicon crystallization. The timesteps in our event-based KMC algorithm are obtained by summing over all possible processes where P(µ → ν) denotes the hopping probability from site µ to site ν according to (7). Diffusion processes to any one of the four nearest-neighbour sites are considered. Material is deposited randomly corresponding to the flux F .

Droplet formation
The first stage in our KMC simulation, as in the experiment, is the formation of the droplets. To specifically simulate the behaviour of a liquid, we have to introduce a new parameter. Depending on the surface energies of the involved materials, a liquid droplet with a specific contact angle forms on the substrate. For a three-phase system with a solid substrate surface, the droplet formation is described by Young's equation, where σ 1 , σ 2 , and σ 21 are the surface tension between solid and vapour phase, between liquid and vapour phase, and between solid and liquid, respectively. These parameters determine the contact angle of the cap-shaped liquid droplet. That is the macroscopic point of view. In our mesoscopic KMC simulation of three-dimensional growth, it is essential to take account of the surface energy in an appropriate way. As a result of the surface tension, an atom on the surface of the droplet behaves differently according to its position: to form the droplet, the atoms at the base are forced to move upward, as indicated by arrows in figure 4. If the atoms are near the pole, they tend to move downward. Thus we have preferred directions of diffusion. In our KMC model this effect is simulated by different energy barriers, so a simple way to consider this effect on an atomic scale is to introduce a barrier for edge diffusion. For a solid, such a barrier is the Schwöbel barrier [12]. But in contrast to the solid state, here we have to take into account the directional dependence in the liquid case. Therefore we introduce an energy barrier E O and set it equal to E O up for upward diffusion, and equal to We restrict our model to constant barriers E O up and E O down neglecting all position dependence of the surface tension. But the main effects like the initiation of three-dimensional growth can be described appropriately. Within our cubic model, this new parameter cannot reproduce the spherical macroscopic morphology of the liquid, since cubic symmetry leads to square-shaped rather than circular bases of the equilibrium droplets, but the essential dependence of the droplet size upon the growth parameters is correctly described. The complete diffusion barrier in the Arrhenius rate (7) is now given by E D = E s + nE b + E O where E s and E b are the binding energies to the substrate and to the n nearest neighbours, respectively.

Crystallization
In the crystallization process we have different conditions. A solid (silicon) crystallizes in a liquid droplet. In our KMC model we confine the deposition and diffusion processes of silicon geometrically to sites within the volume of the liquid droplet, and consider only a single cap-shaped droplet. For the energy barrier of the crystalline growth kinetics we use E D = E s + nE b + E SW with the Schwöbel barrier E SW . There are experimental indications that a temperature gradient inside the droplet, and radial diffusion of silicon within the droplet may be of relevance. Therefore we modify the Arrhenius rates (7) in the following way: We add a position-dependent radial gradient of the temperature T , depending linearly upon the radial coordinate ρ with proportionality factor c T ( figure 5): T Figure 5. Scheme of radial temperature dependence with positive and negative gradients inside the droplet. Left, c T > 0, right, c T < 0. In addition to the temperature gradient we introduce radial diffusion of the silicon atoms inside a liquid indium droplet towards the centre of the spherical droplet. To create a simple but effective model, we simulate this diffusion in terms of a modification of the deposition process, which so far consisted of vertical deposition of the Si atoms upon the planar substrate surface within the droplet with a rate F . Now we assume that the Si atoms impinging upon the spherical surface of the droplets are deflected towards the centre of the sphere and follow a radial path until they reach a position on the substrate surface. This simplified model is valid if the radial diffusion is so fast that the radially diffusing atoms do not interact before they are deposited.

Droplet formation
Here we present the simulations of indium droplet formation on the substrate surface, which depend upon the growth parameters like temperature T , flux to the surface F , and surface tension parameter E O . In addition we observe the time evolution of the growth kinetics. The binding energies between nearest-neighbour atoms E b and to the substrate E s are kept fixed (table 1). First, we will simulate the initial formation of platelets as fundaments of growing droplets in analogy with submonolayer islands [13], setting E O = 0. With increasing temperature, the island size increases while the number of islands decreases, as can be seen in figure 6. This is in agreement with our experimental results as well as with KMC simulations in heteroepitaxial systems in the kinetically controlled regime [5,8]. The diffusivity increases with increasing temperature, and thus the probability of atoms getting attached to existing islands is enhanced. Therefore fewer and larger islands are formed at the same monolayer coverage. If we vary the flux F (in monolayers (ML) s −1 ) to the surface we observe a change in the island size ( figure 7). The flux is increased from F = 0.05 to 0.45 ML s −1 . Note that, at the same time, the coverage increases from c = 0.05 to 0.45, since the simulation time is kept constant. Higher deposition rates F lead to larger islands, because more material is deposited, and although the number of nucleation sites also grows with increasing F , the increasing coverage supplies enough material to obtain larger islands. An inspection of the temporal evolution shows that the islands grow larger in time, while their number remains approximately constant (figure 8).   This indicates that our KMC simulation is still in the kinetic regime, and the thermodynamic equilibrium size has not yet been reached [5]. Because of the high diffusivity, the deposited atoms tend to get attached to existing islands rather than forming new nucleation sites.
In addition to the experimentally adjustable parameters, i.e., temperature, flux, and growth time, we will now discuss the influence of the surface tension parameter E O . As shown in figure 9, the three-dimensional growth strongly depends upon the parameter E O . By increasing the energy barrier E O down and keeping E O up at a constant value, we find a threshold for the initiation of growth in higher monolayers. For low E O down the growth is restricted to the first monolayer. This is consistent with our expectation that the surface tension parameter induces three-dimensional growth. Like the Schwöbel barrier in solid growth systems, this surface tension parameter controls three-dimensional growth in a liquid system.

Crystallization
Next, we simulate the formation of silicon crystallites inside the liquid indium droplets. The parameters used in the simulations of crystallization inside droplets are presented in table 2. The influence of a radial temperature gradient inside the droplets is shown in figure 10. As can be seen in figure 10, a radial temperature dependence with negative gradient, i.e. temperature decreasing away from the centre of the droplet, favours the formation of islands at the border of a droplet. By choosing a positive gradient the formation of a three-dimensional structure is limited to the centre of the droplet. In regions with lower temperature we have reduced diffusion, and therefore agglomerations of atoms are more stable kinetically than in regions of higher temperature.
These simulations cannot explain the growth of crystallites restricted to the centre of the droplets without atoms sticking at the borders (figure 10). In order to explain the experimentally observed pyramidal microcrystallites centred inside the droplets, we have to include the radial diffusion, simulated by radially deflected deposition. Depending on the Schwöbel barrier, we can observe a change in shape from a truncated pyramid for lower E SW to a full pyramid as shown in figure 11. We deduce that the radial diffusion process inside the droplets is the decisive factor leading to pyramidal microstructures close to the centre of the droplets. The morphology is strongly dependent upon the Schwöbel barrier. That is in line with results from previous work, where the influence of a Schwöbel barrier upon epitaxial growth is discussed [14].

Conclusions
In this work we have presented a KMC scheme which is able to model a novel method for selective growth of microcrystals. It is based upon a two-stage Monte Carlo simulation of the formation of indium liquid droplets and the subsequent growth of silicon microcrystals inside these droplets.
Our simulations predict the dependence of the size and shape of the microstructures upon the essential growth parameters and conditions, which might help to optimize them. From Monte Carlo simulations of the indium droplets we have found a strong dependence of their size and number upon the growth temperature, the flux to the substrate, and the growth time. The size of the droplets, which eventually limits the size of the microcrystals, is smallest for low temperature, and for low flux and coverage. We have also introduced a parameter which simulates the droplet surface tension and which is necessary to correctly describe the onset of three-dimensional growth.
From Monte Carlo simulations of silicon crystallization inside the droplets we have obtained an understanding of the conditions under which the microcrystals grow at the centre of the droplets. We have demonstrated that a radial temperature gradient inside the droplets influences the areas and the distribution of mass of the three-dimensional silicon microstructures. However, the most essential parameter for pyramidal growth close to the centre of the droplets has been shown to be radial diffusion of the silicon atoms inside the liquid indium droplets.