The Effect of Disordered Substrate on Crystallization in 2D

TL;DR

This study investigates how amorphous, disordered substrates influence 2D crystal growth, revealing that disorder can sometimes enhance crystallization and identifying optimal growth temperatures through Monte Carlo simulations.

Contribution

It introduces a Monte Carlo simulation approach to analyze the impact of substrate disorder on 2D crystallization, highlighting conditions where disorder improves crystal quality.

Findings

Disorder affects local connectivity and transition rates in crystal growth.

Optimal growth temperature exists for both ordered and disordered substrates.

Disorder can beneficially influence crystallization under certain conditions.

Abstract

In this work the effect of amorphous substrate on crystallization is addressed. By performing Monte-Carlo simulations of solid on solid models we explore the effect of the disorder on crystal growth. The disorder is introduced via local geometry of the lattice, where local connectivity and transition rates are varied from site to site. A comparison to an ordered lattice is accomplished and for both, ordered and disordered substrates, an optimal growth temperature is observed. Moreover, we find that under specific conditions the disordered substrate may have a beneficial effect on crystal growth, i.e., better crystallization as a direct consequence of the presence of disorder.

Equations

$$p_{x_i,x_j}=\frac{e^{\beta \Delta E_{x_i,x_j}}}{{\sum }_j{e}^{\beta \Delta E_{x_i,x_j}}},$$

$$p_{x_i}=e^{-\beta E_{x_i}}.$$

$$E_{S_i}=\left(i_{th}cellperimeter\right)/a$$

$$E_{In{t}_i}=\frac{{\sum }_j{b}_{x_j}}{{\sum }_{j}1}$$

$$E_{In{t}_i}=\frac{{\sum }_j{b}_{x_j}{f}_{x_i,x_j}}{{\sum }_j{f}_{x_i,x_j}}$$

$$NWD=\frac{{\sum }_i{\sum }_j{b}_{x_j}{f}_{x_i,x_j}}{{\sum }_i{\sum }_j{f}_{x_i,x_j}}$$

$$b_{x_j}=\{\begin{array}{ll}1,& \text{mbox}jisfull\\ 0,& \text{mbox}jisempty\end{array}$$

Full text

I introduction

Controlling self-organization of two dimensional crystals is very instrumental due to their potential applications in low-dimensional semiconductors and optoelectronic devices Novoselov and Neto (2012); Putkonen et al. (2005); Wu et al. (2019); Chen et al. (2011a). For example, self-organization of a quantum dot Nurminen et al. (2000), graphene sheet Tetlow et al. (2014) or transition metal dichalcogenides Wu et al. (2019) affects their optical, electrical, magnetic and mechanical properties. Two-dimensional crystals can be also assembled in three-dimensional hetero-structures that do not exist in nature and have tailored properties Novoselov and Neto (2012).

This explains the considerable interest in studying metallic and semiconductor crystal growth. A promising technological development is to grow two-dimensional crystals on a substrate. The most studied two-dimensional crystal is graphene. Graphene has a wide range of interesting properties such as extremely large charge motilities, unprecedented mechanical strength, remarkable heat conduction at room temperature and a uniform absorption of light across the visible and near-infrared parts Novoselov and Neto (2012). Growth rates of two-dimensional crystals, such as graphene are limited by surface diffusion or attachment, that is, by the rate at which deposited atoms jump towards the growing crystal or the rate at which these atoms attach to the crystal. Different parameters such as flux of atom deposition on the surface, temperature, and the nature of the substrate on which the crystal grows, influence the speed of crystal growth as well as the quality of the obtained crystal Wu et al. (2019); Levi and Kotrla (1997); Tetlow et al. (2014); Jourdain and Bichara (2013). Most theoretical studies that are aimed at providing a better understanding of a crystal growth have used Molecular Dynamics (MD) L. Meng and Ding (2012); Tetlow et al. (2014); Shibuta and Maruyama (2003) or Kinetic Monte Carlo (KMC) methods Tetlow et al. (2014); Nurminen et al. (2000); Zhang et al. (2004). MD simulation (MDS) describes the diffusion of atoms according to the forces acting on them. The MDS approach involves the description of all the atomistic details Tetlow et al. (2014); Tang et al. (1992); Ratsch and Venables (2003); Biehl (2005a) and crystal growth is associated with many local events. In this study we concentrate on KMC simulations, that in general, go beyond atomistic details. These methods are based on coarse-graining of “molecular” growth models. Coarse-grained methods take into consideration “important” processes and neglect other details Tetlow et al. (2014); Battaile and Srolovitz (2002); Chatterjee and Vlachos (2007). For example when describing crystallization on substrate one can look at jumps of absorbed atom (adatom) between neighbouring cites of immovable lattice that represent the substrate. In reality, the substrate rearranges and reacts to adatom movement but as long as the substrate movement has a negligible impact on the adatom energy, this movement can be neglected Biehl (2005a). This kind of KMCs describing crystal growth as adatom moving on a lattice, is called solid on solid (SOS) simulations Chatterjee and Vlachos (2007); Pyziak et al. (2004); Biehl (2005a) . In SOS simulations, adatoms are deposed randomly on a site and then jump to a neighbor site according to a microscopic model (see Fig. 2).

The effects of substrate on nucleation and crystal growth have been previously investigated Nurminen et al. (2000); Chen et al. (2011a) using SOS simulations. These models usually use a square lattice and add energetic barriers to specific sites. Specifically, in Nurminen et al. (2000) a patterned substrate was considered for the nucleation process. In order to introduce a substrate pattern, the lattice is divided into square-shaped domains and the adatom energy is a function of its location on the lattice. Nurminen et al. Nurminen et al. (2000) found that patterned substrate affects the nucleation process on the substrate. In Chen et al. (2011a) point defects were introduced and it was found that this addition, at a certain temperature of the substrate material, improves the nucleation and increases the average crystal size. While these studies show that patterning, or addition of local disorder to the substrate, have a good impact on crystal size (increased average island size) it is not clear what happens when the underlying lattice is not ordered, i.e., amorphous substrate. Moreover, experiments K. (2013) show that substrate composition and surface crystallinity also influence crystal growth. Recently, amorphous substrates such as liquid substrate, have been experimentally used for crystal growth, the lack of a crystallographic substrate have been observed to have a good impact on crystallization, i.e. enlarge crystal size. Zeng et al. (2014); Boeck et al. (2017); Zhang et al. (2018). When taking into account this experimental observation, it becomes interesting to explore the effect of substantial substrate disorder on crystal growth. Of special interest is the case when the substrate geometry is sufficiently altered and can’t represent an ordered lattice anymore, i.e., the case of amorphous substrate.

The effects of temperature on crystal growth has been previously addressed in L. Meng and Ding (2012); Zhang et al. (2004); Chen et al. (2011a); Meixner et al. (2001); Pyziak et al. (2004). In L. Meng and Ding (2012) an optimal growth temperature for graphene growth on Ni is found by molecular dynamics simulations. An optimal temperature, at which the surface roughness is minimized, was recently obtained Zhang et al. (2004) using a 3D KMC simulation. Experiments show that at a relatively low temperature, temperature increase has a beneficial impact on crystal growth K. (2013); Loginova et al. (2008); et al. (2013); Chen et al. (2011b) while, at higher temperatures, an opposite effect is found Rafik et al. (2012).

In this work crystallization is studied on two types of latices; square (ordered) and random (amorphous). We use a vectorizable random lattice (VRL) in order to simulate an amorphous substrate. VRL is a lattice with sites that compose a set of randomly chosen points with uniform distribution Moukarzel and Herrmann (1992a). In Bunimovich and Khlabystova (2001) a VRL was used for the study of Lorentz gas. In Schwarcz et al. (2016) one of us used a VRL in order to simulate a semi-solid substrate (agar substrate) that does not posses a crystal order. Here we present two SOS models. Model A in which the interaction energy between adatoms depends only on the number of neighboring adatoms and model B where the interaction energy between adatoms depends on the local substrate geometry. We compare the results obtained on both substrates, for both models, and the effect of temperature is also presented. Interestingly enough, it is found that not only that there is an optimal temperature, but also the amorphous nature of the substrate has a beneficial effect on crystallization.

This paper is organized as follows. In Sec.  (II) the models and computational methods are presented. The effects of temperature variation on crystallization and comparison of ordered and amorphous substrate are presented in Sec.  (III). In Sec.  (III) we also discuss the possible explanation for the beneficial impact of amorphous substrate. The summary is in Sec. IV.

II Models and methods

In order to describe crystal growth on an amorphous substrate by a SOS model one first needs to create a lattice that represent an amorphous substrate. In order to achieve this goal we use a Vectorizable Random Lattice (VRL), VRL is a lattice that consists of random sites that are uniformly distributed in space  Moukarzel and Herrmann (1992b); Schwarcz et al. (2016). During the construction of a VRL following steps are performed: I: A square lattice of $d\times d$ cells is defined, this lattice is called the reference lattice, as shown in Fig. (1) (a). II: A random point is chosen in each cell of the reference lattice (with uniform distribution), while keeping a minimum distance $\delta$ between points as presented in Fig. (1) (b). Any two points that are closer than $\delta$, are re-allocated. $\delta$ controls the degree of randomness of the lattice. These points constitute the VRL sites. III: The Voronoi cell for each lattice site is produced Fig. (1) (c). Voronoi cell is defined as a set of all points that are closer to a given lattice site, than to any other lattice site Vieira et al. (2001). The simulations presented in this work were carried out on a reference lattice of $10^{4}$ cells ($d=100$), the edge length of the reference square lattice was set to be $1$ while $\delta$ is set at $1/10$. Special care was taken in order to enable periodic boundary conditions of the VRL. In the following we describe the dynamics of adatoms on top of a VRL and an ordered square lattice.

Adatoms are allowed to move on the top of the disordered lattice, but first we depose adatoms on the substrate i.e. fill part of the lattice cells with adatoms. Therefore during the first 30 simulated time steps of the simulation 100 adatoms are randomly deposited on VRL sites. If the randomly chosen site is already occupied, the new adatom is added to a neighboring site (in the case of random lattice the neighbor site with the bigger common edge is chosen). After this deposition period, adatoms are no further added. During the simulation, each adatom can jump to one of the neighboring sites of the disordered lattice, i.e. sites that have a mutual edge with the current site that contains the adatom. If an adatom jumps to a site adjacent to another adatom it nucleates to form a new crystal. If it jumps to a site adjacent to a crystal it aggregates to the existing crystal. Adatoms that belong to a crystal border can detach from the crystal or change their emplacement in the crystal. From an energetic point of view, neighboring adatoms has interaction energy therefore if an adatom have a neighbor its jumping probability decreases. On the other hand the adatom probability to stay attached to its neighbors increases. Fig.  (2) summarizes the dynamical process on the lattice, for simplicity the lattice is presented as symmetric. We allow at most one adatom per lattice site.

The probability of a randomly chosen adatom to jump is defined as a combination of two probabilities: the probability to leave the original site and the probability to occupy a new site. First, all the potential locations $x_{j}$s that the adatom may reach from its initial site $x_{i}$ are specified. These potential locations are the sites that have a common edge with $x_{i}$. For each of those sites the transition probability $p_{x_{i},x_{j}}$ is calculated

[TABLE]

where $\beta=\frac{1}{k_{b}T}$, $E_{x_{i}}$ is the adatom energy in position $x_{i}$, $\Delta E_{x_{i},x_{j}}=E_{x_{j}}-E_{x_{i}}$, T is the temperature and $k_{B}$ is the Boltzmann constant. A destination site is randomly chosen in accordance to $p_{x_{i},x_{j}}$. The probability, $p_{x_{i}}$, that an adatom will attempt to leave its position $x_{i}$ is defined as

[TABLE]

Where $E_{x_{i}}>0$ is the adatom energy. $E_{x_{i}}$ is defined as a sum of two terms $E_{Int_{i}}$ and $E_{S_{i}}$ that represent the interaction with the substrate ($E_{S_{i}}$), and the interaction energy contribution ($E_{Int_{i}}$) from all of the occupied nearest neighbors.

The substrate binding energy, $E_{S_{i}}$, is defined in direct proportion to the cell boundary size

[TABLE]

where $a$ is a numerical coefficient. $E_{S_{i}}$ is the local energetic well that dictates the diffusion properties of an isolated adatom. It depends on the lattice site geometry/size and is quenched, i.e. remains constant during the simulation. When the lattice is geometrically disordered, $E_{S_{i}}$ is responsible for introduction of this heterogeneity into the dynamics of the adatoms. In the simple case of an ordered lattice, the substrate binding energy $E_{S_{i}}$ is uniform for all the lattice sites and local diffusion is homogeneous without preferred locations. Cases with uniform quenched energetic barriers were studied in Biehl (2005b); Xu et al. (2017). In our model we assume that the substrate irregularity also affects adatom emplacement and that the adatom emplacement, in its turn, affects the interaction energy between adatoms. The definition of the interaction energy, $E_{Int_{i}}$, between the adatom $x_{i}$ and neighbouring adatoms $\{x_{j}\}$ is varied between the two models. In model A, for each cell every neighboring adatom has the same contribution to the interaction energy. In model B the interaction energy between adatoms depends on the local substrate geometry, specifically it depends on their common edge length $f_{x_{i},x_{j}}$. The interaction energy between two neighboring adatoms in Model A is defined as

[TABLE]

where $b_{x_{j}}$ is $1$ if $x_{j}$ is occupied and [math] otherwise, the sumation is over all the neighbouring sites. In model (B) the binding energy between neighboring adatoms is

[TABLE]

where $f_{i,j}$ (as previously mentioned) is the length of the boundary between $x_{i}$ and $x_{j}$. The summation is again over all neighbouring sites. In the specific case when an ordered lattice is used, Model A and B should provide exactly the same results. For a VRL there is expected to be a difference between the behavior of Model A and Model B, due to the fact that the boundaries of a Voronoi cell are random. VRL simulates an amorphous substrate, therefore for a VRLs the substrate-adatom interaction energy $E_{S_{i}}$ varies for different cells. In model A, for each cell there is a constant interaction energy between adatoms and their neighbors. For VRL the amount of neighbours for different cells is not uniform, this in turn cause to variations of $E_{Int_{i}}$ simply due to the fact that the summation in Eq.  (4) varies. In model B, on top of the diversification of Model A, the magnitude of interaction between two adatoms also varies according to $f_{x_{i},x_{j}}$ ( Eq. (5)). The reasoning behind introduction of this specific dependence of interaction energy between adatoms stems from the geometrical properties of the VRL. We want to take into account only the interaction with nearest neighbours. Number of nearest neighbours varies from site to site due to the disordered structure of the VRL. While for Model A the only thing that matters is the occupation of a nearest neighbour by adatom, in Model B the length of the mutual boundary with this specific site (i.e. $f_{i,j}$) is taken into account. The closer the centers of two VRL sites $i$ and $j$ the larger (on average) is their mutual boundary and vice versa.

Both models include three scaling parameters, $a$, $\delta$ and the edge length of the square lattice. Theses parameters can be adjusted to correspond to a given disordered situation and comparability between $E_{S_{i}}$ and $E_{Int_{i}}$ . The order of magnitude of $f_{i,j}$ is dictated by the edge length of the reference square lattice and $\delta$. When the length edge of the reference square lattice is one, increasing $\delta$ from $\epsilon(\to 0)$ to 0.3 leads to decrease of the variance of $f_{i,j}$. The distribution of $f_{i,j}$ is rather complicated and for high enough values of $\delta$ ($\delta\gtrapprox 0.5$) it stops to behave as a Gaussian centerd around its mean. Therefore we chose to use rather small $\delta=0.1$. The numerical coefficient $a$ that appear in Eq. 3 dictates the ratio of the interaction energy of adatoms ($E_{Int_{i}}$) and the interaction energy between the substrate and the adatom ($E_{S_{i}}$).

To summarize, at each time step a particle is chosen in a random fashion, then its probability to jump, $p_{x_{i}}$ , is calculated in accordance to Eq. (2). If the jump move is accepted, the destination site is chosen among all its empty neighbouring cells. An empty neighbour cell with a lower energy state, i.e. a neighbour surrounded by more adatoms, has better chance to receive the jumping adatom, in accordance to $p_{x_{i},x_{j}}$. The probability $p_{x_{i}}$ (Eq. (2)) of an adatom to exit from its original site, match to trap models where the depth of the trap depends on the energy of the original site Bouchaud and Georges (1990); Burov and Barkai (2011); Miyaguchi and Akimoto (2015); Burov (2017); Magdziarz and Szczotka (0181). The probability to reach a new site, as defined in Eq. (1), resembles barrier models Berne et al. (1998); Bouchaud and Georges (1990), where the transition probabilities between neighbouring sites are asymmetric.

III Results

In the simulation, crystallization on a substrate is reproduced by a two-dimension SOS model described in the previous section. The largest and the most defect-free crystals are the most desirable. In order to quantify the quality of a crystal, we need to take into consideration the crystal size and its uniformity. For this purpose we define an order parameter, the normalized weighted density (NWD): a sum over all cell edges that separate two occupied cells, normalized by the sum of all edges for all occupied cells,

[TABLE]

where

[TABLE]

The first summation (over $i$) is over all lattice cells in which there is an adatom, second summation (over $j$) is over the neighbors of cell $x_{i}$. A brief glance at NWD reveals that there is a resemblance between NWD and the interaction energy of all the adatoms. The reason for this resemblance lays in the fact that the adatoms interaction energy has a strong effect on crystal growth. A global perspective of the system shows that NWD increases with the compactness of a crystal: the higher is the number of occupied sites of a given adatom, the higher is NWD. In the case of irregular crystal morphology, such as protrusion, branches and pores, many adatoms have empty neighbors and NWD is expected to be relatively small.

Two reasons can cause NWD to be small; if the adatoms jumping probability $p_{x_{i}}$ is highly restricted, the adatoms will be stuck in an isolated site, if $p_{x_{i}}$ is too high adatoms will not stabilize even on sites with many neighboring adatoms.

III.1 Temperature effects

In order to study the effect of temperature, we vary $\beta$ which appears in Eqs.  (1) and (2). We observe that the effect of temperature on NWD has an inverse U shape, see (Fig. (3 (a,b) and 4 (a,b) ). At sufficently low $\beta$ (high temperatures) ,an increases of $\beta$ improves crystallization until an optimal $\beta$ is reached, further increase of $\beta$ (low temperatures) damages the crystal.

This observation of optimal temperature for crystal growth can be explained by following reasoning. At sufficiently low temperatures the adatoms mobility is restricted due to the fact that they need quite large activation energy in order to leave their current site, i.e. small $p_{x_{i}}$. Increasing of the temperature leads to increasing of $p_{x_{i}}$ that in its turn contributes to better crystallization, i.e., higher NWD. As the temperature is further increased, another aspect of the dynamics must be taken into account. Above a specific temperature threshold the difference between a neighbor site with many neighbors and a less favorable neighbor is small. From the Eq.  (1), we can see that $\beta$ decrease (T increase) causes a decrease in $\Delta E_{x_{i},x_{j}}$ impact, therefore, the probability that a diffusing adatom will reach the site with highest interaction energy declines when the temperature is increased. This behavior is reproduced in our model (Figs. 3(a,b) , 4(a,b) ) and was already observed in Nishinaga (2015). An optimal temperature for crystal growth was also observed experimentally by Rafik et al. (2012) and reproduced theoretically L. Meng and Ding (2012); Zhang et al. (2004); Chen et al. (2011a).

III.2 Effect of the Disorder

The NWD is also used in order to study the effect of the disorder, i.e., amorphous lattice described by the VRL. For both definitions of the interaction energy, Model A and Model B, the comparison between ordered (square) substrate lattice and VRL with $\delta=0.1$ was performed ( Figs. 3(a,b) and 4(a,b)). In Figs. (3 (a,b)) and (4(a,b)) the behavior of NWD for ordered and disordered cases is presented. Surprisingly enough, the effect of the disorder on crystal-growth is non-negative and can also be positive.

Specifically, for Model B, the NWD parameter is higher for the amorphous case, for any given temperature. The growth of the crystal with time, as described by NWD, is slow in both cases (ordered and disordered) and after sufficient temporal period appear to grow logorithmically (Fig. 5(a-d)). The growth on the disordered substrate appears to be a little bit faster, for Model B, when compared to the ordered substrate. When Model A is considered the situation is somewhat more complicated. In Fig. (3 (a)), there is a small preference for the ordered substrates for most temperatures. But for longer times, Fig. (4(a)),this seems to change and for sufficiently low temperatures the disordered substrate starts to become more beneficial, as compared to the squared lattice. Indeed, when comparing the growth with time of NWD Fig. (5 (a,b)), we see crossover between the ordered and disordered substrates. In both cases the growth is logarithmic but the growth for the disordered case is faster.

We must note that our results describe the stage when the crystal is still growing and local deformations are still in place. In Fig. (5 (a,b)), even after $\approx 10^{6}$ simulation time steps, the steady state, i.e. constant NWD, is not reached.

The observed benefit of introduction of disorder is rooted in the definitions of transition probabilities $p_{x_{i}}$, $p_{x_{i},x_{j}}$ and the local interaction energies $E_{S_{i}}$ and $E_{Int_{i}}$. The interaction with the substrate, $E_{S_{i}}$, is the same for all lattice points in the case of ordered substrate. This is due to the fact that the cell perimeter is constant. When dealing with a VRL, each cell perimeter is different. Only on average the perimeter equals to the perimeter of the ordered lattice. Consequently, on average, the interaction with the substrate is the same for the ordered and the disordered cases (Model A and Model B). If one considers the situation when an adatom is isolated (step 1 in Fig. 2), the probability to perform a jump is similar (on average) for the different models. We can state that $E_{S_{i}}$ is the minimal ”depth” of local energy trap. On top of $E_{S_{i}}$ additional quantity $E_{Int_{i}}$, that depends on occupation of neighbour sites, is added. $E_{Int_{i}}$ can obtain discrete set of $4$ equally spaced values in the case of ordered square lattice. This is graphically displayed in Fig. 6a, the energetic spectrum of a trap is composed of $5$ different sates. When dealing with the disordered case, the situation is different for Model A and Model B. First consider Model A, where $E_{Int_{i}}$ is defined by Eq. (4). Here the $E_{Int_{i}}$ also attains a set of discrete values, but the number of these values varies from site to site. For the case of $\delta=1/10$ the average number of neighbours is $6$. This means that on average, the spectrum of energetic trap in Model A is composed of $7$ different (and equally spaced) states (Fig. 6b). In Model B, according to Eq. (5), there is not only $6$ distinct values (on average) but also they are not equally spaced. Fig. 6c describes the energetic spectrum of traps in Model B. In this situation one can expect to encounter with many situations where there are several energetic levels bunched closely to $E_{S_{i}}$.

This perspective of randomization of ”band gaps” and their number can help intuitively understand the benefit of introduction of disorder. At some stage the homogenization of the crystal is due to local rearrangements of the adatoms (Fig. 2 step 4). In order for a rearrangement to occur the adatom must escape the local energetic trap. The benefit of energetic spectrum of the disordered case is due to existence of energetic states in the vicinity of $E_{S_{i}}$, as compared to the ordered case. This line of thought also supports the observation that the NWD value for Model B is higher as compared to Model A. Existence of ”bunched” energetic states in the vicinity of $E_{S_{i}}$ can facilitate local transformations of adatoms and homogenize the crystal. We observe an effect that is small to moderate. Nonetheless, it becomes evident that disorder can support crystal growth by facilitation of new energetic pathways.

IV Summary

The purpose of this work was to explore the effect of disordered substrate on crystal growth. In order to achieve this goal KMC simulations of SOS models (with disordered lattices) were performed. The geometrical disorder of the lattice influenced the interaction energy between adatoms moving on the lattice and affected the crystal growth process. Numerically computed behavior of normalized weight density (NWD) parameter shows not only that there is a preferred temperature for crystallization but also that the presence of geometrical disorder is beneficial. We suggest that the disorder affects the local energetic states associated with each lattice site. It creates energetic states in the close vicinity of the energy associated with free particles. Those states, that are not present in the case of ordered substrate, can produce additional (and more preferable) energetic pathways to local reformation and faster crystallization. It will be interesting to explore this effect further, especially since it is known that fluctuations of energetic levels can lead to effects such as stochastic resonance Gammaitoni et al. (1998); Burov and Gitterman (2016).