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.

## I. INTRODUCTION

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}

## II. DESCRIPTION OF THE METHODS

### A. Experimental details

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}

### B. Computational details: DFT

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.

## III. COMPUTATIONAL RESULTS: ADSORPTION ENERGY AND DIFFUSION BARRIER

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.

. | 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, *E _{d}*, 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.

## IV. EXPERIMENTAL RESULTS: ISLAND CHARACTERISTICS

### A. Island shapes

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.

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.

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.

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.

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.

### B. Island densities and size distributions

The island density, *N _{isl}*, 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.

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.

These observations are analyzed in Sec. V.

## V. THEORETICAL ANALYSIS OF ISLAND FORMATION FOR DEPOSITION OF DY/GRAHPITE

### A. Island densities

First we assess the island nucleation mechanism by comparing the observed value of *N _{isl}* 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 nm

^{2}. Multiple experiments for deposition of 0.12-0.15 ML provide quite consistent estimates for

*N*(see Fig. 6). Averaging these data yields

_{isl}*N*= (2.4 ± 0.3) × 10

_{isl}^{−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 *E _{d}* ≈ 48 meV, obtained from DFT analysis described in Sec. III. Assuming an attempt frequency of

*ν*= 10

^{12.5}s

^{−1}yields a rate of

*h*=

*ν*exp[ −

*E*/(

_{d}*k*)] for adatom hopping associated with terrace diffusion. A key parameter for nucleation is

_{B}T*h*/

*F*. Using the experimental values of

*T*= 300 K and

*F*= 0.013 ML/s, one has

*h*= 10

^{11.7}s

^{−1}and

*h*/

*F*= 10

^{13.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

*N*∼

_{isl}*θ*

^{(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 E_{rev} = (2/3)E_{d} + E_{NN}, where E_{NN} > 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 × 10^{12} s^{−1}, one requires that E_{rev} ≈ 0.55 eV or above. Thus, for Dy on HOPG with *E _{d}* ≈ 0.05 eV, one requires that E

_{NN}≈ 0.52 eV or above. A crude estimate of E

_{NN}as 1/6 of the bulk cohesive energy yields a value of 0.51 eV, and it is plausible that the true E

_{NN}for a pair of Dy on HOPG is stronger.

For *i* = 1, one has that *N _{isl}* ∼ (

*h*/

*F*)

^{−1/3}for fixed coverage, and simulations indicate that

*N*≈ 2.18 × 10

_{isl}^{−5}per site for

*θ*= 0.13 ML.

^{22,23}Thus, the experimental value of

*N*, 1.26 × 10

_{isl}^{−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 10^{13} s^{−1}, cannot be excluded, but one generally expects low prefactors for low diffusion barriers (as in this system). Alternatively, the lower experimental *N _{isl}* 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

*N*, we note the following. If monomers and dimers have comparable mobility, then

_{isl}*N*∼ (

_{isl}*h*/

*F*)

^{−2/5}

^{22,24}which implies a reduction in

*N*from the classic

_{isl}*i*= 1 value by a factor of

*r*= (

*h*/

*F*)

^{−1/15}= 0.12 for the above

*h*/

*F*= 10

^{13.6}, i.e., the reduction is too strong. For facile dissociation of dimers but stable trimers corresponding to

*i*= 2, one has that

*N*∼ (

_{isl}*h*/

*F*)

^{−1/2}, reducing

*N*for the classic

_{isl}*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 mobility

^{24}(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, *N _{isl}* =

*c*(

*h*/

*F*)

^{−0.351}

*θ*

^{0.325}. Choosing

*c*to recover the experimental

*N*at 0.13 ML and recognizing that

_{isl}*F*=

*θ*/10 in experiment, one obtains

*N*≈ 5.01 × 10

_{isl}^{−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” regime

^{22}). 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

*N*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 ∼

_{isl}*θ*/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.

### B. Island size distribution

Next, we analyze experimental observations for the island size distribution, *N _{s}*, which gives the density of islands of

*s*atoms, so that $Nisl=\u2211s>1Ns$ and $\theta =\u2211ssNs$. If

*s*≈

_{av}*θ*/

*N*denotes the average island size (measured in atoms), then the island size distribution is naturally written in scaled form,

_{isl}*N*≈ (

_{s}*N*/

_{isl}*s*)

_{av}*f*(

*s*/

*s*,

_{av}*θ*), where the scaling or shape function

*f*(

*x*) satisfies ∫

_{x>0}dx

*f*(

*x*) = 1 and ∫

_{x>0}dx

*xf*(

*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.

### C. Layer distributions

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

*P*=

_{j}*θ*−

_{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

*P*and film roughness

_{j}*W*satisfies

*W*

^{2}=

*d*

^{2}

*θ*, 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

*W*

^{2}≈ (1 − 2

*α*)

*θd*

^{2}.

^{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, *P _{j}*, 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.

. | . | Model . | ||
---|---|---|---|---|

. | Experiment . | α^{∗} = α = 0
. | α^{∗} = α = 0.65
. | α^{∗} = 0.69 α = 0.85
. |

P_{0} | 0.631 | 0.670 | 0.634 | 0.631 |

P_{1} | 0.276 | 0.185 | 0.278 | 0.281 |

P_{2} | 0.0904 | 0.096 | 0.072 | 0.081 |

P_{3} | 0.0024 | 0.036 | 0.015 | 0.006 |

. | . | Model . | ||
---|---|---|---|---|

. | Experiment . | α^{∗} = α = 0
. | α^{∗} = α = 0.65
. | α^{∗} = 0.69 α = 0.85
. |

P_{0} | 0.631 | 0.670 | 0.634 | 0.631 |

P_{1} | 0.276 | 0.185 | 0.278 | 0.281 |

P_{2} | 0.0904 | 0.096 | 0.072 | 0.081 |

P_{3} | 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 *P*_{4} = 0) consistent with experiment.

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 *P _{j}* is coupled only to

*P*

_{k<j}so they can be recursively solved exactly to obtain

As shown in Table II, the resultant values of *P _{j}* agree roughly with experimental values, but the fact that

*P*

_{0}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 *P*_{4} = 0. Numerical integration of the associated rate equations indicates that good agreement with experiment comes from choosing *α* = 0.65, with the resultant values of *P _{j}* 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.

## VI. DISCUSSION

One major conclusion from this work is that Dy island nucleation is homogeneous. This is based upon the measured value of *N _{isl}* 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 *N _{isl}* = 2 × 10

^{−2}nm

^{−2}at 0.4 ML, which is two orders of magnitude higher than our value of

*N*at 0.13 ML. While a somewhat higher value of

_{isl}*E*= 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

_{d}*N*. Instead, we propose that the explanation lies in the influence of the number of carbon layers,

_{isl}*n*, on nucleation, growth, and stability of metal nanoclusters. There is now broad evidence for this dependence, with

*N*generally decreasing as

_{isl}*n*increases.

^{6}This has been demonstrated in studies of Au, Ag, and Pd nanoclusters on

*n*-layer graphene sheets supported on SiO

_{2}/Si substrates. This trend, decreasing

*N*with increasing

_{isl}*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.

## VII. CONCLUSIONS

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.

## Acknowledgments

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.

### APPENDIX A: ADDITIONAL NOTES ON LAYER OCCUPATIONS

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 *W*^{2} = (1 − 2*α*) *θd*^{2} quoted in the text. However, reliable results for smoother growth with larger *α* must account for the modified equations, noting that clearly this expression for *W*^{2} cannot apply for *α* > 1/2.

### APPENDIX B: DATA ACCESS

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

## REFERENCES

**15**, Chap. 99, 61–189 (1991).