Self-assembly of a drop pattern from a two-dimensional grid of nanometric metallic filaments
Ingrith Cuellar, Pablo D. Ravazzoli, Javier A. Diez, Alejandro G. González, Nicholas A. Roberts, Jason D. Fowlkes, Philip D. Rack, Lou Kondic
SSelf–assembly of a drop pattern from a two-dimensional grid of nanometric metallicfilaments
Ingrith Cuellar, Pablo D. Ravazzoli, Javier A. Diez, and Alejandro G. Gonz´alez ∗ Instituto de F´ısica Arroyo Seco, Universidad Nacional del Centro de la Provincia deBuenos Aires and CIFICEN-CONICET-CICPBA, Pinto 399, 7000, Tandil, Argentina
Nicholas A. Roberts
Mechanical and Aerospace Engineering, Utah State University, Logan, Utah 84322, USA
Jason D. Fowlkes and Philip D. Rack
Center for Nanophase Materials Sciences, Oak Ridge National Laboratory, Oak Ridge, Tennessee 37381, USA andDepartment of Materials Science & Engineering,University of Tennessee, Knoxville, Tennessee 37996, USA
Lou Kondic
Department of Mathematical Sciences, New Jersey Institute of Technology, Newark, New Jersey 07102, USA
We report experiments, modeling and numerical simulations of the self–assembly of particle pat-terns obtained from a nanometric metallic square grid. Initially, nickel filaments of rectangularcross section are patterned on a SiO flat surface, and then they are melted by laser irradiation with ∼
20 ns pulses. During this time, the liquefied metal dewets the substrate, leading to a linear arrayof drops along each side of the squares. The experimental data provides a series of SEM images ofthe resultant morphology as a function of the number of laser pulses or cumulative liquid lifetime.These data are analyzed in terms of fluid mechanical models that account for mass conservationand consider flow evolution with the aim to predict the final number of drops resulting from eachside of the square. The aspect ratio, δ , between the square sides’ lengths and their widths is anessential parameter of the problem. Our models allow us to predict the δ –intervals within whicha certain final number of drops are expected. The comparison with experimental data shows agood agreement with the model that explicitly considers the Stokes flow developed in the filamentsneck region that lead to breakup points. Also, numerical simulations, that solve the Navier–Stokesequations along with slip boundary condition at the contact lines, are implemented to describe thedynamics of the problem. I. INTRODUCTION
Controlling the placement and size of metallic nanostructures is crucial in many applications [1]. For instance, thesurface plasmon resonance among metallic drops depends on a coordination between their size and spacing [2, 3]. Thus,the inclusion of this type of nanoparticles into photovoltaic devices has led to increased efficiency [4, 5]. Moreover,in the field of biodiagnostics and sensing, functionalized Au nanometric drops bind to specific DNA markers thuspermitting binding detection [6]. In general, the potential applications of organized metallic nanostructures are wide–ranging and include Raman spectroscopy [7, 8], catalysis [9], photonics [10] and spintronics [11]. A methodology togenerate and organize structures at the nanoscale is to take advantage of the natural tendency of materials to theself–assemble [12, 13]. By combining the fact that liquid metals have low viscosity and high surface energy withthe current highly developed nanoscale lithography techniques, we have a platform to study the governing liquid–state dewetting dynamics [14, 15] such as liquid instabilities [16] with the goal of directing the assembly of precise,coordinated nanostructures in one [17] and two [18] dimensions.In this work, we focus on the formation of a two–dimensional drop pattern starting from the pulsed laser–induceddewetting (PliD) [18] of a square grid of Ni strips on a SiO coated silicon wafer. To investigate the behavior ofthese melted square grids, we employ well established nanofabrication techniques and PliD. With this methodologyit is experimentally possible to precisely control the initial far–from-equilibrium geometry and the liquid lifetime viananosecond laser melting.Fluid dynamics is used to rationalize the experimental data, because the evolution (and instability) of the metal ∗ Electronic address: [email protected] a r X i v : . [ phy s i c s . f l u - dyn ] M a y shape occurs in liquid state, and we therefore focus on its analysis from the fluid dynamical point of view. Forsimplicity, we consider the liquid metal as a Newtonian fluid, and ignore the effects that evolving metal temperaturehas on the material properties. The models that are discussed are based on earlier ones developed in the contextof experiments involving the evolution of grids made of silicon oil filaments. While the scale of the experimentsconsidered in the current work is considerably smaller, we will see that the main modeling approaches developed forthe films of millimetric thickness are useful to describe the results on nanoscale as well.The metal geometry analyzed here, while related to the ones studied in previous works involving nanoscale metalfilms [17–19] provides new and interesting challenges and open questions. The process of breakup of an original gridinto filaments is of interest on its own. For example, one could ask whether a drop will form at the intersection points,or whether a dry spot will be present there? Does the answer depend on the deposited film filament? Once theindependent filaments form, can their evolution be described based on the stability analysis of an infinite cylinder?And finally, to which degree could the results be explained based on fluid mechanical models and simulations ofNavier-Stokes equations?We will discuss in the following section the complex procedure involving heating, melting, and consequent solidifi-cation of liquid metal filaments. While significant amount of previous work [16, 20] suggests that focusing on simpleisothermal Newtonian formulations leads to a reasonable agreement with experiments, it is not a priori clear thatsuch an approach is appropriate for a rather complex geometry of metal grids/meshes on the nanoscale. For example,some works suggest that thermal gradients leading to Marangoni effect may be relevant [21–24], although recent workfocusing on liquid filaments suggests that in this geometry they could be safely ignored [25].This paper is organized as follows. In the experimental section we give a brief outline of the setup. In followingsection, results are presented and we use three models to describe the instability observed in experiments. We startfrom the conceptually simplest linear stability analysis, and proceed to consider progressively more complicated modelsbased on mass conservation (MCM) and on a fluid dynamic model (FDM) for the evolution of the breakups. Next, wereport Navier–Stokes numerical simulations using an appropriate geometry based on the experiments and we discusssome special effects observed at the grid corners. Finally, we summarize the results and consider future perspectives. II. EXPERIMENTAL SECTION
Electron beam lithography followed by direct current magnetron sputtering is used to define square grids of Nistrips on Si wafers (coated with a 100 nm thickness SiO layer), which were later melted by nanosecond laser pulses.In this section, we describe the details of the experimental procedure. Electron beam lithography:
Focused electron beam exposure at 100 keV and 2 nA was conducted using a JEOL 9300electron beam lithography system on poly(methylmethacrylate) (PMMA, positive tone electron sensitive resist 495-A4provided by Shipley) in order to define the strips that will form the square grid. The PMMA resist was previouslyspin coated on a 100 mm diameter substrate rotated at 4000 rpm during 45 s. The spin coating process was followedby a 2 min, 180 ◦ C hot plate bake. An electron beam dose of 1000 µ C cm − was required in order to completely exposethe electron resist yielding well defined the thin film strips of the grid. A 495-A4 resist development was carried outin a 1:3 methyl isobutyl ketone (MIBK)–isopropyl alcohol (IPA) solution during 100 s followed by an IPA rinse inorder to expose the strips in the resist down to the underlying SiO layer. Any residual electron resist was removedby exposing the system to an oxygen plasma generated in a reactive ion etcher for 8 s (100 W capacitively coupledplasma, 10 cm min − O flow rate and a pressure setting of 150 mTorr). DC magnetron sputtering deposition:
An AJA International 200 DC magnetron sputtering system was used todeposit the Ni thin film strips. The process was carried on with a constant power deposition mode at 30 W at achamber pressure of 3 mTorr Ar, which was maintained using a gas flow rate of 25 cm min − . The sputter rate ofNi was 5 . − for a target–to–substrate distance of 5 cm. A wet, metal lift–off procedure consisting of theimmersion of the substrate chip in acetone for 1 min was used to dissolve unexposed resist. Thus, the lift–off of theunwanted metal layer surrounding the Ni thin film strip features was achieved. Subsequently, the substrate chip wasrinsed in acetone, afterwards in isopropyl alcohol, and finally blown dry using N gas to remove any remaining debrisfrom the substrate. No specific treatments were done to remove the native Ni oxide prior to laser irradiation, as noobvious influences of the native oxide have been observed in the assembly dynamics. Nanosecond, ultraviolet pulsed laser irradiation:
A Lambda Physik LPX–305 i , KrF excimer laser (248 nm wave-length) was used to irradiate and melt the Ni square grid. During irradiation, the substrate surface was normal tothe incident laser pulse. As a result, the top surface of the strips as well as the surrounding substrate surface wereirradiated. The incident beam size was on the order of ∼ , significantly larger than the grid area ( ∼ µ m ), andthus irradiated the grids in a uniform way. The pulse width of the laser beam was ∼
18 ns (FWHM). A beam fluenceof (200 ±
10) mJ cm − was used to melt the grid, and focusing of the output beam was required in order to achievethis fluence. All samples reported in this work were irradiated with (at most) 30 laser pulses.The Ni square grids were patterned with rectangular cross section strips of width w g = (162 ±
6) nm, and length L g (internal side of the grid squares) in the range (600 , h g = 5, 10, and 20 nm ( ± III. RESULTS AND DISCUSSION
A typical example of the initial state is shown in Figure 1a, which corresponds to h g = 10 nm and L g = 1587 nm.At the end of the first pulse, a single drop appears at the vertices, while shorter and narrower filaments with smallbulges at the ends are formed along the sides of the squares (see Fig. 1b). This structure is a consequence of aliquid–like behavior of the strips due to melting. For subsequent pulses, axial retractions from both ends (shortening)of the remaining filaments are observed in an iterative fashion. The resulting further bulges with the correspondingnew bridges lead to the formation of a certain number of drops along the sides of the squares. Figure 1c shows thepattern obtained after 5 pulses, when the evolution has almost finished.Once the initial strip has been melted, its rectangular cross section evolves into a cylindrical cap shape by parallelcontact line retractions (dewetting) thus leading to a narrower filament of width w = 2 (cid:114) h g w g θ − sin θ cos θ sin θ, (1)where θ is the static contact angle. For the present Ni/SiO system, we have θ = 69 ◦ ± ◦ [17]. Thus, we have w = (2 ± . (cid:112) h g w g . (2)A comparison between the measured widths and the calculated values given by Eq. (1) shows a very good agreement.For instance, for h g = 10 nm, we have measured an average width of (78 . ± .
6) nm, while the calculated value is w = (80 . ±
7) nm. (a) Initial state (b) 1 st pulse (c) 5 th pulse FIG. 1: (a) Initial square grid of Ni strips with rectangular cross section. The as-deposited metal thickness is h g = (10 ±
1) nm.The inner square sides are L g = 1587 nm long. (b) After the first laser pulse, the strips decrease their width by dewetting andthe ends detach from the vertices, leaving drops there. (c) Drops pattern resulting from the breakup of the filaments after 5pulses. Figure 2 shows the patterns observed after 5 pulses for three values of L g , with h g = 10 nm. Clearly, the numberof drops, n , along each side, decreases with L g (note that n does not include the drops at the vertices). Even if thereis a dominant value of n in each case, some dispersion of n is observed. We expect that the main reason for thisdispersion is the experimental noise leading to edge roughness of the strips. Figures 2(a) and (b) show that after5 pulses there are still few filaments which have not yet finished their breakup. Additional pulses will lead to fullparticle formation, but since there are only few we ignore them, and consider in our analysis only the filaments thathave finished evolving.Figure 3 shows histograms that indicate how many times a given number of drops, n , (resulting from each detachedfilament) is observed for a given grid, characterized by the nondimensional parameter δ = L g w . (3) (a) L g = 1387 nm ( δ = 17 .
21) (b) L g = 987 nm ( δ = 12 .
25) (c) L g = 606 nm ( δ = 7 . FIG. 2: (a) Drop patterns for h g = 10 nm and w g = 162 nm after 5 pulses for different initial lengths, L g . For instance, for the grid in Fig. 1, which corresponds to δ = 19 .
69 and h g = 10 nm, we have three bars, namely for 2,3 and 4 drops (Fig. 3b). However, the case with 3 drops is far more frequent (40 counts) than those for 2 and 4 drops(7 and 4 counts, respectively). So, the modal number for this grid is n = 3. Note that each pair ( δ, h g ) representsan experiment such as that in Fig. 1, so that 14 experiments are summarized in Fig. 3. For each one, we analyze allfilaments that have finished evolving; the total number can be found by adding the counts of the corresponding bars.For instance, the above mentioned case includes 51 filaments/edges.Figure 4 shows the modal number obtained from the experimental results presented in Fig. 3 as δ and h g are varied.In this figure, each symbol corresponds to an experiment: for example, the experiment discussed in the precedingparagraph is represented by the point ( h g , δ ) = (10 , .
69) for n = 3. The error bars correspond to uncertainty in w (see Eq. (2)). Note that each experimental point in Fig. 4 corresponds to the modal value of drops from a largenumber of filaments (between 50 and 300, depending on L g ). The various lines shown in Fig. 4 show the predictionsof various models that will be discussed later in the text. (a) h g = 5 nm (b) h g = 10 nm (c) h g = 20 nm FIG. 3: Number of filaments (count) that yield a certain number of drops as the thickness, h g , and the aspect ratio, δ = L g /w are varied. A. Linear Stability Analysis (LSA)
As a first attempt to estimate the emerging spatial scales due to breakups of the filaments, we consider the linearstability analysis (LSA) of an infinitely long filament. There are various approaches in the literature to such ananalysis; see [15, 26] for elaborate discussions. In particular, Fig. 5 from [15] compares several existing modelsshowing that the differences between them are mostly modest for contact angles smaller than π/
2. For simplicity, weconsider only the results obtained by carrying out LSA within long-wave (lubrication) theory, despite the fact that
E x p e r i m e n t a l d a t a L S A M C M F D M n d h g [ n m ] D D D D D D h g [ n m ] h g [ n m ] FIG. 4: Each symbol corresponds to an experiment carried out with a given value of thickness, h g , and the aspect ratio, δ = L g /w . The experiments are grouped according to the modal number of drops, n , as determined from the histograms inFig. 3. The error bars correspond to the uncertainty in w , see Eq. (2). The horizontal dashed lines are the predictions obtainedfrom linear stability analysis (LSA), showing (dimensionless) wavelength of maximum growth, λ m /w ( h g ), see Eq. (5). Thedot–dashed curves are obtained by the Mass Conservation Model (MCM), see Eq. (15). The solid horizontal lines correspondto the Fluid Dynamic Model (FDM), see Eq. (19). the contact angles in the present problem are not small. Such LSA yields the following expression for the critical(marginal) wavenumber, q c (see Eq. (27) in [15]), q c tanh( q c /
2) tanh(1 /
2) = 1 , (4)where q c = 2 πw/λ c , and λ c is the critical wavelength. The solution of this equation is q c = 2 . q RPc = 1. Even if it is usual to consider RP criterion as a firstrough approximation, it lacks other essential features of the problem, such as the contact line physics, and so Eq. (4)is more appropriate. Then, the distance between drops for the varicose unstable mode is given by the most unstablewavelength, λ m = √ λ c = 3 . w . Since all the grid sides are of the length L g , one expects that n is related to howmany times λ m fits into L g . Note that the experiments show that a drop is present at each corner, and then thesepositions must correspond to maximums of the perturbation, which in turn restricts the admissible Fourier modes ofthe LSA. Since n accounts for the number of internal drops (in between both corner drops separated by L g + w g ), wehave the condition ( n + 1) λ m ≤ L g + w g < ( n + 2) λ m . (5)Note that the upper bound for n drops corresponds to the lower bound for n + 1 drops.Figure 4 shows the predicted limits for δ as dashed lines for n = 1, 2 and 3. Clearly, the comparison with theexperiment shows that the straightforward LSA yields a too narrow range of values of δ for given n , leaving someexperimental points out of it. The agreement is lacking particularly for smaller values of δ , as expected since theLSA theory as presented assumes infinite filament length. Moreover, LSA is inconsistent with the the experimentalevolution of the system: the breakups do not occur simultaneously, but in a cascade process starting from the filamentends. Therefore, a model that takes into account these features of the problem is required. We proceed with discussingtwo of such models. B. Mass conservation model (MCM)
Here we consider the fluid mechanical description recently reported by Cuellar et al. [27] in the context of microfluidicexperiments carried out with the silicon oil grids. Although the scale of the experiments considered in [27] is different,visual similarity of the instabilities to the ones considered in the present paper suggests that the instability mechanismmay be similar. The macroscopic experiments from [27] provide however significantly more detailed information aboutgrid evolution, that is useful in explaining the present experimental results for which such detailed information is notavailable.Figure 5 illustrates the instability mechanism discussed extensively in [27]; here we provide a brief overview of themain features. The first step (between Fig. 5a and b), in which the rectangular cross section of the filament changesto a cylindrical one, has been discussed previously (see e.g. Eq. (1)). Afterwards, drops start developing at eachgrid intersection, so that cross–like structures are formed and bridge regions appear at their arms (see Fig. 5c). Thisprocess leads to bridge ruptures and the formation of detached filaments of length L i (see Fig. 5d). The length ofthese bridges, L a , is approximately equal to those of the arms of a cross that dewets to form a corner drop. Then, wecan write L i = L g − L a . (6) FIG. 5: Scheme illustrating the time evolution of a portion of square grid. See the text for the discussion of the various stagesof breakup process.
The detached filaments retract axially and bulged regions start forming at the ends (the stages between Fig. 5dand Fig. 5e). These bulges stop after having retracted a distance L d , so that the new filament length is L = L i − L d . (7)By comparing Fig. 5e with Fig. 1b, we expect that this is the stage achieved at the end of the first pulse. Thisexpectation is supported by the fact that the positions of the bulges in Fig. 1b are coincident with the two side dropsclose to those at the corners (see Fig. 1c). This is the case for most of the detached filaments observed in Fig. 1b.The bulges are connected to the filament by means of additional bridges, whose length is denoted by L b . The ruptureof these bridges leads to the formation of additional drops (see Fig. 5f), and the breakup process continues until onlyfinal drops remain, as also observed in [28].Clearly, the bridge breakup process is the key feature. When the bulges have achieved the equilibrium as in Fig. 5e,we assume a balance between the capillary pressure (due to the longitudinal and transverse curvatures) in the bulge,and that in the filament. Since both the bulge and the detached drop adopt approximately the shape of a sphericalcap (no visible hysteresis effects are present here, in contrast to the microscopic experiments [27], the balance yieldsthe value of the bulge size as (see Eqs. (10) and (11) in [27]) L h = 2 w, (8)so that the drop volume is V drop = π w (2 + cos θ ) sec θ θ . (9)Note that the filament of length L i can be thought as consisting of portions of length d , each one leading to theformation of a single drop. Thus, d can be calculated as the ratio between the volume drop, V drop , and the crosssection of the filament, A : d = V drop A = 8 π w θθ − cos θ sin θ sin θ θ . (10)In order to obtain n drops from a filament of length L i , we have nd ≤ L i < ( n + 1) d. (11)However, the value of L i is not readily available from the experiments, and needs to be estimated. We note first thatthe corner drops result from dewetting of a part of the grid at the intersections (see Fig. 5c). This part correspondsto a cross whose arms have length L a , so that its volume is V cross = A (4 L a + w g ). Assuming that the corner dropshave similar size as those resulting from the filament breakup, we can write V cross = V drop , and obtain L a = d − w g . (12)By using this result in Eq. (6), Eq. (11) gives (cid:18) n + 12 (cid:19) dw − w g w ≤ δ ≡ L i w < (cid:18) n + 32 (cid:19) dw − w g w . (13)In particular, for θ = 69 ◦ Eq. (10) yields d = (5 ± . w, (14)and using Eq. (2) to include the dependence of w on h g , we have5 (cid:18) n + 12 (cid:19) − (cid:114) w g h g ≤ δ < (cid:18) n + 32 (cid:19) − (cid:114) w g h g . (15)The expressions in Eq. (15) are plotted in Fig. 4 as MCM–curves δ versus h g for given n . In general, the boundsgiven by MCM have a better agreement with the experiments than LSA. w L b L h A B C (a) L l L s A B C C' (b)
FIG. 6: (a) Sketch of the head and bridge regions showing the parameters used in the model. (b) Sketch of the longitudinalsection showing the interpretation of both roots for L b , namely, L s and L l . C. Fluid dynamical model (FDM)
Next, we consider a model that includes analysis of the flow during a breakup [27]. Figure 6 illustrates differentstages that can be observed during filament evolution. First, the bulge at the filament end stops its axial retractionwhen it reaches a certain size, at which the bulge is at equilibrium with the unperturbed filament, see Fig. 6a. Letus call by A the static point where the bulge is connected to the filament. The connecting region (bridge) maydevelop a small disturbance in the form of a neck. As a bridge narrows, an axial Stokes flow develops there (in thenext subsection we confirm that inertial effects are not of relevance here). This flow is due to the dynamic balancebetween the viscous forces and pressure difference between the bulge and the depressed center of the bridge (A andB in Fig. 6a, respectively, where L b is the distance between them). This pressure difference occurs because of thedistinctive curvatures (longitudinal and transversal) at the bridge center B, and those at its ends (A and C, onlytransversal), where the curvatures are equal to that of an unperturbed filament. Note that the pressure at C is thesame to that at any point in the rest of the filament (e.g., C’ in Fig. 6b). Since both the head and the unperturbedremaining filament are at equilibrium, points A, C and C’ have the same pressure (note that A and C do not need tobe symmetric with respect to B). As the pressure at point B is different from that at the points A and C, there is anoutflow from B that further depletes the neck region leading to an eventual breakup. Requiring a balance betweenthe resulting Stokes flow and the pressure differences between B and the points A or C, one finds that there aretwo possible equilibrium distances from the center of the neck (B) to the points with the unperturbed pressure of astraight filament (say A and C’).As discussed in [27], there are two positive values of L b : one for a short bridge, L s , and another for a long one, L l ,namely L s = 0 . w, L l = 4 . w, (16)(both calculated for θ = 69 ◦ ). The smaller root, L s , corresponds to the distance between A and B. In order tounderstand the larger one, note that (as the flow develops in the neck, the point C moves away from B towards thefilament. Simultaneously, a new bulge starts to form (dashed red line in Fig. 6b). When the breakup occurs at B, thebulge dewets and grows (dot–dashed green line in Fig. 6b). Finally, C stops at C’, which is an equivalent point to A,because the new bulge (dotted blue line in Fig. 6b) is identical to the former one since it has reached the curvaturesneeded to be at equilibrium with the filament. Consequently, the distance between the fixed point B (breakup point)and C’ (where the static bulge and filament meet) corresponds to the second root, L l .Based on this interpretation of the second root, L l , we can write (see Fig. 5e) L l ≈ L d + L h , (17)as confirmed by the experiments in [27]. Therefore, the characteristic length of a filament needed for the formation ofa single drop is L = L l + L s = 4 . w . Figure 7a illustrates the introduced quantities. Note that L is conceptuallyequivalent to the length d in Eq. (10) (see also Eq. (14)). Although L s is derived from the FDM model for thebreakup of a single filament, it is close to the value found for the bridges (cross arms) that occur at the intersectionof perpendicular filaments, namely L a that was obtained in the MCM. It is noteworthy that while one model focuseson the mass conservation and the other one on the dynamic effects, they yield similar results.Within the dynamical model, if L i = L ≡ L , there is the possibility of generating two drops when the smallbridge between the two heads formed from both ends of the filament is long enough to allow for a breakup at adistance L s from each static bulge (see Fig. 7b). Following a similar reasoning, a general formula for the limits of L i that allow for the formation of n drops can be written as: L n ≡ nL , n = 1 , , . . . (18) L = L t = t > L
1t > t= (a)(b) L l L s L l L s L s L l FIG. 7: Sketches of the filament showing the parameters used in the model to define the limiting lengths of the filaments thatyield: (a) one drop, and (b) two drops.
Note however that these limits are only lower limits for the existence of a certain number of drops, not the upperones. For example, when L i is slightly below L , there is the possibility that both heads coalesce into a single drop.Then, the upper limit of one drop can be estimated as L . Regarding the upper limit for more drops, this coalescenceprocess could occur on both sides of the remaining bridge, and therefore its maximum length should be 2 L . Then,the upper limit for the formation of n drops can be written as ( n + 2) L for n ≥
2. In order to compare this modelwith the experimental data in Fig. 4, we use Eqs. (3) and (6) to define D k = kL + 2 L a w , k = 1 , , . . . (19)where L a = 0 . w, (20)for θ = 69 ◦ (see Eqs. (2), (12) and (14)). Although the model is based on rather rough approximations, the predictedlimits agree very well with the experimental data (see FDM lines in Fig. 4). These limits are horizontal lines becauseEq. (19) does not depend on w ( h g ), in contrast to LSA and MCM. Note that FDM predicts overlapping δ –intervalsfor the existence of a certain number of drops. For instance, for D < δ < D it is possible to find either 2 or 3 drops,as observed in the experiments.Summarizing, we have used three models to predict the number of drops and compared the predictions to experi-ments. The most accurate seems to be FDM, which takes into account the flow in the bridge regions. MCM, whichis simpler, provides less accurate results, while the LSA model is the least accurate. A substantial difference betweenFDM and both LSA and MCM is that the former takes into account the actual sequence of events that lead to thefinal droplet configuration, such as the axial dewetting, the bridge formation and breakup. This iterated sequencepropagates from both filament ends towards the center. On the other hand, both LSA and MCM assume simultaneousevolution of unstable varicose modes and breakups.The models discussed so far have not considered the time scales involved in breakup process. In the followingsection, we present numerical simulations that allow to discuss this time scale. D. Numerical simulations
In this section we discuss time dependent numerical simulations of the dewetting and breakup processes [29, 30]. Wewill see that the results of these simulations are consistent with both the models and experimental results discussed sofar. Furthermore, we will see that the time scales emerging from the simulations are consistent with the experimentalones, for reasonable values of the slip length that is used to define fluid/solid boundary condition. Although precisecomparison of the time scales between experiments and simulations is difficult due to the fact that only limited amountof information is available from the experiments, we find this consistency encouraging.For efficiency of the computations, we consider only a single square of the grid. Since the initial dewetting processthat evolves a strip from rectangular to circular cross-section is very fast, we do not simulate this process. Instead,we consider that at t = 0, the unit cell of the grid is formed by four filaments of a cylindrical cap shape, and of thelength L cyl = L g + w g (see Fig. 5b).The time evolution is obtained by numerically solving the dimensionless Navier-Stokes equation Re (cid:20) ∂(cid:126)v∂t + ( (cid:126)v · (cid:126) ∇ ) (cid:126)v (cid:21) = − (cid:126) ∇ p + ∇ (cid:126)v, (21)where Re = ργw/µ is the Reynolds number. Here, the scales for the position (cid:126)x = ( x, y, z ), time t , velocity (cid:126)v =( u, v, w ), and pressure p , are the w , t c = µw/γ , U = γ/µ , and γ/w , respectively. The fluid parameters for the meltedNi are: density ρ = 7 .
905 g/cm , viscosity µ = 0 . γ = 1780 dyn/cm. For the consideredexperiments we have w = 80 nm, so that t c = 0 . Re = 53. Note that this value of Re is a consequence of thescaling used for the characteristic velocity, U , which yields a capillary number Ca = µU/γ = 1. This is done since thevalue of U is not known a priori . As it will soon be seen (see e.g. the slopes in Fig. 11), the maximum dimensionlessflow velocities are of the order of 10 − , so that the actual Re and Ca are much smaller than the above numbers.The normal stress at the free surface accounts for the Laplace pressure in the formΣ n = − (cid:16) (cid:126) ∇ τ · ˆ n (cid:17) ˆ n, (22)where ˆ n and ˆ τ are the surface normal and tangential vectors. Since the surrounding fluid (e.g., air) is passive, weassume that the tangential stress is zero at this surface, i.e. Σ τ = 0.Regarding the boundary condition at the contact line, we set there a fixed contact angle, θ . As it is commonly donefor the problems involving moving contact lines, we relax the no slip boundary condition at the substrate through theNavier formulation (see e.g. [31]), v x,y = (cid:96) ∂v x,y ∂z at z = 0, (23)0where (cid:96) is the slip length. Based on the experimental comparison from the previous work [16] considering evolutionof liquid metals, we use the value of slip length (cid:96) = 20 nm; this choice is discussed further below.We use a Finite Element technique in a domain which deforms with the moving fluid interface by using the Arbi-trary Lagrangian-Eulerian (ALE) formulation [32–35]. The interface displacement is smoothly propagated throughoutthe domain mesh using the Winslow smoothing algorithm, which consists of mapping an isotropic grid in computa-tional space onto an arbitrary domain in physical space, and it is usually more effective than the Laplace smoothingapproach [36–38]. The main advantage of this technique is that the fluid interface is and remains sharp [39], whileits main drawback is that the mesh connectivity must remain the same, which precludes achieving situations with atopology change (e.g., when the filaments break up). The default mesh used is unstructured, and consists typicallyof 13 × triangular elements and 4 × tetrahedral elements.Since the problem is symmetric with respect to the axis of each filament, we consider only the interior of the square,and apply symmetry boundary conditions along its four sides. Figure 8 shows the time evolution of this square forthe parameters as in Fig. 1 where δ = 19 .
69. The three stages correspond to: (a) the breakups that lead to the cornerdrops and the formation of the detached filaments with the bulges at their ends, (b) the first breakups of these bulges,and (c) the final configuration with four drops along the grid side.Figure 9 shows the thickness profile along the symmetry line of one of the sides of the square. Since completebreakup cannot be simulated using the present numerical method, a remnant film remains between the drops. Notethat this particular case leads to four drops, consistently with some of the experimental outcomes (see Figs. 3b and1c), although most of the time this geometry leads to three drops (modal number n = 3). We have also observed inthe simulations that cases with slightly smaller δ lead to three drops. Moreover, a comparison of the simulations withFDM shows that this particular case is in the overlapping interval between n = 3 and n = 4.Figure 10 shows the results for the parameters corresponding to the experiments from Fig. 2, illustrating how adecrease of δ also leads to increased number of drops in the simulations. Moreover, in these cases the obtained valueof n is fully consistent with the experiments. (a) t = 40 , t c = 8 ns (b) t = 80 , t c = 16 ns (c) t = 105 , t c = 21 ns FIG. 8: Time evolution of the fluid thickness for the parameters of Fig. 1 using (cid:96) = 20 nm. Here, we have w = 80 nm and L cyl = 1749 nm. The lengths are in units of w . Regarding the choice of the slip length, (cid:96) = 20 nm, we note that the main expected influence of the value of (cid:96) is onthe time scale of the evolution. To discuss this issue further, we carry out simulations with (cid:96) in the range [1 ,
40] nm,and record the position of the maximum height in the bulge as a function of time. Figure 11 shows the correspondingresults, together with the resulting number of drops. Not only the time scale of the problem is affected by (cid:96) (smaller (cid:96) implies slower evolution), but also the final value of x h changes and, eventually, the resulting number of drops aswell. While precise comparison of the time scales between experiments and simulations is not possible since we do notknow exactly when the evolution stops in the experiments, it is encouraging to find comparable time scales betweenexperiments and simulations for a reasonable value of slip length.We are also in position to compare the simulation results with the models considered so far. Figure 11 shows(dashed line) the position of the maximum thickness at the bulge (see Fig. 5e) x maxh = L a + L d + L h , (24)where L a is given by Eq. (20). We can estimate x maxh by resorting to the FDM. According to Eqs. (16) and (17) we find x maxh ≈ . w , shown as the dashed line in Fig. 11. As a consequence, the values of (cid:96) that lead to a dewetting distance1 h / w t=0t=1 ns t= 8 ns x / w h / w t=16 ns x / w t= 21 ns x h (t) (a) (b)(c) (d) FIG. 9: Thickness, h , along a filament ( x -coordinate is defined along the symmetry line of a filament), at different times forthe grid shown in Fig. 8 with (cid:96) = 20 nm: (a) t = 0 and 1 ns, (b) t = 8 ns, (c) t = 16 ns, (d) t = 21 ns. Note the corner dropsat x = 0 and x/w = L cyl = 21 .
81. The arrow indicates the position of maximum thickness at the bulk, x h ( t ), further discussedin the text. (a) L g = 1387 nm ( δ = 17 .
21) (b) L g = 987 nm ( δ = 12 .
25) (c) L g = 606 nm ( δ = 7 . FIG. 10: Final numerical drop patterns for the experimental cases shown in Fig. 2 using (cid:96) = 20 nm. The results are shown atthe late times such that no further evolution is expected. of the filament end that are in agreement with the model are in the interval (1 ,
20) nm. Moreover, assuming that thebulge is at rest at the end of first pulse ( t ≈
18 ns), we consider that (cid:96) = 20 nm is an appropriate choice to accountfor the experimental data, consistently with the previous works [16] that considered similar type of experiments.
E. Further effects at the intersections
We observe in the experimental pictures that there are some cases where no drop is formed at the vertices (seeFig. 12). This anomalous effect is more frequent for smaller values of h g , e.g. h g = 5 nm in Fig. 2c. One explanationfor such behavior is an increased importance of the initial irregularities of the strip thickness, leading to instabilitiesof the free surface rather than those related to the contact line. Consequently, the position of the bridges could bealtered by other mechanisms, which could be more of a local character and less related to the symmetry of the system.Such anomalous behavior is particularly common for δ = 14 .
65 for h g = 5 nm and n = 1, and it is therefore notsurprising that this particular data point in Fig. 4 seems to be an outlier which does not agree with the proposedmodels. Careful analysis of the data shows that for ≈
37% of the vertices, a corner drop is missing for this particulargeometry.We can rationalize this effect by noting that there is a volume difference in the vertex region between the originalintersection of two strips with rectangular transversal section and the assumed cross with cylindrical cap arms after2 x h / w l = 1 nm l = 5 nm l =10 nm l =20 nm l =40 nmn=4n=4n=4n=4n=2 FIG. 11: Time evolution of the position of the point of maximum height in the bulge, x h , as a function of time for severalvalues of (cid:96) . The x –coordinate is measured from the filament intersection and n stands for the number of drops formed alongeach side of the square. The horizontal dashed line stands for x maxh as predicted by FDM.FIG. 12: Closeup of a SEM for a grid with h g = 5 nm and L g = 1556 nm ( δ = 24 . the fast initial dewetting stage (see Figs. 5a and b). In fact, the volume of the original intersection region, V = h g w g ,must be compared with the volume, V cross , of the cross region with cylindrical transversal section and width w (seeFig. 13a). Thus, the relative variation can be calculated as∆ VV = V cross − V V = w h g w g (cid:2) w − w g ) cot θ + 2 w cot θ + (3 w g θ sin θ − w ) csc θ (cid:3) − , (25)which is plotted in Fig. 13b as a function of h g for θ = 69 ◦ , w g = 160 nm and w as given by Eq. (2). Note that thisdifference can be as large as ≈ . h g = 5 nm, while it reduces significantly for larger values of h g , such as h g = 10or 20 nm.This volume deficit in the experiments may be the reason why the corner drop is frequently missing for small valuesof h g . Such a deficit implies the formation of either necks at the cross arms or a depression at the cross center. Ingeneral, the first option is more likely to happen, but the probability of the second one increases as h g decreases, since∆ V is so large that neck formation is not enough to compensate for it. For larger h g , this effect appears to be lessrelevant, since no missing corner drops are observed for h g = 10 and 20 nm, and only necks in the arms are formed. IV. SUMMARY AND CONCLUSIONS
In this work we report and analyze a series of experiments focusing on the formation of two–dimensional droppatterns by carrying out pulsed laser–induced dewetting (PliD) of a square grid of Ni strips on silicon wafers. By3 w g w (a) g (nm)0.40.50.60.70.80.91 ∆ V / V (b) FIG. 13: (a) Sketch showing the intersection region for the original grid (dashed lines) of width w g with rectangular crosssection of thickness h g , and the (liquefied) cylindrical cap filaments of width w . V corresponds to the square with thick lines,and V cross to the colored cross region. The dashed diagonals stand for the intersections of the cylindrical surfaces. (b) Relativevariation of volume in the intersection region as a function of h g (see Eq. (25)) when comparing the original volume V withthe assumed cylindrical cap arms as depicted in (a). means of well established nanofabrication techniques we are able to precisely control the initial far from equilibriumgeometry and the liquid lifetime via nanosecond laser melting. The results are presented as a series of snapshots(SEM’s) which show the grid evolution as the number of pulses is increased, and they are interpreted in terms of fluidmechanical models of increasing complexity intended to predict the number of drops that will result from the breakup.The models predictions are given by the curves in Fig. 4, which are compared with the data from 14 experiments fordifferent thicknesses h g .The advantage of this type of experiments is that they allow to study not only the two–dimensional structure ofthe grid as a whole, but also two other fundamental problems, namely, the formation of the corner drops as well asthe dewetting and breakup of short filaments. The modeling of these two phenomena has been combined with theanalysis of the experimental grid patterns. Moreover, the whole grid structure provides a large number of intersectionsand filaments (50 or more) under identical conditions, which is very convenient to verify repeatability and performstatistical analysis.The most basic approach is to use the results of the LSA for an infinitely long filament under long–wave approxi-mation. However, this attempt seems to be too crude for the present problem since its predictions do not comparewell with the data. We expect that the lack of agreement comes from the assumption of infinitely long filament, andnot from the use of long–wave approximation which is known to produce accurate results in this particular contextof filament breakup even for large contact angles [15]. A better approximation is obtained by resorting to a detailedmass conservation formulation that assumes that all drops can be represented by spherical caps. Finally, we obtainan even better agreement with the experiments by applying a fluid dynamical model (FDM), which was previouslysuccessful to account for similar experiments on microscopic scale [27], and that takes into account the dynamics ofthe filament breakup.The time evolution of the grid is also numerically simulated by solving the full Navier–Stokes equation assuming afixed contact angle and a given slip length, (cid:96) . In general, the numerical results regarding the final number of dropsalong each side of the grid agree with both experiments and model. By analyzing the position of the maximum at thebulge, we obtain that (cid:96) should be around 20 nm to obtain times scales comparable to those in the experiments. Thisvalue of slip length is consistent with earlier work [16], that also considered the evolution of liquid metals of nanoscalethickness.While the agreement between relatively simple models, simulations of Navier-Stokes equations, and experiments, ispromising, we note that additional effects could be relevant in the context of dewetting of liquid metal filaments (andother geometries) on nanoscale, such as thermal effects in the metal and substrate [40], as well as the phase changeprocesses. Studies that will include some of these effects are left for the future work.4 Acknowledgments
I. Cuellar and P. Ravazzoli acknowledge post–graduate student fellowships from Consejo Nacional de InvestigacionesCient´ıficas y T´ecnicas (CONICET, Argentina). J. Diez and A. Gonz´alez acknowledge support from Agencia Nacionalde Promoci´on Cient´ıfica y Tecnol´ogica (ANPCyT, Argentina) with grant PICT 1067/2016. P. Rack acknowledgessupport from NSF CBET grant 1603780. The experiments and the lithographic patterning were conducted at theCenter for Nanophase Materials Sciences, which is a DOE Office of Science User Facility. L. Kondic acknowledgessupport by NSF CBET grant 1604351. [1] F. Ruffino and M. G. Grimaldi, Phys. Status Solidi A , 1662 (2015).[2] N. N. H. S. Lal, W. S. Chang, S. Link, and P. Nordlander, Chem. Rev. , 3913 (2011).[3] F. Le, D. W. Brandl, Y. A. Urzhumov, H. Wang, J. K. N. J. Halas, J. Aizpurua, and P. Nordlander, ACS Nano , 707(2008).[4] H. Atwater and A. Polman, Nat. Materials , 205 (2010).[5] J. L. Wu, F. C. Chen, Y. S. Hsiao, F. C. Chien, C. H. K. P. L. Chen, M. H. Huang, and C. S. Hsu, ACS Nano , 959(2011).[6] N. L. Rosi and C. A. Mirkin, Chem. Rev , 1547 (2005).[7] J. N. Anker, W. P. Hall, O. Lyandres, N. C. Shah, J. Zhao, and R. P. V. Duyne, Nat. Mater. , 442 (2008).[8] T. Vo-Dinh, TrAC, Trends Anal. Chem. , 557 (1998).[9] P. Christopher, H. L. Xin, and S. Linic, Nat.Chem , 467 (2011).[10] E. Ozbay, Science , 189 (2006).[11] S. A. Wolf, D. D. Awschalom, R. A. Buhrman, J. M. Daughton, S. von Molnar, M. L. Roukes, A. Y. Chtchelkanova, andD. M. Treger, Science , 1488 (2001).[12] T.-S. L. C. Honisch, A. Heuer, U. Thiele, and S. V. Gurevich, Langmuir , 10618 (2015).[13] J. Koplik, T. S. Lo, M. Rauscher, and S. Dietrich, Phys. Fluids , 032104 (2006).[14] J. Lian, L. Wang, X. Sun, Q. Yu, , and R. C. Ewing, Nano Lett. , 1047 (2006).[15] J. Diez, A. G. Gonz´alez, and L. Kondic, Phys. Fluids , 082105 (2009).[16] Y. Wu, J. D. Fowlkes, N. A. Roberts, J. A. Diez, L. Kondic, A. G. Gonz´alez, and P. D. Rack, Langmuir , 13314 (2011).[17] J. D. Fowlkes, N. A. Roberts, Y. Wu, J. A. Diez, A. G. Gonz´alez, C. Hartnett, K. Mahady, S. Afkhami, L. Kondic, andP. D. Rack, Nano Lett. , 774 (2014).[18] N. A. Roberts, J. D. Fowlkes, K. Mahady, S.Afkhami, L. Kondic, and P. D. Rack, ACS Appl. Mater. Interfaces , 4450(2013).[19] C. A. Hartnett, K. Mahady, J. D. Fowlkes, S. Afkhami, L. Kondic, and P. D. Rack, Langmuir , 13609 (2015).[20] A. G. Gonz´alez, J. A. Diez, Y. Wu, J. D. Fowlkes, P. D. Rack, and L. Kondic, Langmuir , 9378 (2013).[21] A. Oron and Y. Peles, Phys. Fluids , 537 (1998).[22] A. Oron, Phys. Fluids , 29 (2000).[23] J. Trice, D. Thomas, C. Favazza, R. Sureshkumar, and R. Kalyanaraman, Phys. Rev. B , 235439 (2007).[24] J. Trice, D. Thomas, C. Favazza, R. Sureshkumar, and R. Kalyanaraman, Phys. Rev. Lett. , 017802 (2008).[25] I. Seric, S. Afkhami, and L. Kondic, Phys. Fluids , 012109 (2018).[26] K. Sekimoto, R. Oguma, and K. Kawasaki, Ann. Phys. , 359 (1987).[27] I. Cuellar, P. D. Ravazzoli, J. A. Diez, and A. G. Gonz´alez, Phys. Fluids , 102103 (2017).[28] A. G. Gonz´alez, J. Diez, R. Gratton, and J. Gomba, Europhys. Lett. , 44001 (2007).[29] G. Ghigliotti, C. Zhou, and J. J. Feng, Phys. Fluids , 072102 (2013).[30] R. M. S. M. Schulkes, J. Fluid Mech. , 277 (1996).[31] P. J. Haley and M. J. Miksis, J. Fluid Mech. , 57 (1991).[32] T. J. R. Hughes, W. K. Liu, and T. K. Zimmermann, Comput. Methods Appl. Mech. Eng. , 329 (1981).[33] J. Donea, S. Giuliani, and J. P. Halleux, Comput. Methods Appl. Mech. Eng. , 689 (1982).[34] K. N. Christodoulou and L. E. Scriven, Comput. Methods Appl. Mech. Eng. , 39 (1992).[35] C. W. Hirt, A. A. Amsden, and J. L. Cook, Comput. Methods Appl. Mech. Eng. , 203 (1997).[36] A. M. Winslow, J. Comput. Phys. , 149 (1966).[37] P. M. Knupp, Eng. Comput. , 263 (1999).[38] A. A. Charakhchyan and S. A. Ivanenko, J. Comp. Phys. , 385 (1997).[39] T. E. Tezduyar, Comput. Methods Appl. Mech. Eng. , 2983 (2006).[40] V. Ajaev and D. Willis, Phys. Fluids15