Surface Relaxation of High-Aspect-Ratio Nanostructures: Theory and

Feb 22, 2010 - Surface-diffusion-driven decay of high-aspect-ratio gratings: Existence of morphologically related classes. Marcos A. Madrid , Roberto ...
0 downloads 4 Views 4MB Size
Download PDF
J. Phys. Chem. C 2010, 114, 4603–4610

4603

Surface Relaxation of High-Aspect-Ratio Nanostructures: Theory and Experiments Marcos Federico Castez,* Paula Cecilia dos Santos Claro, Patricia Laura Schilardi, Gustavo Andreasen, and Roberto Carlos Salvarezza Instituto de InVestigaciones Fisicoquı´micas Teo´ricas y Aplicadas (INIFTA), Facultad de Ciencias Exactas, UniVersidad Nacional de La Plata-CONICET, Casilla de Correo 16, Sucursal 4, (1900) La Plata, Argentina ReceiVed: NoVember 24, 2009

We discuss some characteristic phenomena that take place when high-aspect-ratio metallic nanopatterns decay by surface diffusion. We compare results obtained by means of a continuous 2D model, a kinetic 3D Monte Carlo model, and also experimental results in which nanostructured Au surfaces decay during an annealing. In all of these cases, we found structures that are very stable against surface diffusion currents in which the surface area is much higher than the one associated with a compact material. According to our models, such structures usually are complex surfaces exhibing voids, tunnels, channels, etc., and they are very different from the compact shapes expected after the annealing of low-aspect-ratio metallic nanopatterns. Introduction Over the past few years the importance of controlling and changing shapes in the nanoworld has been emphasized.1–6 Thermal treatments are commonly used to change the microstructure of polycrystalline materials; however, they can also be used in a controlled way to modify morphological features at the nanoscale. High-aspect-ratio (HAR) nanostructures are employed as platforms for functionalization and sensing of different types of molecules and biomolecules because they provide, due to the high real surface area, amplification of different types of signals used to detect the presence, and to follow the transport of species in different environments. In particular, HAR nanostructured metallic substrates are emerging as good candidates for supporting optically and electrochemically active molecules or biomimetic systems such as phospholipid membranes allowing surface enhanced Raman scattering (SERS), and electrical detection of different analytes.7–9 As already discussed, a crucial point in the use of these nanostructures is their stability under different conditions of operation. In particular Au, Cu, and Ag nanostructures are susceptible to aggregation and sintering even at room temperature due to the high mobility of metallic adatoms.10 The behavior of HAR nanostructures is difficult to predict under the frame of the usual theoretical models of surface relaxation because most of them use the so-called “small slopes approximation” which is clearly inapplicable for these systems.11 Recently we have presented a multiapproach modeling of surface diffusion driven decay of three-dimensional high-aspect-ratio nanostructures by using different strategies,12 and several geometrical properties of these HAR systems during their kinetic evolution were pointed out. In this work we have studied the structural changes and the surface area decay originated by the surface-diffusion-driven process by means of numerical simulations and tested these results with those resulting from experimental evidence concerning the surface relaxation of high-area metallic nanostructured gold substrates. This system consists of a random array of 10-20 nm in diameter columns with variable height that can be adjusted by the growth conditions between few nano* E-mail: [email protected].

meters to the submicrometer range (typically 100 to 1000 nm). Nanostructured gold substrates have been used for electrochemical and Raman detection of submonolayer amounts of adsorbed molecules.13 We show that our simulational results qualitatively agree with the behavior of surface relaxation in this complex system. Simulational Methods We have implemented two complementary strategies to simulate diffusive processes that take place during the annealing of metals: the first one is based on the continuous mesoscopic theory of surface diffusion, whereas the second one, that is intrinsically discrete in nature, is an atomistic kinetic Monte Carlo (KMC) model. The continuous modeling directly deals with the numerical integration of the surface diffusion flow of surfaces, namely the Mullins’ equation14,15

Vn ) -K∆sC

(1)

where Vn is the outward signed normal velocity, ∆s is the intrinsic surface laplacian, and C is the local curvature. Temperature enters into the simulation through the dependence of the coefficient K on temperature. We have parametrized such dependence using experimental results reported for the annealing of polycrystalline Cu.1 The surface diffusion flow defined by eq 1 is conservative, meaning that the total volume remains constant. Moreover, under certain hypotheses concerning the smoothness of the initial condition, it can be proved that this surface diffusion flow makes the total surface area decrease monotonically.16 Finding a numerical procedure to model the motion of curves and surfaces is an intrinsically difficult problem, mainly due to the fact that eq 1 involves a fourth order derivative. Several studies have been reported in the past few decades, dealing either with numerical issues16,17 or with the existence and uniqueness of a certain class of solutions18,19 or even directly applied to surface science related problems.20 We have employed a numerical vectorial model, in which marker particles are used to define a grid to integrate (in a finite differences scheme) eq 1 in the two-dimensional case. The idea behind the model’s

10.1021/jp911176y  2010 American Chemical Society Published on Web 02/22/2010

4604

J. Phys. Chem. C, Vol. 114, No. 10, 2010

Castez et al.

implementation is the following: Let us take a sampling of the initial curve at certain points (“marker points”). The piecewise linear function passing through these marker points is taken as an approximation to the “real” curve and is used to compute the intrinsic laplacian of the curvature and afterward making the marker points evolve in one time step. To avoid well-known problems that marker points algorithms have, e.g., the numerical instabilities that appear when two or more marker points become too close or too far during the kinetic evolution, we have implemented a regriding procedure that gives us a way to keep the distance between marker points bounded. This regriding procedure simply consists in dynamically adding or extracting marker points whenever it is necessary to keep the intermarker distance within the specified bounds, and gives us a chance to study cases that can not be accounted by the ordinary static griding as, for instance, HAR curves, necking processes, etc. Further details on the model’s implementation can be found in ref 2. On the other hand, we have implemented a KMC model21,22 in which particles are able to perform diffusion events between nearest neighbors in an FCC lattice. Transition probabilities P(cini f cfin) from the configuration cini into cfin were taken according to the transition state theory, i.e.

[

P(cini f cfin) ) ν exp -

Eact(cini f cfin) kBT

]

(2)

where ν is the attempt frequency set equal to 5 × 1012 Hz, Eact(cini f cfin) is the activation energy for the transition, kB is the Boltzmann constant, and T is the absolute temperature. Throughout this paper, we have chosen a simple model for the activation energies, in which a harmonic dependence of the energy as a function of the reaction degree is assumed, and the activation energies are simply obtained by considering the crossing of two of such harmonic potentials displaced in one lattice-constant unit. A more detailed description of such model can be found in references,21,23,24 and here we shall restrict ourselves to give the final expression for the activation energy, namely

act Einiffin )

(

β Efin - Eini 1 + 2 β 2

)

2

(3)

where Eini and Efin are the energies of the initial and final states, respectively. Also we called β ) ka2, where k is the force constant of the harmonic wells and a is the lattice constant. Assuming that the energy at a particular state is Eb times the number of occupied nearest neighbors, the energy difference Efin - Eini can be expressed in terms of the variation of the number of occupied neighbors (∆n) caused by a diffusing particle, i.e., Efin - Eini ) Eb∆n. Throughout this work, we have taken a ) 0.3 nm, β ) 1 eV, and Eb ) -0.1 eV. That values were not fit to any particular metal, but they attempt to give a rough estimation of typical values found in real systems. Despite its simplicity, this model will give us a correct qualitative picture about the geometrical properties during the decay of nanopatterns. In fact, as it was pointed out in ref 25, morphologic properties of metallic patterns evolving by surface diffusion currents have shown a rather universal behavior, irrespective of the particular details of the activation energy model. Therefore, most of our conclusions in this work can be

reproduced by using different (physically reasonable) activation energy models. Theoretical Results under the Small-Slopes Approximation. Already in the pioneering work by Mullins,14,15 most attention has been drawn to the case of interfaces in the smallslopes approximation (this was mostly due to the resultant simplicity of the equations rather than to physical motivations). In fact, for a single-valued interface, i.e., one that can be written as a function z(x, y, t), the surface diffusion flow can be described (under the small-slopes approximation) by the following linear partial differential equation:11,26

∂z(x, y, t) ) -K∇4z ∂t

(4)

Due to its linear character, a general solution of eq 4 for an arbitrary small-sloped initial interface can be constructed by studying its Fourier modes and, after that, applying the linear superposition principle. For the sake of simplicity, let us consider the time evolution of a single Fourier mode that only depends on x, i.e., zini(x,y) ) A sin(kxx), where kx ) 2π/λ. It represents a periodic pattern of wavelength λ that varies along the x direction and is invariant along the y direction. The small-slopes condition requires that Akx , 1. It is immediate to verify (performing a simple substitution in eq 4) that such Fourier mode evolves with time in a shape-preserving way, according to the following equation:

z(x, y, t) ) Ae-Kkx t sin(kxx) 4

(5)

From eq 5, the sinusoidal pattern evolves in time exponentially decreasing its amplitude, with a lifetime proportional to λ4. Now, we shall consider a patch of the pattern of sides λ and Ly along the x and y directions, respectively. Let us calculate the surface area S(t) of this patch

S(t) )

∫0L

y

dy

∫0λ dx√1 + ∇2z(x, y, t)

(6)

By combining eqs 5 and 6 and after some elementary calculus, S(t) can be expressed in terms of the complete elliptic integral of the second kind E(k)27

S(t) ) 2Ly



(

λ2 + 4A2(t)E 2 π

1

λ2 1+ (2πA(t))2

)

(7)

where we are calling A(t) ) A exp(-Kkx4t). Equation 7 gives the surface area of the pattern’s patch for an arbitrary amplitude. However, to apply it to the case in which we are interested, namely the small-slopes case, we must expand E(k) in powers of A(t)/λ. Performing such expansion and ignoring higher order terms in A(t)/λ, we obtain the following expression for the surface area S(t), under the small-slopes approximation:

(

S(t) ) Ly λ +

π2 2 -2Kkx4t Ae λ

)

(8)

For the sake of comparison, it is preferable to define a normalized area rather than the absolute surface area S(t), by

Surface Relaxation of HAR Nanostructures

J. Phys. Chem. C, Vol. 114, No. 10, 2010 4605

Figure 1. Left side: Snapshots at successive instants (t1 < t2 < t3) of the interface for the initial condition represented by the rectangular wave corresponding to t1. Right side: (solid line) Time evolution of the normalized surface area (AN) during the simulation. (circles) Theoretical prediction for the time evolution of AN according to the eq 9. Simulated temperature was 400 K.

dividing S(t) by the area of the patch’s projection on the x-y plane, i.e.

AN(t) )

4 A2 S(t) ) 1 + π2 2 e-2Kkx t λLy λ

(9)

Equation 9 says us that the normalized surface area of a lowamplitude sinusoidal pattern decreases exponentially, asymptotically approaching the value 1. By applying superposition ideas, we can generalize this result for an arbitrary, although smallsloped, surface. Thus, the eq 9 contains the essential information supplied by the linear theory of surface diffusion regarding the evolution of the normalized area. Numerical Results Now we will present the simulational results for both modeling strategies, i.e. continuous and KMC models. In the first place, we have performed numerical simulations applying the continuous model in the small slopes case to test at which extent this model agrees with the theoretical expectations given by eqs 5 and 9. On the left side of Figure 1, we can see successive snapshots of the interface evolving from a rectangular pattern as initial condition. [It is important to mention that, although such initial condition should be considered as a singular one from the continuous point of view (due to the presence of points where the function is not differentiable), it is not singular in the framework of our finite differences approach, where no passage to the limit is undertaken.] The initial pattern is a periodic rectangular-shaped pattern with an amplitude (defined as equal to the peak-valley difference) A ) 50 nm and a wavelength λ ) 500 nm. The ratio ε ) A/λ will be named the aspect ratio of such a pattern. The fast filtering of high-frequency Fourier modes predicted by eq 5 turns quickly into sinusoidal the initial rectangular wave (notice that the y scale on the left side of Figure 1 was amplified, so, actually, the initial pattern has a small aspect ratio equal to 0.1). On the right side of Figure 1, we have plotted the time evolution of the normalized surface area (to be exact, we should say “normalized arc-length” rather than “normalized surface area” in the 2D model, although, for the sake of simplicity, we will use this last term in a generalized

Figure 2. Left side: Interface snapshots at successive instants (t1 < t2 < t3 < t4 < t5) for the initial condition represented by the rectangular wave corresponding to t1 (curves are offset for clarity). Right side: Time evolution of the normalized surface area (AN) during the simulation. Simulated temperature was 400 K.

sense) together with the theoretical expectation given by eq 9. Excepting for the short initial transient in which the fast filtering of high-frequency modes transforms the rectangular wave into a sinusoidal wave, the agreement between both curves is evident. When the aspect ratio of the initial pattern is further increased, the actual evolution of the system is far from the linear theory expectations. In fact, superposition ideas are no longer applicable and, even considering a sinusoidal initial profile with an aspect ratio ε ∼1,12 the evolution is no longer shape-preserving: instead, overhangs can be spontaneously developed and, consequently strong deviations from the prediction of eq 9 are expected concerning the time evolution of AN. For instance, if we consider the evolution of a rectangular pattern similar to that considered on the left side of Figure 1 (at time t1) but with a height increased by a factor of 20, morphological and kinetic changes from the results shown in Figure 1 are evident, as can be seen in Figure 2. On the left side of Figure 2, we can see that the initial singlevalued rectangular pattern evolves into a multivalued wavy interface. Such geometrical differences, comparing the cases under and far from the small slopes approximation, are also associated with kinetic differences, as becomes clear when comparing the right sides of Figure 1 and Figure 2. We can see on the right side of 2 that the time evolution of AN strongly disagrees with the exponential dependence expected in the linear case given by the eq 9. In fact, in this case the curve AN(t) vs t shows a fast decreasing at a first transient, followed by a slower decreasing concave downward region and then there is a new inflection point, where the final decreasing-exponential stage starts. Despite these differences, both curves on the right side of Figure 1 and Figure 2 asymptotically approach 1, i.e., in these cases the interface is surface-diffusion-driven to that interface with the minimum surface area that is consistent with the imposed boundary conditions: a flat interface. After a quick examination of Figure 2, particularly concerning the lateral growth that the structures exhibit, it is immediate to ask: What happens in the case in which the distance between the rectangular structures in Figure 2 is decreased? The answer to this question is provided by Figure 3, where we can see how, as a consequence of the lateral growth, two adjacent rectangular structures merge, producing necks (Figure 3b) and, after that, giving place to the formation of voids (Figure 3c) that remain below the external interface (of course, to see voids formation

4606

J. Phys. Chem. C, Vol. 114, No. 10, 2010

Castez et al.

Figure 3. (a-f) Interface snapshots at successive instants for the initial condition represented by the rectangular wave shown in (a). (g) Normalized surface area (AN) values associated to states shown in (a-f). Simulated temperature was 400 K.

is required, from the numerical point of view, the application of a regriding procedure to define physically realistic interfaces, since the Mullins’ eq 1 with no regriding prescription leads to unphysical self-intersected interfaces). The external interface and the voids evolve subsequently in an independent way, asymptotically approaching two particular stationary states:2 the external interface evolves approaching a flat interface and voids evolve into circular shapes (Figure 3d-f). Both circular and linear interfaces are stationary solutions of the 2D surface diffusion flow. It implies that, once the system reaches such a state, no subsequent evolution takes place and the normalized surface-area won’t get its minimum value 1, but it will be “frozen” at a larger value, as can be seen in Figure 3g. Surface diffusion alone is not able to relax the system beyond the situation shown in Figure 3f; therefore, other relaxation mechanisms (e.g., bulk diffusion) will be needed in order to further decrease the surface area. The continuous model, based on the eq 1 only accounts for surface diffusion processes. On the other hand, the kinetic Monte Carlo model described before, despite its simplified character, provides a more realistic description of the system, in the sense that it is intrinsically discrete, takes into account the crystalline structure of a typical metal, and diffusive effects are not limited to the surface since bulk diffusion of vacancies is also present. In Figure 4, we show results obtained from KMC simulations concerning the decay of a high-aspect-ratio rectangular pattern, in the case in which the rectangular structures are well-separated. The morphological aspects of the evolving interface are very similar to those found in the continuous counterpart (see Figure 2). Moreover, time evolution of the normalized surface area in the discrete case (Figure 4 (f)) is also qualitatively very similar to the evolution in the continuous case (see the right side of Figure 2). On the other hand, when the separation between the nanopattern walls is very narrow, we find again a behavior that is a natural extension to the 3D case of the predictions from the continuous 2D model: the interface evolves into a flat external surface but leaving holes in the bulk (Figure 5). It should be noted that, although the linear continuous theory predicts a flat surface as the asymptotic state of the decay process (see eq 5), it is evident that, because of the internal noise,28 discrete systems (both real systems and KMC models) will depart from such

Figure 4. (a-e) Interface snapshots at successive instants for an initial condition represented by the rectangular pattern shown in (a) evolving in the framework of the KMC model. Rectangular structures in the initial condition have 10 × 40 × 22 atoms each, laying on the top of a substrate of 40 × 40 × 3 atoms. (f) Time evolution of the normalized surface area (AN) during the simulation. Simulated temperature was 400 K.

asymptotically flat states. Evidently, discrete systems will have surface defects such as monoatomic vacancies, steps, grain boundaries, etc. in the equilibrium state. Despite these facts, for abuse of language, we will call them “flat surfaces” even in the discrete case, since we are not focusing on the increasing of the surface area because of such equilibrium defects, that is much more small than the other sources for the increasing of the surface area that we are discussing in this paper. The holes in Figure 5, panels b and c, are cylindrical due to the imposed periodical boundary conditions along the axes x and y. Also in this case surface diffusion drives the interface into a stationary state in which the normalized surface area is higher than the corresponding to a compact structure. Moreover, Figure 5d shows that AN have a nearly exponential decay dependence against time, although its saturation value is not 1, as it is under the small-slopes approximation (9). Surface diffusion alone is unable to drop out the interface from this state: to further relax the interface into surfaces with smaller AN, different relaxation mechanisms than surface diffusion (e.g., bulk diffusion) are required. However, at moderate temperatures, bulk diffusion has a slower kinetics than surface diffusion; therefore, a complete decaying of the interface mediated by bulk diffusion will require a longer period of time. While patterns in Figures 4 and 5 were straightforward extensions to the 3D case of 2D situations, in the sense that such patterns were invariant under translations along the y axis, it is interesting to consider “true” 3D patterns. So, in Figure 6 we consider a situation in which the initial condition is a periodic

Surface Relaxation of HAR Nanostructures

J. Phys. Chem. C, Vol. 114, No. 10, 2010 4607

Figure 6. (a-e) Interface snapshots at successive instants for an initial condition represented by the rectangular pattern shown in (a) evolving in the framework of the KMC model. Rectangular structures in the initial condition have 12 × 12 × 52 atoms each, laying on the top of a substrate of 45 × 45 × 3 atoms. (f) Time evolution of the normalized surface area (AN) during the simulation. Simulated temperature was 500 K.

Figure 5. (a-c) Interface snapshots at successive instants for an initial condition represented by the rectangular pattern shown in (a) evolving in the framework of the KMC model. Rectangular structures in the initial condition have 17 × 40 × 52 atoms each, laying on the top of a substrate of 40 × 40 × 3 atoms. (d) Time evolution of the normalized surface area (AN) during the simulation. Simulated temperature was 400 K.

array of HAR nanocolumns: Looking at the intermediate stages during its kinetic evolution (Figure 6b-e), we can see that a complex evolution takes place, where nanocolumns merge together while, at the same time, several tunnels connect the external surface with holes inside. Remarkably, despite that this is a case where we are very far from the assumptions of the linear theory of surface diffusion, an exponential decay of AN as a function of time is still a good approximation to describe the asymptotic behavior of the decay process (Figure 6f). [It is in order to mention that we are not attempting to extract information regarding the short time dynamics of our KMC simulations, since that such behavior is strongly dependent on the rather artificial initial condition imposed. However, as kinetic evolution by itself produces missorientations among nanostructured walls, the first stage of the kinetic evolution is considered as an equilibration step where the system rapidly relaxes. We are only extracting relevant physical conclusions at later times in the kinetic evolution.] Besides the morphological features of the obtained nanostructure, it is important to stress that the saturation value in such exponential decay is once again higher

than that expected under the framework of the linear theory of surface diffusion (9). In fact, as is shown in Figure 6f, the asymptotic value of AN is more than four times bigger than the associated to a compact flat structure. Experimental Methods Nanostructured Au substrates were prepared from polycrystalline Au wires or foils following the procedure described by Arvia et al.29 They were first cleaned with piranha solution, rinsed with Milli Q water, and introduced in a conventional 3-electrode glass cell using a saturated calomel electrode (SCE) and high area platinum foil as reference and counter electrodes, respectively. The polycrystalline wires and foils substrates were then anodized in 0.5 M H2SO4 solution for ta (1 < ta < 15 min) respectively at 2.40 V (vs SCE) in order to form a thick hydrous gold oxide layer. Immediately afterward, the hydrous oxide was electroreduced to metallic Au by applying a potential sweep from 2.40 to - 0.6 V at 0.025 V s-1. The electrode was maintained at the final potential for 4 min in the case of wires and for 7 min when the substrates were Au foils, in order to ensure the total electroreduction of the thick oxide layer. The nanostructured material is a “black metallic gold” consisting of 10-20 nm in diameter and 100-1000 nm in length columns (depending on ta, the anodization time that controls the thickness of the oxide layer) containing large amounts of pores. The real surface area of the nanostructured Au (after the anodization/reduction procedure) was voltammetrically estimated in the same electrolyte (0.5 M H2SO4 solution) by recording a triangular

4608

J. Phys. Chem. C, Vol. 114, No. 10, 2010

Castez et al.

Figure 8. Evolution of the normalized surface area of a gold overlayer with a large initial area, as a function of the immersion time in an acid electrolyte. In the final stage AN decays with a nearly exponential dependence (dotted line).

Figure 7. (a) Evolution of the normalized surface area of a nanostructured gold overlayer, determined by voltammetry, as a function of the immersion time in an acid electrolyte. (b) Detail of the slow concave descending region. Finally the system reaches an exponential decay (dotted line). (c and d) In situ STM images (75 × 75 nm2) of the surface taken in the acid electrolyte immediately afterward oxide reduction (c) and after a relaxation time of 10800 s (d).

potential sweep from - 0.25 to 1.60 V at V ) 0.1 V s-1. The real surface area was then calculated as A ) Q/qmon where Q is the charge involved in the Au oxide monolayer electroreduction after the anodizing procedure and qmon is the charge density related to the electroreduction of a gold oxide monolayer (qmon ) 0.44 mC cm-2). This method, which use the O atom as yardstick, is widely used in electrochemistry to obtain the real surface area of different metals.30 Typical initial surface area values used in this work were 3-30 cm2. The nanostructured Au was then allowed to relax in the electrochemical cell containing 0.5 M H2SO4 solution at open circuit potential and room temperature. The time evolution of the surface area was in situ measured at different times by using by the voltammetric procedure described above for the evaluation of the initial surface area. In situ STM images were taken at open circuit potential with a ECSTM Nanoscope III equipment using its electrochemical cell containing 0.5 M H2SO4, and Apiexon covered Pt-Ir tips. Experimental Results Here we analyze the surface relaxation of a nanostructured gold overlayer as a model system to test our numerical findings. This system consists of a gold overlayer prepared by oxidation at high potentials of a gold substrate in acid solution, which results in the formation of a thick hydrous oxide layer, followed by fast reduction of the hydrous gold oxide to metallic gold.29 Scanning electron microscopy of the overlayer cross section and scanning tunneling microscopy (STM) of the overlayer surface reveal a columnar-like structure formed by nanosized columns.31 The thickness of the metallic overlayer (the column height) can be controlled between 30 and 1000 nm by increasing the anodization time of the gold substrate, the real surface area of the overlayer increasing sharply with the overlayer thickness.30,31 Figure 7c shows a typical STM image of the overlayer gold surface taken in situ in the acid electrolyte immediately afterward oxide reduction. The nanosized columns (average size 10 nm) are separated by pores and channels (dark regions in Figure 7c). In fact, ellipsometric data have shown that the

density value of the electrochemically prepared nanostructured gold is nearly a half than the one corresponding to pure gold.32 Figure 7a shows the evolution of the normalized surface area of the overlayer, determined by voltammetry, as a function of the immersion time (t) in the acid electrolyte. The AN vs t plots shows an initial sharp decay, followed by a slow concave descending region (shown in more detail in Figure 7b), to reach finally an exponential decay (dashed line in Figure 7a). The overall decay plot qualitatively resembles those shown in Figures 4-6. The relaxation process results in the lateral growth of the columns, a fact that is clearly observed in the STM images shown in Figure 7d, taken in situ in the same region depicted in Figure 7c, after a relaxation time t ) 10 800 s. We also observe the decay in height (amplitude) of the columns that is reflected in the decrease in the color contrast (compare panels c and d in Figure 7). It is evident that the surface diffusion driven process turns the surface more compact eliminating small columns by merging and smoothening. Note that in this case the system approaches asymptotically the real surface area expected for a conventional granular surface, i.e., slightly higher than the geometric area of the system. In Figure 8 the AN vs t plot for an overlayer with a large initial area is shown. While the shape of the curve is similar to the ones depicted in Figure 7, the slow concave descending region is even more evident than that shown for the thinner film, in a similar way to what we observed in the continuous model when we increased the pattern aspect ratio (see the right sides of Figure 1 and Figure 2). Asymptotically, the decay curve approaches a saturation value ∼5, which implies that the system is unable to relax completely to a flat surface leaving the real surface area more than four times greater than the geometric area. Monoenergetic positron beam results point out that a strong positron trapping by defects exists as a consequence of the large amounts of pores in the columnar structure.33 This means that surface diffusion is not able to further relax the system from this limit in agreement with the results of our simulations. From the analysis of Figures 7 and 8, it can be concluded that the final state clearly depends on the initial excess of area that in turn depends on the column height and pore density. Discussion and Conclusions Summarizing, in the first place we have performed numerical simulations based on the continuous theory of surface diffusion. After a discussion about the exact solvable linear case (i.e., under the small-slopes approximation), we focused on the case of high-

Surface Relaxation of HAR Nanostructures aspect-ratio nanostructures. From these simulations we have shown that, depending on the characteristic separation length between two adjacent nanofeatures, we can conclude that surface diffusion can lead to a complete decay of the nanostructure or, when this distance is small enough, it can lead to the formation of voids distributed in the bulk of the solid volume. Such voids adopt circular shapes in the 2D case, and they do not longer evolve by surface diffusion currents; that is, they represent stationary state structures. Other relaxation mechanisms (e.g., bulk diffusion) are needed to further relax those structures into other ones with even lower energy. The results obtained under the continuous approach were complemented by kinetic Monte Carlo simulations on threedimensional nanostructures. These simulations were in qualitative agreement with the results found in the continuous case, in the sense that, after an intermediate regime and a final exponential decay, the interface adopts a shape that does not evolve or, if it evolves at all, its kinetic evolution is extremely slow when compared to the typical time-scales expected from the continuous theory of surface diffusion. These metaestable structures typically present voids and channels under the external surface, showing a richer variety of situations than those found in the two-dimensional case. In fact, depending on the imposed boundary conditions, such voids and channels can adopt spherical, cylindrical, or even more complex shapes. Nevertheless, it is expected that in real systems with no long-range symmetries\ a disordered mixture of several kinds of shapes should coexist. Our experiments regarding HAR nanostructured Au substrates evolving by surface diffusion have shown that such structures are unable to fully relax into a flat surface, staying in metaestable states during long periods of time (on the order of days). Although the lost of surface area, grain boundary formation, and remaining channels are known facts of surface relaxation of nanostructured gold, the final real surface area of these conventional granular systems does not exceed a factor 2 with respect to the geometric area.29 The final factor higher than 4 reached for the thick metal overlayers can only be associated to a more complex and disordered final state, that is accounted in our simulations. The kinetic evolution of the normalized surface area is in very good qualitative agreement with our numerical simulations. Such metaestable states are characterized by a normalized surface area that is several times the corresponding to a flat surface with no voids nor tunnels underneath, in agreement with our simulational results for both the 2D continuous model and the 3D KMC model. Although we cannot get cross-sectional views of the experimental systems, because the small sizes of the nanocolumns and due to the destructive character of cross sections, the time dependence of the measured total area, the information obtained from top views of the nanostructures, and the structure information after decaying given by the monoenergetic positron beam measurements are in qualitative agreement with our simulational results. It is worth discussing briefly the time scales involved in both our simulations and our experiments. For example, KMC data in Figure 6 have a relaxation characteristic time on the order of 10-7 s. On the other hand, in the case of the results using the continuous 2D modeling shown in Figure 2, the time scale is on the order of 106 s, whereas the typical relaxation time in our experiment was on the order 104 s. Such huge variations in time scales are partially due to the fact that we have not parametrized specially such models to fit with the values of our nanostructured Au used in the experiments. Moreover, to reduce the computing time involved in our simulations, we performed them using

J. Phys. Chem. C, Vol. 114, No. 10, 2010 4609 temperature values of 400 and 500 K, whereas our experiments were performed at ambient temperature, which implies a modification in time scales in a range of 2-4 orders of magnitude. However, the most important variation in time-scales is not rooted in the values of related physical coefficients, but it has a geometrical source: In fact, it is well-known that the characteristic periods of time in surface diffusion processes change, when the overall system’ size is rescaled in an amount b, as tfb4t. This implies that if we consider systems in a scale length in the order of the nanometers (as it is the case for our KMC simulations) and systems in a scale length in the order of micrometer, then the characteristic periods of time in the bigger system will be increased in a factor 1012. Among the simplifications (intrinsic to any modeling) implicit in our KMC simulations, it is worth mentioning that we have ignored grain boundaries and also the fact that HAR initial patterns in the real system certainly are expected to be structures exhibiting a wide distribution of geometrical characteristics much more complex than our periodic simulated patterns. However, it is important to stress that, despite these simplifications, our simulations give a quite good qualitative description about the behavior of our experimental systems. Such agreement is due to the fact that the effects and phenomena that we discussed have a geometrical origin, whose essential aspects can be correctly described by our phenomenological modeling. Very interesting experimental examples supporting the qualitative correctness of our morphological description (concerning both the continuous and discrete approaches) regarding the surface diffusion evolution of HAR nanostructures are provided by recent experiments on thermal annealing HAR semiconductor devices.34,35 Acknowledgment. This work has been accomplished as part of projects PICT 06-621 and PICT 05-32906 of ANPCyT, PIP 112-200801-00362 of CONICET and 11/X532 of UNLP (Argentina). Authors acknowledge CONICET (Argentina National Research Council) and UNLP (Universidad Nacional de La Plata) for their support. References and Notes (1) Andreasen, G.; Schilardi, P. L.; Azzaroni, O.; Salvarezza, R. C. Langmuir 2002, 18, 10430. (2) Castez, M. F.; Salvarezza, R. C.; Solari, H. G. Phys. ReV. E 2006, 73, 011607-1-011607-12. (3) Rieth, M. Nano-Engineering in Science and Technology; World Scientific: Singapore, 2003. (4) Son, C.-S.; Kim, T.; Wang, X.-L.; Ogura, M. J. Cryst. Growth 2000, 221, 201–207. (5) Murty, M. V. R. Phys. ReV. B 2000, 62, 17004. (6) Ballestad, A.; Tiedje, T. Phys. ReV. B 2006, 74, 153405. (7) Wu, Y.; Livneh, T.; Zhang, Y. X.; Cheng, G.; Wang, J.; Tang, J.; Moskovits, M.; Stucky, G. D. Nano Lett. 2004, 4 (12), 2337–2342. (8) Chaney, S. B.; Shanmukh, S.; Dluhy, R. A.; Zhao, Y.-P. Appl. Phys. Lett. 2005, 87, 031908. (9) Song, Y.-Y.; Gao, Z.-D.; Kelly, J. J.; Xia, X.-H. Electrochem. Solid State Lett. 2005, 8, C148–C150. (10) Zangwill, A. Physics at Surfaces; Cambridge University Press: Cambridge, 1988. (11) Barabasi, A. L.; Stanley, H. E. Fractal concepts in surface growth; Cambridge University Press: Cambridge, 1995. (12) Castez, M. F.; Salvarezza, R. C. Appl. Phys. Lett. 2009, 94, 053103. (13) Tognalli, N. G.; Fainstein, A.; Vericat, C.; Vela, M. E.; Salvarezza, R. C. J. Phys. Chem. C 2008, 112, 3741–3746. (14) Mullins, W. W. J. Appl. Phys. 1957, 28, 333. (15) Mullins, W. W. J. Appl. Phys. 1959, 30, 77. (16) Mayer, U. F. Computat. Appl. Math. 2001, 20, 361–379. (17) Chopp, D. L.; Sethian, J. A. Interfaces Free Boundaries 1999, 1, 1–18. (18) Escher, J.; Mayer, U. F.; Simonett, G. SIAM J. Math. Anal 1998, 29, 1419–1433. (19) Giga, Y.; Ito, K. Comm. Appl. Anal. 1998, 2, 393–405.

4610

J. Phys. Chem. C, Vol. 114, No. 10, 2010

(20) Coleman, B. D.; Falk, R. S.; Moakher, M. Physica D 1995, 89, 123–135. (21) Kang, H. C.; Weinberg, W. H. J. Chem. Phys. 1989, 90, 2824. (22) Fichthorn, K. A.; Weinberg, W. H. J. Chem. Phys. 1991, 95, 1090. (23) Kang, H. C.; Weinberg, W. H. Phys. ReV. B 1988, 38, 11543– 11546. (24) Castez, M. F.; Albano, E. V. J. Phys. Chem. C 2007, 111, 4606– 4613. (25) Castez, M. F.; Albano, E. V. Phys. ReV. E 2008, 78, 031601. (26) Wolf, D. E.; Villain, J. Europhys. Lett. 1990, 13, 389–394. (27) Abramowitz, M.; Stegun, I. A. Handbook of Mathematical Functions; Dover Publications: New York, 1970. (28) Kampen, N. G. V. Stochastic Processes in Physics and Chemistry; Elsevier Science B.V.: Amsterdam, The Netherlands, 1992. (29) Vazquez, L.; Bartolome, A.; Baro, A. M.; Alonso, C.; Salvarezza, R. C.; Arvia, A. J. Surf. Sci. 1989, 215, 171–189.

Castez et al. (30) Salvarezza, R. C.; Arvia, A. J. A Modern Approach to Surface Roughness Applied to Electrochemical Systems. In Modern Aspects of Electrochemistry; Bockris, J. O., Conway, B. E., White, R. E., Eds.; Plenum Press: New York, 1996; Vol. 28, pp 289-373. (31) Alonso, C.; Salvarezza, R. C.; Vara, J. M.; Arvia, A. J.; Vazquez, L.; Bartolome, A.; Baro, A. M. J. Electrochem. Soc. 1990, 137, 2161– 2165. (32) Vela, M. E.; Zerbino, J. O.; Arvia, A. J. Thin Solid Films 1993, 233, 82–85. (33) Macchi, C.; Somoza, A.; Mariazzi, S.; Brusa, R. S.; Vericat, C.; Vela, M. E.; Salvarezza, R. C. Phys. Status Solidi (C) 2009, 6, 2585–2588. (34) Lee, M.-C. M.; Chiu, W.-C.; Yang, T.-M.; Chen, C.-H. Appl. Phys. Lett. 2007, 91, 191114. (35) Sudoh, K.; Iwasaki, H.; Hiruta, R.; Kuribayashi, H.; Shimizu, R. J. Appl. Phys. 2009, 105, 083536.

JP911176Y