We develop a theory of surface ripples that can emerge spontaneously when an amorphous thin film is grown by oblique-incidence sputter deposition. For simplicity, we consider the case in which two diametrically opposed, broad atomic beams are simultaneously incident on the substrate and focus on the angles of incidence just above the threshold angle for ripple formation. At early times, the ripples are roughly sinusoidal in form, but as time passes, they become increasingly faceted as their wavelength and amplitude grow. The facet slopes and the form of the ripple crests and troughs are found analytically at long times, and these predictions are in good agreement with our simulations. Ultimately, the ripples become highly ordered in the transverse direction and few dislocations remain. However, because the ripple wavelength and the amplitude vary in the longitudinal direction, the ripples are not perfectly ordered, even at long times.
In a series of fascinating experiments, Voronov and co-workers deposited multilayer Mo/Si thin films by oblique-incidence sputter deposition.1 They discovered that the surface of the growing thin film spontaneously developed a rippled topography if the angle between the beam and the substrate normal exceeded critical value . The ripples were oriented so that, on average, their crests were perpendicular to the incident atomic flux. Similar results have recently been obtained by another group.2 These findings are not just of academic interest: multilayer coatings consisting of alternating layers of high and low Z-materials have high x-ray reflectivities and are important in extreme ultraviolet lithography, x-ray imaging, high resolution x-ray spectroscopy, and x-ray astronomy. Moreover, in the fabrication of many advanced optical devices, multilayers must be deposited on strongly curved or nanostructured substrates, and this inevitably leads to oblique-incidence deposition. When multilayer blazed gratings for extreme ultraviolet radiation or soft x rays are produced, for example, multilayers are deposited on substrates with sawtooth topographies.3,4
In the experiments of Voronov et al., the energies of the atoms arriving at the surface of the growing multilayer coating were on the order of 10 eV. This is below the sputter yield threshold for both Mo and Si, and so essentially no sputtering of atoms that had already been incorporated into the multilayer thin film occurred. Although the impinging atoms produced no sputtering, momentum transfer from the incident atoms to atoms near the solid surface would have led to inelastic displacement of atoms already resident in the solid, a phenomenon that is known as mass redistribution (MR).5–7 Molecular dynamics simulations of the deposition of diamond-like carbon thin films, for example, show that arriving carbon atoms with energies as low as 30 eV produce craters.6 MR has a smoothing effect on a solid surface for less than a critical value, but it has a destabilizing influence for higher angles of incidence.5–7 Voronov et al. argued convincingly that MR was responsible for the ripple formation they observed in their experiments.1
In this paper, we will develop a theory of surface ripples that emerge when an amorphous thin film is deposited by oblique-incidence sputter deposition. For simplicity, we will consider the deposition of a thin film composed of a single atomic species, rather than studying multilayer deposition as Voronov et al. did. We will also simplify the problem by analyzing the case in which two diametrically opposed, broad atomic beams are simultaneously incident on the substrate and focus our attention on just above the threshold angle for ripple formation . This makes the problem more tractable by ensuring that the surface slope remains small as time passes. As a result, no overhangs develop as the film grows and no part of the surface is shadowed from the incident atomic fluxes. In contrast, overhangs formed and shadowing occurred in the single-beam experiments conducted by Voronov et al. Shadowing is a nonlocal phenomenon that is particularly difficult to deal with theoretically.
Interestingly, we find that the ripples that form on the surface of the growing film develop a sawtooth form as time passes. The ripple crests and troughs are both rounded, but, in general, the extent to which the crests and the troughs are rounded differs. Thus, the surface ripples display vertical asymmetry, i.e., the appearance of the ripples changes if they are reflected about the substrate plane. At early times, the ripples are rounded and roughly sinusoidal in form, but as time passes, their wavelength and amplitude grow and they become more faceted. At long times, the ripples become highly ordered in the transverse direction and few dislocations are present. The ripple wavelength and amplitude vary in the longitudinal direction, however, and so the ripples are not perfectly ordered.
This paper is organized as follows. In Sec. II, we begin by considering the oblique-incidence sputter deposition of an amorphous thin film by a single atomic beam before moving on to the case in which two diametrically opposed beams are incident on the surface in Sec. III. An in-depth analysis of the dual-beam equation of motion (EOM) is carried out in Sec. IV. We discuss our results and place them in context in Sec. V and conclude in Sec. VI.
II. SINGLE-BEAM EQUATION OF MOTION
Consider the oblique-incidence sputter deposition of an amorphous thin film. We take the substrate surface to be nominally planar and place the origin on it. We orient the axis so that it is perpendicular to the substrate surface and so that it points away from the substrate. A single broad beam of atoms is now directed onto the substrate, and a thin film begins to accumulate on it, as shown in Fig. 1(a). For simplicity, we assume that the incident atoms all have the same energy and momentum. The angle that the beam direction makes with the axis will be called the angle of incidence and will be denoted by . We orient the axis so that it is parallel to the projection of the beam direction onto the plane. The incident atomic flux is then , where .
We assume that the impinging atoms are all incorporated into the growing film, which is taken to have uniform density. The incident atoms will be taken to have an energy low enough that essentially no sputtering takes place. They are assumed to be sufficiently energetic that appreciable MR occurs, however. In contrast, in molecular beam epitaxy, the incident atoms or molecules have low enough energies that MR is completely negligible.
We will employ a continuum description of the surface dynamics in which the position of an arbitrary point on the film’s surface is given by , where is the height of the point above the plane at time . The surface height is obtained by coarse-graining the detailed microscopic surface configuration and is assumed to be a smoothly varying function of its arguments , , and . We assume that no overhangs are initially present or develop as the film grows, so that is a single-valued function of and for all times . We will also assume that the slope of the surface remains small enough that shadowing of the incident atomic flux does not occur.
If there were no currents on the surface of the solid, the normal velocity of a point on the surface would simply be , where is the atomic volume and
is the outward-pointing surface normal. Therefore, the EOM is in this case
where the subscripts on denote partial derivatives. Equation (2) is of course easy to solve for arbitrary initial surface heights .
We now move on to consider the case in which there are surface currents. Let be the atomic current on the surface of the solid. is everywhere tangent to the solid surface. The equation of motion for the surface is
The surface current depends on the angle of incidence . Its value at time also depends in principle on the form of the entire surface at that time. In other words,
for and . We can simplify the EOM obtained by combining Eqs. (3) and (4) a little by setting , i.e., by transforming to a frame of reference that moves upward with the spatially averaged surface velocity. We obtain
Equation (5) completely specifies the dynamics of the surface but, in general, it is an exceedingly complicated nonlinear partial differential equation. In addition, how and depend on all of the derivatives , , , , is not known for any atomic species, energy or angle of incidence. Fortunately, the surface dynamics become much simpler just above the threshold for pattern formation, i.e., when is small and positive. We seek solutions to Eq. (5) of the form
, , and are the so-called “slow” variables and , , and are the corresponding “fast” variables. Heuristically speaking, Eqs. (6) and (7) say that for just above the critical angle , the height of the surface disturbance varies slowly in space and time. In addition, the spatial variation in the direction is more gradual than in the direction because ripples with their wave vector pointing in the direction develop for . An a posteriori justification for adopting the scaling ansatz given by Eqs. (6) and (7) will be obtained once we have arrived at an EOM that is well-behaved in the limit for the case in which diametrically opposed beams are incident on the substrate.
Let denote the derivative of with respect to its th argument after the semicolon evaluated for . We define in an analogous fashion. The way in which quantities like and are defined should now be apparent. We expand Eq. (8) in powers of and retain terms up to fourth order in . The EOM must be invariant under transformation and so the coefficients of terms that do not have this invariance must have coefficients equal to zero. We obtain
The constants , , , , , , and that appear in Eq. (9) all depend on .
The term that appears in Eq. (9) can be eliminated by transforming to a moving frame of reference. The resulting equation can then be recast in terms of dimensionless quantities, which yields an EOM that depends on three dimensionless parameters. This EOM does not have a well defined limit, however. Physically, the reason for this is that the scaling ansatz given by Eqs. (6) and (7) implies that is of order and so is small, but if a single beam is incident on the film, surface slopes of order one develop.1 This means that the scaling ansatz is invalid for the single beam case.
III. DUAL-BEAM EQUATION OF MOTION
We now turn our attention to the problem in which there are two diametrically opposed atomic beams, each with flux . These beams both have the polar angle but one has azimuthal angle and the other has an azimuthal angle of , as shown in Fig. 1(b). If only the beam that is incident from the left were present, the EOM would be Eq. (9). Conversely, if only the beam that is incident from the right were present, the EOM would be Eq. (9) with replaced by , i.e.,
Note that in addition to being invariant under transformation , the dual-beam EOM (11) is invariant under , as it must be.
Equation (11) holds for just above the critical angle , i.e., for small . Experiments performed with a single incident atomic beam show that ripples with their wave vector in the direction form for , while the surface remains flat for .1,2 We assume that this would also be the case if the experiments had been carried out with diametrically opposed beams. This assumption implies that is negative for , zero for , and positive for , at least provided that is not too large. Accordingly, for close to the critical angle, , where is a positive constant that does not depend on . On the other hand, is negative for because it has been assumed that ripples with their wave vector pointing in the direction do not develop for just below . Equation (11) may now be written
We must have , since otherwise arbitrarily short wavelengths are linearly unstable and the continuum description breaks down. In addition, cannot be positive since if it were, the slope of the surface would grow without bound. Equation (12) becomes
where , , , , and . The coefficients , , and in Eq. (13) are positive while is nonnegative. In addition, by replacing by if necessary, we can arrange for to be nonnegative.
The term in Eq. (13) describes the destabilizing effect of MR above the critical angle. Thermally activated surface diffusion could produce the term ,8 but we suspect that MR could also contribute to this term. MR is smoothing for disturbances in the transverse direction for all angles of incidence.7 This is the physical origin of the term . The term in Eq. (13) is familiar from the theory of mounding instability that can occur during molecular beam epitaxy. There the term results from the Ehrlich–Schwoebel effect and can lead to the formation of a faceted surface.9 In our problem, the slope dependence of the surface current produced by MR leads to the presence of this term. The term , on the other hand, is the so-called conserved Kuramoto–Sivashinsky (CKS) nonlinearity. Although this term was first encountered in molecular beam epitaxy,9,10 it is also believed to play a role in the deposition of amorphous thin films.11 The CKS nonlinearity tends to produce coarsening of the surface patterns, i.e., the characteristic lateral and vertical length scales increase with time. It also breaks the symmetry that would be present if were zero. Since there is vacuum above the surface and solid below, there is no reason that such a symmetry should exist.
Our EOM (13) conserves mass, i.e., is a constant of the motion. This is a result of transforming to a frame of reference that moves upward with the spatially averaged surface velocity of the growing thin film.
How slowly must the height of the surface disturbance vary with for Eq. (13) to be valid? Let be the characteristic size of the collision cascade, or, equivalently, the characteristic lateral length scale of the crater that the impact of a single atom produces on the surface of the growing thin film. The characteristic lateral length scale (or length scales) of the surface ripple must be large compared to for Eq. (13) to be a good approximation. If the surface disturbance is sinusoidal, it has a single characteristic lateral length scale—its wavelength. On the other hand, if the surface disturbance is periodic but the width of the crests differs from the width of the troughs , then the characteristic lateral length scales and must both be large compared to for Eq. (13) to apply.
IV. ANALYSIS OF THE DUAL-BEAM EQUATION OF MOTION
We introduce the dimensionless quantities
Making these substitutions in Eq. (13) and dropping the tildes results in an equation of motion with a single dimensionless parameter ,
where is defined by the relation . This parameter is a measure of the relative strength of the quadratic and cubic nonlinearities.
Because the substrate surface is nominally flat, Eq. (15) can be linearized at early times. The amplitude of a sinusoidal ripple with wave vector increases exponentially in time with the rate
[The ripple amplitude decays exponentially if is negative.] It follows that ripples with wave number and with their wave vector along the direction emerge shortly after the growth of the thin film begins.
If the surface height does not depend on the transverse coordinate , Eq. (15) reduces to
which we will refer to as the one-dimensional (1D) EOM. For the case , we set . Differentiating Eq. (17) with respect to yields
which is the 1D Cahn-Hilliard (CH) equation.12 This leads us to the conclusion that for long times the solution to Eq. (17) with tends to a state in which most of the surface has a slope nearly equal to one of the two selected values . The slope changes rapidly in interfacial regions between adjacent intervals where the slope is nearly constant. As some regions of nearly constant slope contract and then disappear and others expand, the pattern coarsens. For , on the other hand, Eq. (17) is the CKS equation, which has been studied as a model of amorphous thin film growth11 and the step meandering instability on a crystal surface.13 A family of periodic steady-state solutions to the CKS equation exists;11 these consist of concave, nearly parabolic segments that meet at “kinks.” These kinks are not discontinuities in but are instead relatively narrow regions where is negative. For generic nominally flat initial conditions, a nearly periodic pattern with the linearly selected wave number emerges at early times, but at longer times coarsening occurs: kinks merge and the average size of the parabolic segments grows in time.13
Let us consider the behavior of solutions to the two-dimensional (2D) EOM (15) with on the domain in which and and apply periodic boundary conditions. We introduce the effective free energy
will be referred to as the effective potential. Equation (15) can be written
where denotes the variational derivative of with respect to surface height . Equation (21) implies that , i.e., the effective free energy can never increase. The dynamics therefore tends to minimize the value of . The effective potential has minima at . Therefore, the surface will tend toward a state in which most of the surface has a gradient nearly equal to , i.e., the surface will facet. A flat facet with one of the two selected gradients has a free energy equal to zero. Adjacent facets are separated by “edges” that have a positive free energy per unit length. changes rapidly but not discontinuously as an edge is traversed.
The EOM (15) cannot be written in the variational form (21) for . Nevertheless, we will see that the surface tends to facet for . The selected gradients will be shown to be equal for , and so their common magnitude increases with . We will also see that there are no selected slopes for . Our first step toward these conclusions will be to study steady-state solutions to Eq. (15) that depend only on the longitudinal coordinate . This will be done in Sec. IV A.
We have an additional motivation for studying these steady-state solutions: our numerical integrations of the EOM (15) show that as time passes and the patterns coarsen, they grow close to steady-state solutions of increasingly long wavelength. Steady-state solutions to Eq. (15) of infinitely long wavelength are therefore of particular interest.
A. Steady-state solutions
We will seek steady-state solutions to Eq. (17) in this subsection. Setting in Eq. (17), recalling that , and integrating the resulting ordinary differential equation (ODE) with respect to , we obtain
where is a constant of integration. We will look for steady-state solutions to Eq. (17) that are even functions of , which means that must be an odd function of . Applying this condition to Eq. (22), it follows that . Then, making the replacements and in Eq. (22), we obtain
where the dots denote time derivatives. For , Eq. (23) is the EOM for an undamped, undriven Duffing oscillator.
For , Eq. (23) can be written
Equation (24) is formally identical to Newton’s second law for a particle with unit mass moving in potential . The potential has global maxima at and a local minimum at . Equation (24) has the first integral
where the constant of integration can be thought of as the energy of the particle. For , the particle oscillates about the origin between the classical fixed points at that are given by the equation . For , Eq. (24) has the solutions
where is an arbitrary real constant. In terms of original variables, Eq. (27) is
where is an arbitrary constant. Equation (28) is the well-known solution to the 1D CH equation that describes the interface between regions in two different phases. Here, it gives the form of the surface between two regions in which the slope is constant.
For , we set in Eq. (17), which yields
Integrating this equation twice gives , where and are constants of integration. We once again seek steady-state solutions to Eq. (17) that are even functions of , and so . Moreover, by setting and then dropping the bar, we obtain
Making the replacements and gives
where the effective potential is given by
has a global minimum at , a local maximum at , and tends to infinity for . Equation (32) has the first integral
where the constant can be thought of as the energy of a particle of unit mass moving in potential .
For , Eq. (32) has the solution
where is an arbitrary constant. In terms of the original variables, Eq. (35) is
which is a parabolic surface profile. Naturally, if we translate this surface profile up or down, we still have a steady-state solution to Eq. (17).
For , the particle oscillates between the classical fixed points and that are given by the equation for and 1. We choose to be the left classical fixed point, so that and as . The period of the oscillations increases as tends to zero from below and becomes infinite in the limit that . For close to , we have . Therefore, choosing the zero of time so that , we obtain
where is an arbitrary constant. Equation (38) gives the form of the surface close to one of the wave crests. The radius of curvature at the crest is , which tends to zero in the limit that , i.e., the limit in which the wavelength of the surface ripple tends to infinity. Thus, the ripple crests become increasingly sharp as the wavelength of the periodic steady-state solution grows large. An analogous analysis shows that the radius of curvature at the bottom of a trough tends to 2 as the ripple wavelength tends to infinity. In the limit that the wavelength becomes infinite, the surface takes on the parabolic form given by Eq. (36).
Equation (23) cannot be written in the form of Newton’s second law for . To make progress for any , we set and . This gives us a system of two first-order ODEs,
A temporally periodic solution to this system is equivalent to a spatially periodic, steady-state solution to Eq. (17).
The system of equations given by Eqs. (39) and (40) defines a flow in the phase space with coordinates and , as illustrated in Fig. 2 for three values of . The flow has three fixed points for which appear at = and . For , there is a single fixed point at . The phase plane has a region in which all of the trajectories are closed; outside the trajectories are unbounded. The boundary of is the separatrix.
Linearizing the systems (39) and (40) about , we find that the trajectories about this fixed point are approximately circular: , where is a small, positive constant. The corresponding solutions to Eq. (17) have the form , where the amplitude of the sinusoidal surface disturbance is small and is an arbitrary constant.
Inside the separatrix, there are nested closed orbits. The period of these orbits is finite, but if we look at orbits that are progressively closer to the separatrix, the period diverges. It takes infinitely long to traverse either of the heteroclinic trajectories that travel from one of the fixed points to the other.
In terms of the original variables, Eq. (41) is
where is an arbitrary constant. The slope tends to the constant values with magnitude for and it switches sign at . Note that for , Eq. (45) reduces to the result we obtained earlier, Eq. (28). Integrating Eq. (45) with respect to , we obtain
Equation (44) is
when expressed in terms of original variables.
Equation (46) describes the curved surface that separates two regions in which the slope differs and takes on the constant values . If the surface is convex (concave), the width of the curved region is (). The widths and are both equal to for . If , then and hence the concave curved region is wider than the convex curved region, as seen in Fig. 3. This lack of symmetry under transformation stems from the presence of CKS nonlinearity in Eq. (17), which breaks the symmetry. is a decreasing function of and tends to zero for , whereas is increasing and it tends to infinity in this limit.
For , one of the parabolic arcs that make up the separatrix reduces to and the other goes off to infinity. Thus, for , the separatrix is the line . Equation (39) shows that along this line and so . In terms of the original variables, this is . Integrating this, we obtain
where is an arbitrary constant. This is precisely the parabolic surface profile we found earlier for the case .
B. Simulations in 1D
We integrated the 1D EOM (17) numerically for a selection of values using the fourth-order Runge–Kutta exponential time-differencing (ETD4RK) method described by Cox and Matthews.14 Simulations were carried out on the interval , and periodic boundary conditions were employed. The initial condition was spatial white noise with an amplitude of 0.001. The linear terms in the EOM were calculated in Fourier space and the nonlinear terms were computed in real space using central finite differencing accurate to fourth order in the grid spacing.
As time passes, the pattern coarsens and the radii of curvature of the peaks decrease. In our numerical work, it was important to ensure that the chosen grid spacing was always relatively small compared to the width of the peaks. The peak widths are on the order of half the linearly selected wavelength at early times. The average peak width decreases as time passes, however, and tends to at long times. We therefore chose the spacing for the spatial grid to be approximately equal to . Because vanishes for , we were not able to carry out reliable simulations for values larger than about 0.8.
Figures 4(a)–4(c) show the results of a simulation of Eq. (17) with . At early times, the ripple wavelength is close to the linearly selected wavelength , and there is little evidence of a selected slope. The ripple pattern becomes progressively more faceted as time passes. The average width of the facets grows and the number of interfacial regions that separate facets decreases. This leads to a reduction in the effective free energy for the surface,
as shown in Fig. 5. As expected, the slope distribution has two pronounced peaks at at long times, as shown in Fig. 6(a) for . The wavelength and amplitude of the ripples vary with position—see Fig. 4(c) for example. This disorder in the pattern is not surprising, since the band of linearly unstable wave numbers is not narrow: Eq. (16) shows that it extends from to .
Even though Eq. (17) cannot be written in the variational form (21) for , the behavior we found for is in some respects similar to what we observed for . Figures 4(d)–4(i) show the time evolution of the surface for and . These figures show that the pattern coarsens with time, just as it does for . In addition, the slope distribution once again has pronounced peaks at long times [see Figs. 6(b) and 6(c)], showing that the surface facets.
For , the solution to Eq. (17) grows increasingly close to the steady-state solutions with infinite wavelength as the time grows long. This is illustrated in Fig. 7. In this figure, is plotted vs at time for three different values of . Also included in the plots are the separatrices given by Eq. (47). The agreement is excellent.
The similarity between the form of the solutions to Eq. (17) at long times and the steady-state solutions is further reflected in the agreement between the predicted and measured selected slopes. In our simulations, the selected slope was taken to be the positive value of where the peak in the slope distribution appears. This is compared to the predicted value in Fig. 8. Once again, the agreement is impressive.
In Sec. IV A, we found that the peaks on the surface are sharper than the troughs for . This is seen in Figs. 4(f) and 4(i), for example. In contrast, for , the radii of curvature of the peaks and troughs are the same, as illustrated in Fig. 4(c).
During the time evolution of the surface, there are long intervals in which little change occurs. These are punctuated by brief periods in which the surface changes rapidly and the number of wavelengths in the spatial domain decreases. For , this was already evident in Fig. 5. A free energy cannot be defined for , but we can plot as a function of time to provide an illustration of the punctuated time evolution, as we have done in Fig. 9 for three representative values of . Each spike in these plot coincides with a decrease in the number of wavelengths in the domain by one.
To compute the mean wavelength in the simulations, we found the average inter-peak and inter-trough separation distances, averaged them over 100 runs, and plotted the results as a function of time for various values of . We also used the root-mean-square deviation of the surface height from its mean value as a measure of the interface width. The results are shown in Fig. 10.
The 1D CH equation produces coarsening, but at long times the coarsening is very slow—it is logarithmic.15 For the CKS equation, in contrast, and at long times.13 Accordingly, we would expect that increasing would increase the rate that coarsening occurs. Our simulations show that the wavelength and the interface width do indeed grow more rapidly when is increased, as seen in Fig. 10. Our simulations do not give compelling evidence for a regime in which and display power-law scaling. In fact, the coarsening may be logarithmic for , just as it is for . This type of logarithmic behavior would be very difficult to discern in simulations.
C. Simulations in 2D
We carried out simulations of the 2D EOM (15) for a selection of values, again using the ETD4RK method. As in the 1D case, we adopted periodic boundary conditions and used spatial white noise with an amplitude of 0.001 as the initial condition. The grid spacing in the direction was chosen as before. Because there is no instability in the direction, a coarser grid spacing could be used along that direction with a negligible loss of accuracy.
The behavior in 2D is in most respects similar to that in 1D. Simulations for three different values of are shown in Fig. 11. For , the surface is faceted at sufficiently long times and the slope distributions have sharp peaks, as shown in Fig. 12. The selected slopes are once again close to the predicted value —see Fig. 8.
Cuts of the surface along the direction are shown in Figs. 12(d)–12(f) for . As before, the amplitude and wavelength of the surface ripples vary with , and, for , the peaks are narrower than the troughs.
At long times, the surface is again close to the steady-state solutions with infinite wavelength. Evidence for this appears in Fig. 13, where is plotted vs for three values of at time and the separatrices given by Eq. (47) are also shown. Additional evidence is given in Fig. 14, where is plotted vs for three simulations with different values of . In each case, there is little change in the surface past time .
In 2D, a phenomenon that is impossible in 1D occurs—annihilation of dislocations of opposite sign. An example of such an event appears in Fig. 15. This leads to faster coarsening in 2D than in 1D, as Fig. 10 shows. After an initial transient, does not display the same type of high, sharp spikes as does, because wavelengths in the ripple pattern can be removed by dislocation propagation and annihilation in 2D. At long times, no dislocations remain, and there is very little variation of the surface height in the transverse direction (see Fig. 11). However, as already noted, the amplitude and wavelength vary with , and so the patterns are not perfectly ordered.
D. Selected slopes and interfacial widths in the original units
In terms of the original, unscaled variables, Eq. (13) is
where , and Eq. (28) is
Equation (51) shows that the unscaled selected slope is . Interestingly, is independent of the coefficient of the CKS term . It also tends to zero as as approaches from above.
The widths of the curved interfacial regions between adjacent facets are in the original units. These widths diverge as as approaches from above. Because depends on and in turn depends on , the widths do depend on . In particular, and coincide for . The width of a convex curved region is a decreasing function of and tends to zero for . In contrast, the width of a concave interfacial region increases with and tends to infinity for .
As we have seen, the case in which diametrically opposed atomic beams are simultaneously incident on the growing thin film is much simpler than the single beam case because in the former case, the EOM must be invariant under transformation . In an experiment, it would actually not be necessary for two atomic beams to be used. Instead, the sample could be rotated periodically through increments about the -axis while a single atomic beam was obliquely incident upon it. If the time between rotations were made sufficiently small, the effect would be essentially the same as if diametrically opposed beams were concurrently incident on the film.
In addition to carrying out experimental work on thin film deposition by a single obliquely incident atomic beam, Voronov and co-workers introduced a theory of ripple formation that they observed.1 In their theory, the Carter–Vishnyakov model of MR5 was adopted. The Carter–Vishnyakov model, however, is unlikely to be correct for arbitrary local angles of incidence. In addition, Voronov et al. neglected the effect of the CKS term and considered only surface morphologies that do not depend on transverse coordinate .
Surface ripples have also been observed in oblique-incidence epitaxial growth of the (100) surfaces of metals.16–20 In contrast to the problem studied in this paper, the energies of the incident atoms are too low to induce appreciable mass redistribution in oblique-incidence epitaxy, and shadowing effects appear to be responsible for ripple formation.19,20
Equations that are closely related to the EOM we derived for the surface of the growing thin film with two incident beams (50) have been encountered in other physical contexts. For example, EOMs that describe the epitaxial growth of single crystal thin films and that include the Ehrlich–Schwoebel and CKS nonlinearities have been studied by Golubović and co-workers.9 Our EOM (50) is a special case of the EOM that they used to model the epitaxial growth of a (110) surface,21–24 although they did not analyze the case in which and there is no instability in the transverse direction.
When a (100) GaAs surface is maintained at a temperature in excess of its recrystallization temperature and is bombarded with a normally incident argon ion beam, highly ordered ripples form.25–27 Ou et al.25 modeled the time evolution of the GaAs surface using a variant of Eq. (50) in which the anisotropic terms and are replaced by their isotropic counterparts and . The term in our EOM (50) suppresses variations in the surface height in the transverse direction, and so the differences in the time evolution produced by Eq. (50) and the equation studied by Ou et al. are expected to be modest.
In the problem studied by Ou et al., the damage done to the crystal structure by the incident ions is rapidly annealed away, and the material remains essentially crystalline. The facets that develop on the GaAs surface are low-index crystallographic planes and so are a manifestation of the underlying crystal structure. The same is true of the facets that form in epitaxial growth. In contrast, in our problem, the deposited material is amorphous and the facets that develop are not low-index crystal planes. Instead, the facets in our problem represent a compromise between the destabilizing effect of MR for angles of incidence greater than and the leading-order nonlinear effect of MR, which gives rise to the term proportional to in Eq. (50).
In Ref. 28, it was shown that Eq. (50) is the EOM that governs the surface of an amorphous solid that is bombarded with diametrically opposed, obliquely incident noble gas ion beams if the energy of the incident ions is below the sputtering threshold. The method used to derive Eq. (13) in Ref. 28 differs from the method employed here and applies only if thermally activated surface diffusion is negligible and the ion flux is sufficiently small. In addition, the steady-state solutions and their relationship to non-steady solutions of Eq. (13) were not studied in Ref. 28.
In this paper, we developed a theory of surface ripples that emerge when an amorphous thin film is deposited by oblique-incidence sputter deposition. For simplicity, we considered the case in which two diametrically opposed, broad atomic beams are simultaneously incident on the substrate and focused our attention on angles of incidence just above the threshold angle for ripple formation . At early times, the ripples are roughly sinusoidal in form, but as time passes, they become progressively more faceted and their wavelength and amplitude grow. The magnitude of the facet slopes is proportional to for close to . The form of the surface between adjacent facets was found analytically for long times, and the theoretical prediction is in good agreement with our simulations. In general, the extent to which the ripple crests and troughs are rounded differs. At long times, the ripples become highly ordered in the transverse direction and few dislocations remain. However, because the ripple wavelength and amplitude vary in the longitudinal direction, the ripples are not perfectly ordered.
We are grateful to Kevin M. Loew, Daniel A. Pearson, Patrick D. Shipman, and Dmitriy L. Voronov for valuable discussions. This work was supported by Grant No. DMS-1814941 awarded by the U.S. National Science Foundation.
The data that support the findings of this study are available from the corresponding author upon reasonable request.