We have studied nucleation and growth of Dy islands on the basal plane of graphite at 300 K using scanning tunneling microscopy, density functional theory (DFT) in a form that includes van der Waals interactions, and analytic theory. The interaction of atomic Dy with graphite is strong, while the diffusion barrier is small. Experiment shows that at 300 K, the density of nucleated islands is close to the value predicted for homogeneous nucleation, using critical nucleus size of 1 and the DFT-derived diffusion barrier. Homogeneous nucleation is also supported by the monomodal shape of the island size distributions. Comparison with the published island density of Dy on graphene shows that the value is about two orders of magnitude smaller on graphite, which can be attributed to more effective charge screening in graphite. The base of each island is 3 atomic layers high and atomically ordered, forming a coincidence lattice with the graphite. Islands resist coalescence, probably due to multiple rotational orientations associated with the coincidence lattice. Upper levels grow as discernible single-atom layers. Analysis of the level populations reveals significant downward interlayer transport, which facilitates growth of the base. This island shape is metastable, since more compact three-dimensional islands form at elevated growth temperature.

The interaction of magnetic materials with graphitic surfaces, from single-layer graphene to few-layer graphene to bulk graphite, holds possibilities for designing and controlling magnetic properties at the nanoscale. For graphene, in particular, it has been proposed that the carbon sheet may play a valuable role in new spintronics devices, or as a support or modifier of magnetic materials.1–5 In most cases, it is necessary to achieve a good contact between the magnetic material and the carbonaceous substrate, in the form of a thin, essentially two-dimensional (2D) metal layer.1 However, this has proven difficult since the adsorption energy of a metal on a graphitic surface is often much less than its cohesive energy, a factor which promotes three-dimensional (3D) growth of metals.6 

In terms of magnetic properties, rare earths such as Dy attract special interest because of their high magnetic moment. This paper is a study of Dy islands that form when Dy is deposited on bulk graphite at room temperature. It sheds light on the fundamental energetics and mechanisms involved in adsorption, diffusion, nucleation, and growth, which may contribute to an ability to predict and manipulate Dy nanostructures, and the conditions under which they form.

There are two prior sets of studies of Dy’s interaction with related substrates. In one set, the substrate was graphene grown on SiC(0001).7–9 There, it was shown that deposition of Dy at 660 K produced compact three-dimensional (3D) islands, with (mainly) triangular shapes indicative of fcc atomic structure, rather than the bulk hcp structure.7 Smaller, more irregular islands were formed at 300 K and were very stable against coarsening when annealed at 600 K.8,9 In the second set, the substrate was highly oriented pyrolytic graphite (HOPG).10 A 20 nm film of Dy was deposited, corresponding to about 70 atomic layers if distributed uniformly. Annealing to 1200 K produced a carbide-like surface species, consistent with the high propensity for rare earths to form carbides.11 

The experiments were performed in an ultrahigh vacuum (UHV) chamber with base pressure 5 × 10−11 mbar. The graphite samples were HOPG of grade ZYA in experiments at low coverage, i.e., below 0.8 monolayer (ML) of Dy, or ZYB in experiments at higher coverage. The clean graphite surface was prepared by tape-cleavage in air, followed by transfer into the UHV chamber where all subsequent experiments took place. Samples were heated in the manipulator to either 800 K for 60 min (low coverage experiments) or 1000 K for 20 min (high coverage experiments) in UHV to remove contaminants and then transferred to the scanning tunneling microscopy (STM) stage for subsequent Dy deposition and STM imaging. Dy was deposited via physical vapor deposition from a Mantis QUAD-EV-C Mini e-beam evaporator, using a Mo crucible lined with pyrolytic boron nitride. For heating, the crucible was biased at +2 kV with respect to an electron filament mounted parallel to and near the top of the evaporation target, with total power typically 25 W.

In a study of Cu/graphite, we showed that this evaporator, under similar conditions, produced a significant fraction of metal ions that damaged the graphite surface and influenced nucleation. We also showed that this could be circumvented by heating the evaporator, then shutting off the high voltage and/or filament current during deposition, with only a small drop in flux over times up to 10 s.12 For all experiments reported herein, we employed this protocol. The deposition time was constant at 10 s; coverage was varied by adjusting power and hence flux.

In STM, tunneling parameters were typically in the range −1.0 V to −1.2 V tip bias and 0.15 nA–0.25 nA tunneling current for low coverage experiments, and 0.70 V–1.3 V and 0.13 nA–0.77 nA for higher coverage experiments. The STM tip was electrochemically etched W. All STM images were planed for data analysis. Other details of STM experiments and data analysis were the same as reported elsewhere.12 

First-principles energy calculations with non-local van der Waals correction were performed based on density functional theory (DFT) using VASP.13,14 The exchange and correlation energy functional adapted the opt-B88 scheme developed by Klimeš et al.15–17 This functional has been verified to describe accurately the energy and other properties of graphite and metals.15,17 The electron-ion interaction was described by the projector augmented wave method.18 The energy cutoff for the plane wave basis set used in the calculation was 600 eV. The calculated fundamental properties of graphite agreed well with the corresponding experimental values.19 The calculated lattice constants of hexagonally close-packed Dy were 0.357 nm for a and 0.558 nm for c, respectively, in good agreement with experimental data of 0.359 nm for a and 0.565 nm for c.20 The calculated cohesive energy was 3.58 eV/atom, somewhat larger than the experimental value, 3.04 eV/atom.20 

The graphite substrate was modeled by a slab with a 6 × 6 unit cell in the xy plane and four layers along the (0001) direction, plus enough vacuum (1.57 nm) to avoid interaction between the slab and its images under the periodic boundary condition. Spin polarization and dipole correction were considered in all calculations. A Γ-centered k-point grid of 15 × 15 × 1 was used for Brillouin zone sampling to ensure energetic convergence. During geometric optimization, the bottom three layers were fixed at their bulk positions, while carbon atoms in the top layer and metal adatom were relaxed fully with a force tolerance of 0.1 eV/nm.

Table I shows calculated values of the adsorption energy of a Dy atom at different sites on the graphite surface. The adsorption energy at the favored site is 1.90 eV. This is higher than the value of 1.47 eV calculated for Dy on graphene,21 which reinforces a trend noted elsewhere: For a given metal, its calculated adsorption energy on graphite exceeds that on graphene.19 The adsorption energy is about 2/3 the cohesive energy of bulk Dy, 3.04 eV. The adsorption energy falls in the range of chemisorption and implies an essentially infinite surface lifetime for a Dy atom at 300 K. By all measures, the interaction of Dy with graphite is reasonably strong.

TABLE I.

DFT results for a Dy atom adsorbed on graphite. A schematic of atoms in the (0001) graphite plane is shown in top left, where large circles represent top-layer C atoms and smaller circles represent second-layer C atoms. Four high-symmetry adsorption sites are identified: two top (T) sites, which differ in occupation of the second-layer site; hollow (H); and bridge (B).

TαTβBridge (B)Hollow (H)
Adsorption energy (eV) −1.85 −1.83 −1.83 −1.90 
Height above carbon atom plane (nm) 0.241 0.242 0.239 0.214 
TαTβBridge (B)Hollow (H)
Adsorption energy (eV) −1.85 −1.83 −1.83 −1.90 
Height above carbon atom plane (nm) 0.241 0.242 0.239 0.214 

The Dy atom is most stable at the hollow (H) site, in the center of the hexagon of carbon atoms. (See inset in Table I for identification of adsorption sites.) This is the same as the site predicted for Dy and other rare earths on free-standing graphene.9 The diffusion barrier, Ed, can be approximated as the difference in adsorption energies along a path that links minimum-energy sites. In this case that path leads from the H site, over Tα and back to H. The energy difference between the H and Tα sites is 0.048 eV.

Fig. 1 shows representative STM images of Dy islands on graphite terraces, for coverages spanning 0.12–2.5 ML. The islands consistently have a flat base with an irregular footprint. The flat base is decorated with smaller upper features, many of which are found near the island edges. The base is 0.87 ± 0.02 nm high, based on a sample size, N, of 36 islands. This height is invariant with island size. To first order, it can be compared with the 0.282 nm spacing between close-packed planes in bulk Dy, which yields a thickness of 3 atomic layers. Thus, the base can be described as quasi two-dimensional (2D), which is not unreasonable in light of the relatively high adsorption energy found in Sec. III.

FIG. 1.

STM images of Dy islands on HOPG terraces, 100 nm × 100 nm. Coverages are: (a) 0.12 ML; (b) 0.28 ML; (c) 1.2 ML; (d) 2.5 ML. Noise is reduced in (d) using a fast Fourier transform (FFT) filter.

FIG. 1.

STM images of Dy islands on HOPG terraces, 100 nm × 100 nm. Coverages are: (a) 0.12 ML; (b) 0.28 ML; (c) 1.2 ML; (d) 2.5 ML. Noise is reduced in (d) using a fast Fourier transform (FFT) filter.

Close modal

In many cases a honeycomb pattern with amplitude 0.04 ± 0.01 nm (N = 144 islands) can be observed on the base, as shown in Fig. 2. This is a moiré pattern between the Dy film and the graphite support, which indicates that the Dy film is atomically well-ordered. The modulation periodicity is 1.23 nm. In moiré patterns such as this, a rotation between the unit cell vectors of the overlayer and the substrate, along with possible distortion of the overlayer, brings the two lattices into coincidence. Unfortunately, atomic-scale resolution of the graphite support could not be obtained in the presence of Dy islands, so the rotational angle is unknown. However, the existence of a coincidence lattice provides further evidence of significant interaction between Dy and the carbon substrate.

FIG. 2.

STM images of Dy islands, 20 nm × 20 nm, at (a) 0.28 ML and (b) 0.12 ML. The enlargement shows the honeycomb pattern visible on the base of (a), 7.5 nm × 3.2 nm. The enlargement is shown after FFT filtering.

FIG. 2.

STM images of Dy islands, 20 nm × 20 nm, at (a) 0.28 ML and (b) 0.12 ML. The enlargement shows the honeycomb pattern visible on the base of (a), 7.5 nm × 3.2 nm. The enlargement is shown after FFT filtering.

Close modal

At high Dy coverage, Dy islands become very close but often retain their individual identity, as illustrated in Fig. 3. For example, the oval in Fig. 3(a) encompasses 3 islands separated by thin gaps, indicating a barrier to coalescence. We attribute this to the different epitaxial orientations of the Dy base, associated with the moiré pattern. Only when orientations happen to match can adjacent islands merge easily.

FIG. 3.

STM images of Dy islands at 2.5 ML. Image sizes are 250 nm × 250 nm. The white oval in (a) shows three islands that are very close but still separated.

FIG. 3.

STM images of Dy islands at 2.5 ML. Image sizes are 250 nm × 250 nm. The white oval in (a) shows three islands that are very close but still separated.

Close modal

The bases often have taller dot-like features around their edges. The height of these dots corresponds to 1 Dy atomic layer above the base. Examples are visible in Figs. 1 and 4. Larger, layer-like features can also cover part of the islands’ interiors. These are 1.13 ± 0.02 and 1.42 ± 0.02 nm above the graphite. Figure 4 provides examples of these features, including line profiles that illustrate heights. The fact that each layer is higher by 0.28 nm indicates that single Dy layers are populated successively, above the base. With increasing Dy coverage, it becomes more probable to find these partially filled layers. The layers can exist near the center of an island, or emanate from an edge. When more than one layer exists atop the base, layers are stacked in a wedding-cake morphology such that each level below the top one is partially exposed.

FIG. 4.

STM images of Dy islands on graphite terraces, with line profiles. Each image has been differentiated to facilitate viewing of small features on different levels. Image size is 60 nm × 60 nm. (a) 0.14 ML ((b)-(d)) 1.2 ML.

FIG. 4.

STM images of Dy islands on graphite terraces, with line profiles. Each image has been differentiated to facilitate viewing of small features on different levels. Image size is 60 nm × 60 nm. (a) 0.14 ML ((b)-(d)) 1.2 ML.

Close modal

For comparison, Fig. 5 shows the result when Dy is deposited at 800 K. Now the Dy islands are much taller and more three-dimensional; the three-layer mesa-like base is no longer evident. This indicates that the quasi-2D shape at 300 K is metastable. Indeed, the equilibrated Wulff shape of crystalline islands should be relatively tall compared to their width for a metal weakly bound to the substrate. The footprint should also be geometric, e.g., hexagonal for hcp crystal structure. Presumably, upward diffusion beyond the first three layers becomes facile at the higher temperature, allowing islands to achieve a more equilibrated shape during growth, than that which can be attained at 300 K.

FIG. 5.

Topographic STM image of Dy islands on a graphite terrace after deposition at 800 K, with associated line profiles. Image size is 250 nm × 250 nm, and tunneling conditions are 1.3 V, 0.15 nA.

FIG. 5.

Topographic STM image of Dy islands on a graphite terrace after deposition at 800 K, with associated line profiles. Image size is 250 nm × 250 nm, and tunneling conditions are 1.3 V, 0.15 nA.

Close modal

The island density, Nisl, is shown as a function of Dy coverage in Fig. 6. The density increases from a value of (2.4 ± 0.3) × 10−4 nm−2 at 0.12-0.15 ML and saturates at 6.0 × 10−4 nm−2 at 1.2 ML.

FIG. 6.

Dy island density as a function of Dy coverage. The curve represents the point island model described in the text, with Nisl (in units of nm−2) = 9.633 × 10−4θ0.676. Error bars reflect the standard error of the mean, calculated from the set of images available within each experiment, and do not reflect variations between different experiments.

FIG. 6.

Dy island density as a function of Dy coverage. The curve represents the point island model described in the text, with Nisl (in units of nm−2) = 9.633 × 10−4θ0.676. Error bars reflect the standard error of the mean, calculated from the set of images available within each experiment, and do not reflect variations between different experiments.

Close modal

The island size distribution is shown in Fig. 7 for two coverages, 0.13 and 1.2 ML. In both cases the distribution is monomodal.

FIG. 7.

(a) 0.13 ML sav = 4827 with bin size = 0.21. The bin size is chosen to correspond to 1000 atoms/bin. The histogram represents 287 islands. (b) 1.2 ML sav = 20 434 with bin size = 0.20. The bin size is chosen to correspond to 4000 atoms/bin. The histogram represents 270 islands.

FIG. 7.

(a) 0.13 ML sav = 4827 with bin size = 0.21. The bin size is chosen to correspond to 1000 atoms/bin. The histogram represents 287 islands. (b) 1.2 ML sav = 20 434 with bin size = 0.20. The bin size is chosen to correspond to 4000 atoms/bin. The histogram represents 270 islands.

Close modal

These observations are analyzed in Sec. V.

First we assess the island nucleation mechanism by comparing the observed value of Nisl and its coverage dependence for deposition at 300 K with theories for homogeneous nucleation. It is useful to note that the area of the surface unit cell for graphite is Ω = 0.052 nm2. Multiple experiments for deposition of 0.12-0.15 ML provide quite consistent estimates for Nisl (see Fig. 6). Averaging these data yields Nisl = (2.4 ± 0.3) × 10−4 nm−2 (or 1.3 × 10−5 per unit cell or per adsorption site) for coverage θ = 0.13 ± 0.1 ML corresponding to a deposition flux of F = 0.013 ML/s.

To analyze this behavior, we use the estimate for the terrace diffusion barrier of Ed ≈ 48 meV, obtained from DFT analysis described in Sec. III. Assuming an attempt frequency of ν = 1012.5 s−1 yields a rate of h = νexp[ − Ed/(kBT)] for adatom hopping associated with terrace diffusion. A key parameter for nucleation is h/F. Using the experimental values of T = 300 K and F = 0.013 ML/s, one has h = 1011.7 s−1 and h/F = 1013.6 s−1. For homogeneous nucleation and growth of islands with critical size i (where islands of size >i atoms are stable), the island density satisfies Nislθ(1−χ)/2(h/F)χ where χ = i/(i + 2) in the regime of low coverage.22 

For this system, it is reasonable to anticipate that homogeneous nucleation and growth of islands is at least close to irreversible, corresponding to critical size i = 1. Theoretical analysis indicates that the temperature for transition to reversible island formation is proportional to Erev = (2/3)Ed + ENN, where ENN > 0 is the strength of the nearest-neighbor interaction.22 More precisely, to ensure a transition temperature above 300 K for a standard deposition flux of 0.005 ML/s and attempt frequency of ν = 5 × 1012 s−1, one requires that Erev ≈ 0.55 eV or above. Thus, for Dy on HOPG with Ed ≈ 0.05 eV, one requires that ENN ≈ 0.52 eV or above. A crude estimate of ENN as 1/6 of the bulk cohesive energy yields a value of 0.51 eV, and it is plausible that the true ENN for a pair of Dy on HOPG is stronger.

For i = 1, one has that Nisl ∼ (h/F)−1/3 for fixed coverage, and simulations indicate that Nisl ≈ 2.18 × 10−5 per site for θ = 0.13 ML.22,23 Thus, the experimental value of Nisl, 1.26 × 10−5 per site, is close to but slightly reduced from classic i = 1 prediction (reduced by a factor of r = 0.58).

The lower experimental value could be achieved by choosing an attempt frequency which is higher by a factor of 5.1. Such a prefactor, just above 1013 s−1, cannot be excluded, but one generally expects low prefactors for low diffusion barriers (as in this system). Alternatively, the lower experimental Nisl could be naturally explained by deviations from classic i = 1 behavior due to either significant mobility of adsorbed Dy dimers (and perhaps other small Dy clusters) or due to the onset of reversibility in island formation. To obtain a sense of the impact of dimer mobility or reversibility in reducing Nisl, we note the following. If monomers and dimers have comparable mobility, then Nisl ∼ (h/F)−2/522,24 which implies a reduction in Nisl from the classic i = 1 value by a factor of r = (h/F)−1/15 = 0.12 for the above h/F = 1013.6, i.e., the reduction is too strong. For facile dissociation of dimers but stable trimers corresponding to i = 2, one has that Nisl ∼ (h/F)−1/2, reducing Nisl for the classic i = 1 value by a factor of r = (h/F)−1/6 = 0.005, which is much too strong. In conclusion, the observed reduction by r = 0.58 relative to classic i = 1 could be induced by some dimer mobility24 (still significantly lower than adatom mobility) or by slight reversibility in island formation.

Another way of accounting for these effects is by introducing a scaling exponent χ = 1/3 + δχ which is slightly above the value of 1/3 for i = 1, so that (h/F)δχ = r = 0.58. This implies that δχ = 0.018 so χ = 0.351 and (1 − χ)/2 = 0.325. As a result, Nisl = c(h/F)−0.351θ0.325. Choosing c to recover the experimental Nisl at 0.13 ML and recognizing that F = θ/10 in experiment, one obtains Nisl ≈ 5.01 × 10−4θ0.676 per site, which should apply in the lower coverage regime where the footprint of the islands covers a small fraction of the substrate area (i.e., the “point island” regime22). This curve is shown in Fig. 6 up to a coverage of 0.4 ML, beyond which one expects that the point island model will significantly overestimate Nisl because it underestimates the capture probability at existing islands. (We emphasize that for the Dy islands with a 3-layer base, the areal coverage or fraction of the substrate covered is only ∼θ/3 in contrast to θ for 2D islands. Thus, nucleation persists to higher coverages in the 3-layer case.) There is a good agreement between the model and the data in this regime.

Next, we analyze experimental observations for the island size distribution, Ns, which gives the density of islands of s atoms, so that Nisl=s>1Ns and θ=ssNs. If savθ/Nisl denotes the average island size (measured in atoms), then the island size distribution is naturally written in scaled form, Ns ≈ (Nisl/sav) f(s/sav, θ), where the scaling or shape function f(x) satisfies ∫x>0 dxf(x) = 1 and ∫x>0 dxxf(x) = 1.22 For homogeneous nucleation at lower coverages around 0.1 ML, f(x) is monomodal with a peak height around x = 1 increasing from about 0.75 for i = 1 to 1.1, 1.3, and 1.5 for i = 2, 3, and 6, respectively (and also f(x) narrows with increasing i). The peak also increases to about 1.1 with high dimer mobility.24 There is a significant population of smaller islands reflected in f(0) ≈ 0.3 for i = 1, with f(0) quickly decreasing with increasing i (e.g., to 0.15 for i = 2). The experimental f(x) shown in Fig. 7(a) for 0.13 ML has a peak height around 1.0 consistent with slight deviations from classic i = 1 behavior. However, there is a significant depletion in the population of small islands relative to classic i = 1 behavior. This feature together with an increase in peak height and shift to smaller x-values has been seen previously and attributed to post-deposition diffusion and coalescence of small clusters.25 Fig. 7(b) shows the experimental f(x) for a much higher coverage of 1.2 ML. The peak height is similar to 0.13 ML, but the peak and the weight of the distribution have shifted to higher (scaled) island sizes. This feature might be expected, as inspection of STM images reveals that island growth is impeded by impingement upon other islands (without coalescence) in this regime. As noted in Sec. IV, we attribute this to epitaxial mismatch between different islands. This would tend to induce a sharper cutoff on the right side of the distribution.

Despite the complexity of multilayer growth in this system (described in more detail below), simplified modeling can produce significant insights, e.g., into the extent of interlayer transport. To provide some context, such models have been traditionally developed for systems with equivalent layers where growth in higher layers proceeds before growth in lower layers is complete. If θj denotes the coverage of layer j, then Pj = θjθj+1 denotes the exposed fraction of layer j. In the extreme case where interlayer transport is absent, one finds a Poisson distribution for Pj and film roughness W satisfies W2 = d2θ, where d is the interlayer spacing.22 This modeling has been refined to account for downward interlayer transport by specifying that a fraction, α, of atoms deposited in each layer reach the next lower layer (so α = 0 recovers the above case).26–28 Then, for rough growth with small α, one can show that W2 ≈ (1 − 2α) θd2.27,28 Although not well-recognized, this result must be modified for larger α [see  Appendix A].

A special feature of the Dy/graphite system is that atoms deposited on the substrate form bases consisting of three atomic layers (rather than single layer islands). This implies facile transport from the substrate upward to the second and third layers. However, on top of these bases, growth occurs in single layers. Thus, it is convenient to refer to the fraction, Pj, of exposed level, j (rather than layer), where j = 0 is the substrate, j = 1 is the top of the 3-layer base, j = 2 is the top of the single layer island on the base, etc. Experimental values for these populations or fractions of exposed levels for coverage θ = 1.2 ML are given in Table II.

TABLE II.

Layer distributions in Dy islands at 1.2 ML. Pj is the fraction of area exposed for each level j, where j = 0 indicates the graphite substrate, j = 1 is the top of the Dy base, j = 2 is the first Dy layer atop the base, and j = 3 is the second Dy layer above the base. The parameters α and α are defined in Fig. 8. The experimental values were obtained by analyzing a total area of 1.38 × 106 nm2, which included 792 Dy islands.

Model
Experimentα = α = 0α = α = 0.65α = 0.69 α = 0.85
P0 0.631 0.670 0.634 0.631 
P1 0.276 0.185 0.278 0.281 
P2 0.0904 0.096 0.072 0.081 
P3 0.0024 0.036 0.015 0.006 
Model
Experimentα = α = 0α = α = 0.65α = 0.69 α = 0.85
P0 0.631 0.670 0.634 0.631 
P1 0.276 0.185 0.278 0.281 
P2 0.0904 0.096 0.072 0.081 
P3 0.0024 0.036 0.015 0.006 

To elucidate this behavior, we develop a theory in the spirit of the more traditional one mentioned above, for multilayer growth with equivalent layers. However, here we naturally specify that all atoms deposited on the graphite substrate are incorporated immediately into 3-layer bases. A fraction α of atoms deposited on top of bases hop down to become incorporated in the base and the remainder stay on top of the base forming single-layer islands. Of atoms deposited in higher layers, a fraction α hop down to be incorporated in the next lower layer, and the remainder stay in the layer where they were deposited, forming or joining single-layer islands in that layer. See the schematic in Fig. 8. With this prescription, the rate equations describing the evolution of the level populations are

Simple modification of these equations is possible if one wishes to exclude population of level 4 (i.e., to enforce P4 = 0) consistent with experiment.

FIG. 8.

Schematic of model for multilayer growth.

FIG. 8.

Schematic of model for multilayer growth.

Close modal

It is natural to first consider the extreme case of no interlayer transport (other than upward transport of atoms deposited on the substrate to form 3-layer bases). In this case α = α = 0, and the above equations partially decouple in the sense that Pj is coupled only to Pk<j so they can be recursively solved exactly to obtain

As shown in Table II, the resultant values of Pj agree roughly with experimental values, but the fact that P0 is higher in the model than in experiment (0.670 vs. 0.631, respectively) is a clear indication that net downward 3 transport with α > 0 is operative in the Dy/graphite system.

To refine the above simplest model, we introduce interlayer transport, initially setting α = α. We also enforce P4 = 0. Numerical integration of the associated rate equations indicates that good agreement with experiment comes from choosing α = 0.65, with the resultant values of Pj shown in Table II. A slightly better fit is achieved by allowing independent values α = 0.69 and α = 0.85. In either case, this analysis indicates that there is significant net downward transport in the Dy/graphite system, which impacts the film height distribution.

The existence of small single-layer features at the edges of islands was noted in Sec. IV. Their presence indicates that nucleation is somewhat more likely near islands’ edges than in islands’ interiors. This could divert some Dy atoms from reaching and crossing over the edge of the base, leading to α < α as suggested by the modeling.

One major conclusion from this work is that Dy island nucleation is homogeneous. This is based upon the measured value of Nisl at 0.13 ML, in comparison with the value predicted from a point-island growth model using the diffusion barrier from DFT and i = 1. The agreement is quite good, with the experimental value being slightly lower than the predicted value. Elsewhere, we have argued that homogeneous nucleation of metals on graphite terraces may be less common than previously thought, with heterogeneous nucleation playing a larger role than expected due to the susceptibility of graphite substrates to inadvertent ion beam damage.12 We have argued that if the measured island density is significantly higher than the predicted value, then heterogeneous nucleation should be considered. If the opposite is true, as in this case, then homogeneous nucleation is viable with some adjustments. Here, we have shown that the degree of adjustment is too extreme if i = 2 or if dimer diffusion is fully operative, though a small contribution from one or both of these mechanisms could explain the discrepancy.

The island size distributions are also consistent with homogeneous nucleation. Their monomodal shapes are quite different from the monotonically decreasing shape found in our previous study of Cu on graphite, where Cu ions in the metal atom flux introduced defects in the graphite surface during deposition. Precautions are taken in the present work to prevent such damage, consistent with the radically different distributions observed.

As noted in Sec. I, Hupalo et al. studied a related system, Dy deposited on monolayer graphene supported on SiC, at 300 K.8 They reported Nisl = 2 × 10−2 nm−2 at 0.4 ML, which is two orders of magnitude higher than our value of Nisl at 0.13 ML. While a somewhat higher value of Ed = 0.125 eV was calculated for that system, this and other differences—in experimental flux and coverage—are insufficient (by far) to explain the difference in Nisl. Instead, we propose that the explanation lies in the influence of the number of carbon layers, n, on nucleation, growth, and stability of metal nanoclusters. There is now broad evidence for this dependence, with Nisl generally decreasing as n increases.6 This has been demonstrated in studies of Au, Ag, and Pd nanoclusters on n-layer graphene sheets supported on SiO2/Si substrates. This trend, decreasing Nisl with increasing n, describes exactly the trend in going from Dy/graphene (n = 1) in the prior work, to Dy/graphite (n → ∞) in this study. One explanation involves charge transfer from the metal to the underlying carbon and consequent Coulomb and dipole repulsion terms, which have significant effect at small n due inefficient screening.29 

Finally, we comment on the morphology of the Dy islands. Experiment indicates that growth is quasi-2D. The islands exhibit a 3-layer base, and single atomic layers grow on top of the base. DFT shows that the interaction between Dy and the graphite is reasonably strong, providing a rationale for the quasi-2D base. Analysis of the level populations reveals significant downward interlayer transport, which facilitates growth of the base. On the other hand, the quasi-2D structure is metastable, based on comparison with the taller, more compact 3D islands that form at 800 K; the most natural explanation is that upward transport (beyond the 3-layer base) is kinetically limited at 300 K.

It is possible that the base is limited to 3 atomic layers because of strain between Dy and graphite lattices. It is well-established that strain influences growth in heteroepitaxial systems. For the growth mode known as Stranski-Krastanov (SK), the film wets the substrate (growth is smooth) up to a critical thickness, beyond which growth is 3D.30,31 In SK systems strain is often manifest in a moiré pattern at or below the critical thickness.32–35 Strictly speaking, SK is a thermodynamic (equilibrium) picture of growth, but the same factors that lead to SK growth can influence structures formed under kinetically limited conditions.

We have shown that experimental measurements of the value of Dy island density, and coverage-dependence of Dy island density at low coverage, are well described by analytic theory for homogeneous nucleation, using energetics provided by DFT. The Dy islands have a quasi-2D 3-layer base which exhibits atomic order coincident with the graphite substrate. This natural tendency toward quasi-2D growth may be useful in situations where good contact between a magnetic metal and a graphitic substrate is desired.

The experimental and DFT components of this work (Secs. III and IV) were conducted by the following authors: E.J.K., H.L., A.L.-R., M.W., Y.Z., C.-Z.W., M.C.T., and P.A.T. Experimental and DFT effort was supported by the US Department of Energy (DOE), Office of Science, Basic Energy Sciences, Materials Sciences and Engineering Division. Research was performed in part at the Ames Laboratory, which is operated by Iowa State University under Contract No. DE-AC02-07CH11358. DFT was performed, in part, with a grant of computer time at the National Energy Research Scientific Computing Centre (NERSC). NRSEC is a DOE Office of Science User Facility supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231. H.L.’s DFT was supported in part by the National Science Foundation of China (NSFC) Grant No. 11575230, and Y.Z.’s participation in the experimental effort was supported by the China Scholarship Council. The modeling and analyses described in Sec. V were performed by J.W.E., with support from NSF Grant No. CHE-1111500.

For equivalent layers with deposition at rate F so θ = Ft and a fraction α of deposited atoms reaching the next lower layer, one has that

Separate equations are needed for j = 0 and j = 1. Neglecting these different equations yields the simple result for W2 = (1 − 2α) θd2 quoted in the text. However, reliable results for smoother growth with larger α must account for the modified equations, noting that clearly this expression for W2 cannot apply for α > 1/2.

Data used to generate this manuscript are available via the DOI: http://dx.doi.org/10.17039/ameslab.dmse.2016.DS2/1240605.

1.
C.
Vo-Van
,
Z.
Kassir-Bodon
,
H.
Yang
,
J.
Coraux
,
J.
Vogel
,
S.
Pizzini
,
P.
Bayle-Guillemaud
,
M.
Chshiev
,
L.
Ranno
,
V.
Guisset
,
P.
David
,
V.
Salvador
, and
O.
Fruchart
, “
Ultrathin epitaxial cobalt films on graphene for spintronic investigations and applications
,”
New J. Phys.
12
,
103040
(
2010
).
2.
B.
Uchoa
,
T. G.
Rappoport
, and
A. H.
Castro Neto
, “
Kondo quantum criticality of magnetic adatoms in graphene
,”
Phys. Rev. Lett.
106
,
016801
(
2011
).
3.
H.
Sevinçli
,
M.
Topsakal
,
E.
Durgun
, and
S.
Ciraci
, “
Electronic and magnetic properties of 3d transition-metal atom adsorbed graphene and graphene nanoribbons
,”
Phys. Rev. B
77
,
195434
(
2008
).
4.
S. M.
Binz
,
M.
Hupalo
,
X.
Liu
,
C. Z.
Wang
,
W.-C.
Lu
,
P. A.
Thiel
,
K. M.
Ho
,
E. H.
Conrad
, and
M. C.
Tringides
, “
High island densities and long range repulsive interactions: Fe on epitaxial graphene
,”
Phys. Rev. Lett.
109
,
026103
(
2012
).
5.
A. N.
Rudenko
,
F. J.
Keil
,
M. I.
Katsnelson
, and
A. I.
Lichtenstein
, “
Adsorption of cobalt on graphene: Electron correlation effects from a quantum chemical perspective
,”
Phys. Rev. B
86
,
075422
(
2012
).
6.
X.
Liu
,
Y.
Han
,
J. W.
Evans
,
A. K.
Engstfeld
,
R. J.
Behm
,
M. C.
Tringides
,
M.
Hupalo
,
H.-Q.
Lin
,
L.
Huang
,
K.-M.
Ho
,
D.
Appy
,
P. A.
Thiel
, and
C.-Z.
Wang
, “
Growth morphology and properties of metals on graphene
,”
Prog. Surf. Sci.
90
,
397
(
2015
).
7.
M. T.
Hershberger
,
M.
Hupalo
,
P. A.
Thiel
, and
M. C.
Tringides
, “
Growth of fcc(111) Dy multi-height islands on 6H-SiC(0001) graphene
,”
J. Phys.: Condens. Matter
25
,
225005
(
2013
).
8.
M.
Hupalo
,
X.
Liu
,
C.-Z.
Wang
,
W.-C.
Lu
,
Y.-X.
Yao
,
K.-M.
Ho
, and
M. C.
Tringides
, “
Metal nanostructure formation on graphene: Weak versus strong bonding
,”
Adv. Mater.
23
,
2082
(
2011
).
9.
X.
Liu
,
C. Z.
Wang
,
M.
Hupalo
,
W. C.
Lu
,
M. C.
Tringides
,
Y. X.
Yao
, and
K. M.
Ho
, “
Metals on graphene: Correlation between adatom adsorption behavior and growth morphology
,”
Phys. Chem. Chem. Phys.
14
,
9157
(
2012
).
10.
S. A.
Gorovikov
,
A. M.
Shikin
,
G. V.
Prudnikova
,
V. K.
Adamchuk
,
S. L.
Molodtsov
,
C.
Laubschat
, and
A. M.
Ionov
, “
Formation of surface intercalation compounds at Gd(Dy)/graphite interfaces under thermal annealing
,”
Surf. Sci.
474
,
98
(
2001
).
11.
G.-Y.
Adachi
,
N.
Imanaka
, and
F.
Zhang
, “
Rare earth carbides
,”
Handbook on the Physics and Chemistry of Rare Earths
15, Chap. 99, 61–189 (1991).
12.
D.
Appy
,
H.
Lei
,
Y.
Han
,
C.-Z.
Wang
,
M. C.
Tringides
,
D.
Shao
,
E. J.
Kwolek
,
J. W.
Evans
, and
P. A.
Thiel
, “
Determining whether metals nucleate homogeneously on graphite: A case study with copper
,”
Phys. Rev. B
90
,
195406
(
2014
).
13.
G.
Kresse
and
J.
Furthmüller
, “
Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set
,”
Comput. Mater. Sci.
6
,
15
(
1996
).
14.
J. P.
Perdew
,
K.
Burke
, and
M.
Ernzerhof
, “
Generalized gradient approximation made simple
,”
Phys. Rev. Lett.
77
,
3865
(
1996
).
15.
J.
Klimeš
,
D. R.
Bowler
, and
A.
Michaelides
, “
Van der Waals density functionals applied to solids
,”
Phys. Rev. B
83
,
195131
(
2011
).
16.
J.
Klimeš
,
D. R.
Bowler
, and
A.
Michaelides
, “
Chemical accuracy for the van der Waals density functional
,”
J. Phys.: Condens. Matter
22
,
022201
(
2010
).
17.
G.
Graziano
,
J.
Klimeš
,
F.
Fernandez-Alonso
, and
A.
Michaelides
, “
Improved description of soft layered materials with van der Waals density functional theory
,”
J. Phys.: Condens. Matter
24
,
424216
(
2012
).
18.
P. E.
Blöchl
, “
Projector augmented-wave method
,”
Phys. Rev. B
50
,
17953
(
1994
).
19.
D.
Appy
,
H.
Lei
,
C.-Z.
Wang
,
M. C.
Tringides
,
D.-J.
Liu
,
J. W.
Evans
, and
P. A.
Thiel
, “
Transition metals on the (0001) surface of graphite: Fundamental aspects of adsorption, diffusion, and morphology
,”
Prog. Surf. Sci.
89
,
219
(
2014
).
20.
C.
Kittel
,
Introduction to Solid State Physics
, 8th ed. (
John Wiley & Sons
,
USA
,
2004
).
21.
X.
Liu
,
M.
Hupalo
,
C.-Z.
Wang
,
W.-C.
Lu
,
P. A.
Thiel
,
K.-M.
Ho
, and
M. C.
Tringides
, “
Growth morphology and thermal stability of metal islands on graphene
,”
Phys. Rev. B
86
,
081414
(
2012
).
22.
J. W.
Evans
,
P. A.
Thiel
, and
M. C.
Bartelt
, “
Morphological evolution during epitaxial thin film growth: Formation of 2D islands and 3D mounds
,”
Surf. Sci. Rep.
61
,
1
(
2006
).
23.
Y.
Han
,
E.
Gaudrie
,
T. J.
Oliveira
, and
J. W.
Evans
,
J. Chem. Phys.
(in Press).
24.
M. C.
Bartelt
,
S.
Günther
,
E.
Kopatzki
,
R. J.
Behm
, and
J. W.
Evans
, “
Island-size distributions in submonolayer epitaxial growth: Influence of the mobility of small clusters
,”
Phys. Rev. B
53
,
4099
(
1996
).
25.
C. R.
Stoldt
,
C. J.
Jenks
,
P. A.
Thiel
,
A. M.
Cadilhe
, and
J. W.
Evans
, “
Smoluchowski ripening of Ag islands on Ag(100)
,”
J. Chem. Phys.
111
,
5157
(
1999
).
26.
P. I.
Cohen
,
G. S.
Petrich
,
P. R.
Pukite
,
G. J.
Whaley
, and
A. S.
Arrott
, “
Birth-death models of epitaxy
,”
Surf. Sci.
216
,
222
(
1989
).
27.
W. C.
Elliott
,
P. F.
Miceli
,
T.
Tse
, and
P. W.
Stephens
,
Kinetic Roughening During Ag Homoepitaxy, Nato Science Series B
, 1st ed. (
Plenum
,
New York
,
1997
), Vol.
360
.
28.
V.
Fournée
,
A. R.
Ross
,
T. A.
Lograsso
,
J. W.
Evans
, and
P. A.
Thiel
, “
Growth of Ag thin films on complex surfaces of quasicrystals and approximant phases
,”
Surf. Sci.
537
,
5
(
2003
).
29.
Z.
Luo
,
L. A.
Somers
,
Y.
Dan
,
T.
Ly
,
N. J.
Kybert
,
E. J.
Mele
, and
A. T. C.
Johnson
, “
Size-selective nanoparticle growth on few-layer graphene films
,”
Nano Lett.
10
,
777
(
2010
).
30.
J. A.
Venables
, “
Nucleation and growth of thin films
,”
Rep. Prog. Phys.
47
,
399
(
1984
).
31.
J.
Zhang
,
M.
Brehm
,
M.
Grydlik
, and
O. G.
Schmidt
, “
Evolution of epitaxial semiconductor nanodots and nanowires from supersaturated wetting layers
,”
Chem. Soc. Rev.
44
,
26
(
2015
).
32.
Y.
Kuo
,
P.
Yen
,
W.
Chen
,
S.
Chen
, and
S.
Yau
, “
In situ scanning tunneling microscopy study of cobalt thin film electrodeposited on Pt(111) electrode
,”
Electrochim. Acta
112
,
831
(
2013
).
33.
E.
Lundgren
,
B.
Stanka
,
M.
Schmid
, and
P.
Varga
, “
Thin films of Co on Pt(111): Strain relaxation and growth
,”
Phys. Rev. B
62
,
2843
(
2000
).
34.
P.
Allongue
,
L.
Cagnon
,
C.
Gomes
,
A.
Gündel
, and
V.
Costa
, “
Electrodeposition of Co and Ni/Au(111) ultrathin layers. Part I: Nucleation and growth mechanisms from in situ STM
,”
Surf. Sci.
557
,
41
(
2004
).
35.
B.
Voigtländer
,
G.
Meyer
, and
N. M.
Amer
, “
Epitaxial growth of thin magnetic cobalt films on Au(111) studied by scanning tunneling microscopy
,”
Phys. Rev. B
44
,
10354
(
1991
).