Introduction

Understanding the growth of molecular materials such as the prototypical molecular semiconductor fullerene C60 (refs 1, 2) on surfaces is an indispensable prerequisite for the rational design of complex nanomaterials from molecular building blocks, as well as for optimizing the performance in thin-film-based applications such as solar cells3,4,5 and organic light-emitting diodes6,7. So far, molecular self-assembly and growth8 has often been described by scaling laws to describe surface roughening and evolving island densities9,10. On a molecular level, a range of studies have elucidated the kinetics of diffusion and nucleation (see, for example, refs 11, 12, 13, 14, 15, 16) and the Ehrlich–Schwoebel barrier for interlayer transport across a molecular step edge11,17,18 (see Fig. 1). In the last decades, the energy barriers for atomic growth have been refined to take into account the local neighbourhood during multilayer growth, for example, by including concerted gliding of islands or by distinguishing between different step-edge orientations19,20,21,22. Yet to date, there is no organic compound for which even the ‘minimal’ set of the three parameters diffusion barrier, lateral binding energy and Ehrlich–Schwoebel barrier have been simultaneously quantified to describe multilayer molecular growth. Therefore, predictive simulations of the rate- and temperature-dependent morphology in molecular multilayer growth have so far been impossible, contrary to the situation for elemental atomic systems23,24,25 and colloids26,27,28. Importantly, C60 exhibits properties in between those of atoms and colloids, which makes it a test case of fundamental relevance. On one hand, its van der Waals diameter of 1 nm29 is closer to atomic dimensions than to the μm length scale of colloidal systems. On the other hand, C60 resembles colloids with its short-range nature of the effective centre-of-mass interactions30, which decay as −1/r9 with r being the centre-of-mass separation stemming from the averaged van der Waals interactions (approximately −1/r6) between the individual carbon interaction sites31. These forces between atomic, molecular or colloidal building blocks are of prime importance for kinetic growth processes, similar to their role in equilibrium phase behaviour and self-assembly32,33. For example, C60 lacks a stable equilibrium liquid phase30, contrary to most elemental atomic systems. C60 is therefore not only relevant for device applications, but also an important, fundamentally unique material bridging atoms and colloids. From the experimental side, a particular challenge in studying C60 growth is that post-growth changes make the interruption of this non-equilibrium process to image different growth stages potentially misleading. It is therefore essential to use in situ real-time techniques.

In this article, we employ the combination of specular X-ray growth oscillations34 with real-time diffuse X-ray scattering35,36 to simultaneously follow the vertical and lateral morphology during growth. Further understanding on a nanoscale level is provided by kinetic Monte Carlo (KMC) simulations of coarse-grained C60 molecules without internal degrees of freedom. Then, the three relevant parameters determined by a fit of the data are the Ehrlich–Schwoebel barrier, the surface diffusion barrier and the lateral binding energy (see Fig. 1). With these parameters alone, we achieve quantitative agreement with the experimental data, enabling us to predict the rate-, temperature- and thickness dependency of the film morphology. Moreover, our analysis demonstrates that the short interaction range of C60 as compared with atoms affects the relative heights of diffusion barrier and binding energy and results in comparatively long diffusion times. However, unlike the colloidal systems, C60 has a true energetic Ehrlich–Schwoebel barrier, rather than the pseudobarrier that colloids display26.

Figure 1: Surface processes in C60 growth.
figure 1

The diffusion barrier ED, binding energy EB and Ehrlich–Schwoebel barrier EES determine island nucleation and interlayer transport in multilayer growth. Included are numerical values determined by fitting the experiment using KMC simulations.

Full size image

Results

Experimental results for the layer-by-layer growth of C60 on mica

For a comprehensive understanding of the processes during growth, the surface morphology has to be measured on the molecular length scale with an experimental time resolution that is fast compared with the minute timescale of the deposition of a monolayer. Interrupting growth to take a series of real-space microscopy images can be problematic, as the kinetics can be altered. For our system of C60 on top of a closed first C60 layer on mica, this route is indeed impossible because of quick dewetting effects characterized by a time constant of ~10 min. Also, in situ low-energy electron microscopy unfortunately—while very successfully used in a range of studies37,38—cannot be applied due to charging effects on mica. Therefore, we use X-ray scattering that can be performed non-invasively during growth and yields time-resolved information about the layer formation. This is extracted through specular reflectivity measurements at the so-called anti-Bragg position of C60 (see Fig. 2a) corresponding to half the Bragg value of the C60(111) reflection. Lateral information is available through simultaneous measurement of the diffuse scattering (grazing incidence small-angle X-ray scattering (GISAXS)), giving information about the island distance (Fig. 2a).

Figure 2: Specular and diffuse X-ray scattering during C60 growth.
figure 2

(a) Scattering geometry: both the specular reflected X-ray beam and the diffuse scattering are detected. The two-dimensional scattering pattern contains both lateral (transfer momentum q‖) and vertical (q⊥) information on the surface morphology. (b) The specular X-ray reflectivity at the anti-Bragg point q⊥=0.38 Å−1 oscillates with increasing molecular exposure (time × growth rate) during growth of C60 on mica indicating layer-by-layer growth (T=60 °C). (c) The diffusely scattered intensity oscillates with the nucleation and coalescence of every layer and exhibits a characteristic peak-splitting Δq‖. The latter corresponds to the inverse average island distance, which changes with film thickness.

Full size image

The time-dependent specular X-ray reflectivity as a function of molecular exposure, which is time × deposition rate, is shown in Fig. 2b for growth at T=60 °C substrate temperature and a deposition rate of f=0.1 ML min−1. The anti-Bragg intensity oscillates with a period of two monolayers (ML) as the X-rays are reflected from consecutive C60 layers and alternately interfere destructively and constructively with an intensity modulation of up to 90%. Here, the diffusely scattered intensity can be neglected in an analysis of the specular reflectivity, as it represents <1% of the total intensity. The oscillations are indicative of a layer-by-layer growth and from the change in oscillation period, a variation of the sticking coefficient is deduced (see Methods). Only after the first three layers, one observes a damping of the oscillations, reflecting the onset of slight roughening. An additional discussion on the anti-Bragg intensity during the growth of the first monolayer of C60 on mica is given in Supplementary Note 1 and illustrated in Supplementary Fig. 1.

While the diffuse scattering is weak, it nevertheless contains important lateral information. Figure 2c shows a map of the diffusely scattered intensity as a function of q‖ and molecular exposure (see Supplementary Fig. 2 for a graph of the diffusely scattered intensity at a molecular exposure of 0.3 nm). In contrast to the anti-Bragg oscillations, the diffusely scattered intensity oscillates with a period of one monolayer. As the first molecules are deposited in a monolayer, the surface roughness and therefore the diffusely scattered intensity rises due to nucleation of islands. Eventually, as the islands coalesce, the roughness and diffuse intensity decrease again, before reaching a minimum for a smooth complete layer. For each C60 layer, the diffusely scattered intensity has two maxima along q‖, because the characteristic average island distance D causes an increase in the diffuse scattered intensity at Δq‖≈±2π/D (refs 39, 40).

From a crystallographic perspective, we find the established41 epitaxial order of C60 on top of mica(001) as confirmed by grazing incidence X-ray diffraction experiments shown in Supplementary Fig. 3 and explained in Supplementary Note 2.

KMC simulations of the growth process

To understand the morphological evolution on a molecular level, we employ KMC simulations, which are capable of describing the entire growth process of (coarse-grained) C60 molecules into a face-centred cubic (fcc) lattice. KMC models the growth as a stochastic process, in which the molecules adsorb with a constant net adsorption rate f=fadsorb−fdesorb. The molecules are treated on a coarse-grained level, that is, we do not take into account any internal (rotational or vibrational) degrees of freedom. This coarse-graining approach is supported by the fact that for the temperatures studied here, C60 rotates freely both in bulk crystals42 and in one-dimensional confinement43. Once adsorbed, a particle at site i then can diffuse to a neighbouring fcc site j via an activated process with Arrhenius-type rate ri,j. We follow the Clarke–Vvedensky bond-counting approach44,45, where the rate is defined as

The pre-factor v0=2kBT/h is chosen in accordance with previous KMC studies for atomic systems46,47,48, consistent with our coarse-grained description of C60 as a sphere. The total energy barrier for molecular hopping consists of a barrier for free diffusion, ED, and contributions determined through the local neighbourhood of the particle. The neighbour binding energy EB contributes with a number of lateral neighbours ni. The sum of ED and niEB then determines the lateral diffusion (si,j=0) and thus, the growth of islands. Other pre-factors to the neighbour binding energy have been suggested in literature19,24, which increase the diffusion rate of particles along island edges. As a consequence, the islands become more compact. In our C60 system, however, the islands are quite compact from the very onset of the growth (see Fig. 3). Therefore, the details of the pre-factor of EB do not significantly influence the results at the parameters considered. If a particle at site i crosses an up- or downward step to reach site j, an additional Ehrlich–Schwoebel contribution EES is added to the total energy barrier (si,j=1). As a result, a particle diffusing onto an island from an edge site with two neighbours has to overcome the activation energy ΔE=ED+2EB+EES, while a particle on the island has to overcome only ΔE=ED+EES to diffuse downwards over the island edge. The step-edge energy barrier used in our simulations is, by construction, an average energy barrier. For this, we recall that our energy barriers are exclusively gained by comparison with experiment, and that the experimental (X-ray scattering) data are intrinsically averaged in lateral direction. Therefore, we did not take into account the orientation of the step edge in this study. The KMC input parameters T (substrate temperature) and f (adsorption rate) are taken directly from experiment. The KMC simulations have been performed from the second layer onwards as we concentrate on the C60–C60 interactions and do not model C60–mica interactions. This strategy is justified, as we know from the experiment that the first C60 layer is completely filled and that there is no lattice strain; thus, we can assume a smooth C60(111) surface as initial surface in simulations. Furthermore, we assume defect-free growth without cavities or overhangs. We also note that we do not take collective diffusion mechanisms into account. Different concepts for collective diffusion have been suggested in the literature, one example being dimer shearing49. More recently, approaches have been suggested for shearing, reptation and concerted gliding of islands50. These phenomena are certainly worth studying in more detail, however, it would not have been possible to simulate the time and length scale required in our study if these effects were included.

Figure 3: Experimental and simulated measures of surface morphology.
figure 3

(a) Island density (inset: 2D island growth regimes as simulated by KMC; scale bar, 100 nm), (b) anti-Bragg growth oscillations and (c) layer coverages are shown as a function of the molecular exposure for a C60 film grown at T=60 °C and f=0.1 ML min−1. Parts b,c include data from an analytical growth model. (d) Maximal island density for the third layer for both a low deposition rate of 0.1 ML min−1 and a high deposition rate of 1 ML min−1 as a function of temperature. The KMC simulations have been performed from the second layer onwards. The confidence interval in a and the error bars in d are calculated from the systematic experimental uncertainties. For the complete morphology evolution during growth for T=60 °C, f=0.1 ML min−1 as well as 40 °C, 0.1 ML min−1 and 60 °C, 1 ML min−1 simulated by KMC, see Supplementary Movies 1–3.

Full size image

Energy barriers for surface processes in C60 growth

For the comparison of experiment and simulations, we use the time-dependent layer coverages from KMC simulations to calculate anti-Bragg oscillations using kinematic scattering theory51 (see Methods). The energy barriers ED, EB and EES (see equation (1)) are then adjusted until both the simulated anti-Bragg oscillations and island densities fit the experiment. Figure 3a,b shows experimental (black dots) and KMC simulation data (red solid line) for the island density and the anti-Bragg intensity for the temperature T=60 °C. The experimental island density is directly extracted from the data in Fig. 2c, using the average island distance D≈2π/Δq‖, assuming a hexagonal island arrangement (see also Supplementary Fig. 4 and Supplementary Note 3 for a comparison with real-space atomic force microscope data). Both experiment and simulation predict that the island density changes markedly during the deposition of each monolayer. Initially, in the nucleation regime, the island density increases. Then, lateral island growth sets in, where the island density stays constant. Finally, the island density drops again as islands merge in the coalescence regime. The inset in Fig. 3a shows the corresponding KMC simulation snapshots for the three growth regimes. In all cases, we observe compact island shapes in the simulations as well as in the experiments. A more detailed comparison of the morphology is given in Supplementary Note 4 and shown in Supplementary Fig. 5. The sequence of growth regimes is observed for the first five layers at each temperature and deposition rate employed.

As it is clearly seen from Fig. 3a,b, there is excellent agreement between the experimental and simulated data regarding the island density and anti-Bragg growth oscillations. The minima and the maxima in the island density, as well as the trend of decreasing density for the different layers (increase in island size), are clearly reproduced. The apparent increase in the island density in the fifth layer, which starts to differ slightly from the true island density, indicates the limits of our data analysis. The analysis takes into account only the islands in a single, currently growing layer, however, due to the roughening of the film, both islands in the simultaneously growing 4th and the 5th layer contribute to the diffuse scattering at that stage. The vertical layer filling and roughening are also highly consistent, as can be seen from the good agreement between experimental and simulated evolution of anti-Bragg intensity in Fig. 3b. As an independent confirmation of the KMC results, we have employed a mean-field analytical model for thin-film growth (see refs 52, 53, 54), the results of which agree with the layer coverages of the KMC simulations, as can be seen in Fig. 3c.

Even beyond the specific experimental parameters chosen in Fig. 3a–c, KMC simulations show a good agreement with the experimental findings for all studied rates (0.1 and 1 ML min−1) and the full experimental temperature range of 40–80 °C (see Supplementary Note 5 and Supplementary Fig. 6 for a comparison of 40 °C and 0.1 ML min−1). This is seen in Fig. 3d, where we compare the experimental and simulated values for the maximum island density in the third monolayer. In accordance with growth theories predicting a scaling of island density with deposition rate/diffusivity10,23, we find that the island density decreases for higher substrate temperature and lower deposition rate by an order of magnitude. Furthermore, KMC simulations correctly predict the change in island density by an order of magnitude when changing deposition rate and temperature. Notably, this comprehensive agreement of temperature-, rate- and time-dependent data was achieved with a physical model of surface processes that contains only three parameters for the nanoscopic energy barriers for diffusion, nucleation and step-edge crossing. The resulting values are ED=(540±40) meV for the diffusion energy, EB=(130±20) meV for the lateral binding energy and EES=(110±20) meV for the step-edge/Ehrlich–Schwoebel barrier (see also Fig. 1). For a more detailed discussion of the mutual correlations between energy parameters, see Supplementary Note 6.

Discussion

It is instructive to compare the self-consistent parameter set obtained in this study to energy values reported earlier. The height of the C60 Ehrlich–Schwoebel barrier (110 meV) is comparable to atomic systems, such as Pt/Pt(111) (80 meV)24 and is close to the value of 100 meV for C60 from recent density functional theory calculations by Goose et al.55 Our value for the binding energy, EB=130 meV, is smaller than that related to the minimum of the pair interaction potential of two C60 molecules, in particular the Girifalco potential, EC60–C60=270 meV, which has been derived theoretically56,57 and has recently been measured in atomic force microscope experiments58. There are several factors contributing to this difference: first, we are considering molecules close to a substrate, which has not been taken into account in refs 56, 57 but has already been shown to weaken the interaction58. Second, we are considering dense and thus strongly correlated systems, not two molecules in vacuum as assumed in refs 56, 57. Third, and maybe most importantly, our value for the binding energy has been obtained such that experimental data are fitted over a range of temperatures. It is well known that effective potentials (and thus binding energies) can strongly depend on the temperature59; thus our value has to be considered as a temperature average. Finally, we stress that our value for EB is very close to an estimate gained from the cohesion energy per neighbour of C60 in its bulk fcc crystal, EC=133 meV (1.6 eV is the total cohesion energy60,61 divided by the 12 bulk lattice neighbours). Regarding our value for the diffusion barrier (ED=540 meV), we note that this is significantly larger than the corresponding value derived from a potential landscape analysis, Epot=168 meV (ref. 62). This is likely due to the fact that in our KMC simulations, we do not consider all energy minima as lattice sites. Thus, the travelled distances across several minima are larger, leading effectively to a larger barrier. In addition, we cannot exclude stacking faults and domains in the epitaxial C60 adlayers, which could contribute to a larger effective diffusion barrier in our calculation as transport across domain boundaries is hindered. A more detailed comparison of our value for the diffusion barrier with values derived from pair potential calculations and molecular dynamics simulations is given in Supplementary Note 7. Without this strategy, the simulation of the full multilayer growth would have been impossible. Furthermore, the same strategy is used in simulations of metallic growth24,63,64 enabling a comparison with these studies.

In addition to the quantities discussed so far, KMC simulations allow us to extract single-particle trajectories and, thus, to study the dynamics on a particle level, which is not yet possible with current experimental techniques. An example of a single C60 particle trajectory (red) on top of a third monolayer island (light blue) is shown in Fig. 4a. Clearly, the Ehrlich–Schwoebel barrier leads to a ‘caging’ of the C60 molecule close to the borders of the island, that is, the standard random walk behaviour is restricted by the step edge of the island.

Figure 4: Particle-resolved dynamics during C60 growth.
figure 4

(a) Trajectory of a single molecule in the 4th ML (T=40 °C and f=1 ML min−1; scale bar, 5 nm). The influence of the Ehrlich–Schwoebel barrier can be clearly seen as a caging of the single C60 molecule on the island. The letters A and B denote the adsorption of one molecule on the surface (A) and the formation of a dimer (B). (b) MSD=‹|r(t)−r(0)|2› of C60 on C60(111), for T=60 °C and f=0.1 ML min−1 as a function of time spent on the surface. Results are averaged over 500 realizations. The particles considered arrive in the 2nd ML after the growth of 1.5 monolayers. For comparison, we show data for a system with atom-like ratio ED/ED+EB=0.34. Note that the quasi-free diffusion of C60 extends substantially further than for atom-like systems, even if scaled by the lattice parameter, signifying the qualitatively different behaviour of C60. (c) Schematic illustration of energy landscape for atoms, colloids and the fullerene C60 near an island step edge: The interaction range of the different materials clearly affects the character of step-edge barrier as one can distinguish between real and a diffusion-mediated pseudobarrier26.

Full size image

Importantly, the particle-resolved dynamics reveal crucial differences in the diffusion behaviour of C60 and atomic systems. For C60 on C60(111), the diffusion barrier ED is relatively large compared with the binding energy EB. Specifically, the ratio R=ED/(ED+EB) is R=0.83. This is significantly larger than in typical atomic systems, such as Pt on Pt(111) where R≈0.29–0.34, or Ag on Ag(111) with R≈0.29–0.39 (refs 23, 24). We suggest that this pronounced difference is related to the relatively short attractive interaction range of C60, as compared with the attraction range of atoms, if normalized to their respective size (for details see Supplementary Fig. 7 and Supplementary Notes 8 and 9). The comparatively large ratio R for C60 has a profound impact on the mobility of the particles. This is shown in Fig. 4b, where we plot the mean-squared displacement, MSD=‹|r(t)−r(0)|2› for particles arriving between islands after the growth of 1.5 monolayers for C60 and for a system with an atom-like ratio R=0.34.

The linear increase with time of the C60 MSD in the very beginning corresponds to free diffusion, depicted in grey, as the molecules perform a random walk on the underlying fcc(111) surface. After a time of about 0.1 ms, encountering an upward island edge as well as interactions with neighbours hinder the diffusion of the molecules, the MSD saturates. Similar sub-diffusive behaviour also occurs in the atom-like system, but at much shorter times. This is because atoms can form new bonds more quickly due to the longer range of atomic interactions and the stronger binding energy. As a result, a C60 molecule is able to explore an area that is nearly two orders of magnitude larger than in the atom-like system before it is immobilized.

The different diffusion behaviour of C60 prompts the question on the nature of the Ehrlich–Schwoebel barrier in comparison with atomic and colloidal growth. Indeed, regarding their narrow interaction range, C60 ‘nanocolloids’ are more similar to colloids than atoms. In colloids, the range of attractive interactions is so small that the reduced coordination associated with an edge is not ‘sensed’. This effectively leads to the vanishing of an energetic barrier at the edge. Instead, one observes a purely diffusive Ehrlich–Schwoebel barrier in colloids, arising from a lower diffusion probability along the geometrically longer path across the step edge26. In contrast, atoms crossing an island edge have to overcome an energetic Ehrlich–Schwoebel barrier, as bonds are missing at the step-edge. For C60, we can estimate an upper bound for a diffusive barrier based on the waiting time of a typical hopping process. Multiplying this time by a geometric factor (see ref. 26), which accounts for the longer path of a step-edge crossing, we obtain a diffusive pseudobarrier of EES,geo=ln(F) kBT<50 meV (see Supplementary Note 10 and Supplementary Fig. 8 for details). This is markedly smaller than the value of 110 meV obtained from the KMC simulations. We thus conclude that the Ehrlich–Schwoebel barrier in C60 surface growth is, at least partially, of energetic character, consistent with the intermediate range of the C60 interactions (which lies between the range of colloidal and atomic interactions). This is schematically shown in the energy landscapes for atoms, colloids and C60 in Fig. 4c.

In conclusion, the present experimental and theoretical study yields, for the first time, a quantitative description of molecular thin-film growth for the important case of C60, as an intermediate between atoms and colloids. We have demonstrated that in situ specular X-ray reflectivity and diffuse GISAXS oscillations are powerful tools for non-invasive real-time studies of the morphological evolution during molecular growth. Relating the experimental data to results from KMC simulations, we have been able to determine a consistent set of energy parameters determining the growth kinetics on the molecular level. This way we can quantitatively predict C60 deposition at different temperatures and rates, including the evolution of island density and surface roughening with film thickness. Thus, our combined analysis provides a detailed understanding of C60 in terms of molecular-scale processes. Moreover, our study sheds new light on various dynamical aspects accompanying the growth. In particular, we show that the colloid-like, short-ranged character of C60 interactions leads to relatively long surface diffusion times before immobilization occurs at existing islands. Nevertheless, the step-edge crossing barrier of C60 differs from colloids in that it is not a pseudo-step-edge barrier arising from lower diffusion probability at a step edge, but a true energetic barrier as observed for atoms. Since C60 features aspects of both atomic and colloidal systems, our findings will help to gain insight into island nucleation and surface growth processes for van der Waals-bound molecules between the scales of atomic and colloidal systems. This quantitative, scale-bridging understanding enables predictive simulations and a rational choice of growth conditions, which, together with molecular design and synthesis, ultimately leads to optimized design of functional materials.

Methods

X-ray surface scattering and thin-film preparation

The X-ray surface scattering experiments during growth were carried out at the MiNaXS beamline P03 (ref. 65) of PETRA III (DESY, Hamburg) at an X-ray wavelength of 0.946 Å. The growth was performed in a portable ulta-high vacuum (UHV) chamber designed for molecular beam deposition, equipped with a Be window for X-ray access, C60 effusion cell and a quartz crystal microbalance, at a base pressure of 10−8 mbar. Fullerene C60 (Sigma Aldrich, >99.5% purity) was thermally deposited on cleaved mica (diameter: 10 mm, Plano GmbH) for two different deposition rates (0.1 and 1 ML min−1) and for three different substrate temperatures (40, 60 and 80 °C) to study rate-, temperature- as well as time- and thickness dependency of the island density and layer coverage. Films were grown repeatedly on the same substrate after heating the mica substrate to ~450 °C, resulting in a clean substrate, as confirmed by specular and diffuse X-ray scattering before every growth run. The high brilliance of the beamline and high dynamic range of the PILATUS 300 K (Dectris) area detector enable a simultaneous measurement of the strong specular X-ray reflectivity and weak diffuse X-ray scattering. An incident angle of αi=1.65°, the so-called anti-Bragg position of C60 corresponding to half the Bragg value of the (111) reflection, was chosen. Here the reflectivity shows time-dependent oscillations during layer growth, which provide information on the vertical layer filling16,53. Lateral information is available through simultaneous measurement of the diffuse scattering (GISAXS), giving information about the island distance39, as a function of the lateral momentum of transfer q‖ at a resolution in q‖ of 0.001 Å−1. We avoided beam damage due to the high photon flux at PETRA III by laterally moving the substrate during the real-time growth experiments and confirmed that pristine and previously exposed spots gave the same scattering pattern in post-growth experiments.

Anti-Bragg intensity and sticking coefficient

The time-dependent anti-Bragg intensity can be calculated in kinematic approximation using

with the layer coverages θn for the nth layer. The substrate amplitude Asub, the substrate phase ϕsub and the molecular form factor f(qz) are determined by maximal, minimal and saturation intensity of the real-time experiment51. The anti-Bragg intensity for the KMC simulations was calculated using equation (2) and the simulated layer coverages shown in Fig. 4c. Furthermore, we have fitted the experimental data according to analytical growth models51,53 to extract the coverage evolution for each layer. In addition, we can extract the sticking coefficient from the anti-Bragg growth oscillations, which is found to decrease during the growth of the first four layers for all studied temperatures. Quantitatively, we find for a temperature of 60 °C that with respect to the growth of the first monolayer, the sticking coefficient decreases by 5% in the 2nd ML, 25% in the 3rd ML and 30% from the 4th layer onwards. This decrease is due to the different mica–C60 and C60–C60 interactions. It is further influenced by a different island density in each layer, which leads to a change in the free diffusion times and aggregation behaviour. In our KMC simulations, which otherwise assume complete condensation, we have accounted for the changing sticking coefficient by scaling the molecular exposure axis accordingly. The same sticking coefficients have also been included in our analytical mean-field modelling.

Time step in KMC simulations

Assuming that exactly one process takes place in one simulated time step, we can define an average time-step length as

This time unit allows us to compare simulated and experimental timescales. The simulation is carried out on a triangular lattice. In this way, the growth process generates a fcc structure in accordance with the C60 bulk crystal (see the studies of Cox et al.22 for a similar simulation strategy for the growth of Ag on Ag(111), and of Heinrichs et al.66 for corresponding theoretical considerations). Starting point of the simulation is a completely filled, defect-free layer of C60 molecules (corresponding to the C60(111) surface). Within the subsequent growth process, we exclude the formation of overhangs. To achieve this, we assume that particles on overhang sites relax instantaneously (with a relaxation probability proportional to the corresponding diffusion rate) until they reach a stable site. Typical simulations involve a lattice with 1,000 × 1,000 unit cells, and they cover a time range up to 4,000 s, corresponding to (1011–1012) events.

Additional information

How to cite this article: Bommel, S. et al. Unravelling the multilayer growth of the fullerene C60 in real-time. Nat. Commun. 5:5388 doi: 10.1038/ncomms6388 (2014).