Molecular dynamics simulations of sodium nanoparticle deposition on magnesium oxide
Yannick Fortouna, Pablo de Vera, Alexey Verkhovtsev, Andrey V. Solov'yov
aa r X i v : . [ phy s i c s . a t m - c l u s ] N ov Molecular dynamics simulations ofsodium nanoparticle deposition on magnesium oxide
Yannick Fortouna , , Pablo de Vera , , ∗ Alexey Verkhovtsev , , and Andrey V. Solov’yov , MBN Research Center, Altenh¨oferallee 3, 60438 Frankfurt am Main, Germany Departamento de F´ısica–Centro de Investigaci´on en ´Optica y Nanof´ısica,Regional Campus of International Excellence “Campus Mare Nostrum”, Universidad de Murcia, 30100 Murcia, Spain On leave from A. F. Ioffe Physical-Technical Institute,Polytekhnicheskaya 26, 194021 St. Petersburg, Russia and Currently at Department of Materials Science and Engineering,University of Ioannina, 45110 Ioannina, Greece
The interaction of mass-selected atomic clusters and nanoparticles with surfaces attracts stronginterest in view of fundamental research and technological applications. Understanding dynamics ofthe deposition process is important for controlling structure and functioning of deposited nanopar-ticles on a substrate, but experimental techniques can usually observe only the final outcome ofthe deposition process. In this paper, the deposition of 4 nm-sized sodium nanoparticles on anexperimentally relevant magnesium oxide substrate is studied by means of classical molecular dy-namics simulations. An empirical force field is derived which accounts for the interaction of highlypolarizable Na atoms with the surface, reproducing the results of previously reported quantummechanics/molecular mechanics simulations. Molecular dynamics simulations permit exploring thedynamics of deposited nanoparticles on long timescales on the order of hundreds of picoseconds, thusenabling the analysis of energy relaxation mechanisms and the evolution of nanoparticle structureup to its thermalization with the substrate. Several nanoparticle characteristics, such as internalstructure, contact angle, and aspect ratio are studied in a broad deposition energy range from thesoft landing to multi-fragmentation regimes.
I. INTRODUCTION
Metal clusters, nanoparticles and nanoalloys have beena subject of intense research over the past decades [1–8].Unique and size-dependent structural, electronic, opticaland magnetic properties of these systems have led to var-ious technological applications. For instance, metal clus-ters and nanoparticles, both monatomic and bimetallic,can be used as junctions in nanoelectronic devices [9] oras elements of photonic crystals [10]. They are also use-ful for energy, environmental and medical applications,e.g. as catalysts [11], contrast agents in medical imag-ing [12], and radiosensitizers in cancer treatment withionizing radiation [13–15].Many of these applications involve the interaction ofclusters with molecular environments or with surfaceswhich serve as a support [4]. Understanding the dynam-ics of cluster deposition is of high relevance for control-ling the structure and properties of supported clusters.Depending on deposition conditions, structure of clus-ters on the substrate can be either preserved or changedsubstantially, or clusters may experience fragmentationincluding possible penetration into the substrate and/ormodification of the latter. The ability to control theseprocesses lays in the core of key experimental techniquesfor the fabrication of thin films and nanodevices, such ascluster ion beam assisted deposition, sputtering, surfacesmoothing or substrate implantation [16, 17]. ∗ [email protected] The shape of clusters and nanoparticles on a surfaceis determined by the interplay of different processes andphenomena, including electron shell closure [18–20], in-teraction of deposited systems with the substrate [21, 22],relaxation of thermal energy remaining after the collisionby means of heat transfer, and atomic rearrangementscaused by collision-induced mechanical stress of phasetransitions [5, 6]. These processes depend on the typeand temperature of the substrate, the size and compo-sition of the deposited cluster/nanoparticle and on thedeposition energy. For small atomic clusters containing N .
200 atoms quantum effects, such as even-odd oscil-lations in cluster abundance spectra and the appearanceof “magic” numbers associated with electron shell clo-sure, play a crucial role in determining the shape of clus-ters on a surface [19]. However, these effects shrink withincreasing the system size up to N ∼ [6]. While clus-ters made up of a few atoms usually keep their structureupon soft landing (i.e. when kinetic energy per atom ismuch smaller than the cluster cohesion energy) [23, 24],large clusters and nanoparticles can experience signifi-cant deformations, such as flattening, surface wetting orepitaxial alignment [16]. Hard deposition at kinetic en-ergies exceeding the cohesion energy can lead to clusterfragmentation.The deposition of metal clusters and nanoparticles onvarious surfaces has been widely studied experimentally(see e.g. Ref. [16] and references therein). However,experimental studies are usually limited to the observa-tion of the final state of the deposition process, and thestudy of the cluster deposition mechanisms thus com-monly relies on theoretical methods [16, 25]. For in-stance, detailed quantum mechanics/molecular mechan-ics (QM/MM) simulations of the deposition of smallsodium clusters on magnesium oxide substrates were re-ported in Refs. [24–27]. These materials were selected tohighlight the importance of both electron and nuclear dy-namics in the process of cluster deposition. Sodium is ahighly polarizable metal having a single valence electron.Magnesium oxide is a hard ionic crystal with a highly cor-rugated potential energy surface which can easily polarizeNa atoms and produce complex interaction patterns.The study of metallic aggregates deposited onto oxidesurfaces is very important for technological applications,especially in the field of catalysis [28]. Deposition, dy-namics and diffusion of metal clusters on MgO films hasparticularly attracted both experimental and theoreticalinterest [28–30], including a number of recent studies [31–34].QM/MM simulations conducted in Refs. [24–27] ex-plored different deposition regimes for small Na and Na clusters including soft landing, hard collision and reflec-tion of the clusters. Structure and dynamics of the clus-ters as well as the mechanisms of energy transfer to thesubstrate were analyzed in detail. However, due to thecomplexity and computational cost of QM/MM calcu-lations, this analysis could only be carried out for verysmall clusters, short simulation times (on the order ofseveral picoseconds) and zero temperature.Contrary to ab initio approaches, classical moleculardynamics (MD) permits simulating deposition of muchbigger systems at finite temperature and over signifi-cantly longer timescales, thus allowing to explicitly ac-count for the process of energy relaxation [22, 35–39].Empirical force fields employed in MD simulations canbe tuned to effectively reproduce the collision dynamicsobserved in QM/MM calculations on much shorter timescales.This paper reports a detailed theoretical analysis ofthe deposition of Na nanoparticles ( ∼ II. COMPUTATIONAL METHODOLOGY
In the performed classical MD simulations the cou-pled Langevin equations for all atoms in the system weresolved numerically by means of the leapfrog algorithm[40]. Simulations of the deposition of a 4 nm diame-ter sodium nanoparticle on MgO(001) surface were per-formed for systems pre-equilibrated at 77 K and 300 K.These temperatures allow the comparison of the depo-sition dynamics for a solid nanoparticle and a liquiddroplet, having into account that the melting temper-ature of Na is slightly below 300 K [41].All simulations were performed by means of MBN Ex-plorer [42], a software package for the advanced multi-scale modeling of complex molecular structure and dy-namics. The MgO substrate and the sodium nanoparti-cle were constructed by means of its dedicated graphicaluser interface, MBN Studio [43]. This software was usedalso to prepare all other necessary input files and to an-alyze simulation results. In the following sections, theconstruction and preparation of each part of the systemis described, together with the potential used for eachinteratomic interaction.
A. Magnesium oxide substrate
The interactions involving Mg and O atoms were de-scribed based on the nuclear contribution to the empiricalpotential which was defined in earlier QM/MM calcula-tions [24, 25]. Details of this potential, being a combi-nation of exponential and power potentials, as well asits parameters can be found in Table 4 of Ref. [25]. Allatoms in the substrate carried partial charges of ± | e | (with e being the elementary charge) and thus interactedalso through the Coulomb potential. The electrostatic in-teractions were treated by means of the Ewald algorithmimplemented in MBN Explorer [40].A face-centered cubic structure of MgO with a latticeparameter of 4.212 ˚A was employed to create substrateswith the size of 12 . × .
21 nm and 24 . × .
43 nm in the x - y directions, simulated using periodic boundaryconditions. Following Ref. [44] the substrate was formedby seven atomic layers in the z -direction normal to thesurface. The height of the constructed substrate (12.6˚A) exceeded significantly the range of Na–Mg and Na–Ointeratomic interactions (see next subsection).The structure of MgO was optimized using the veloc-ity quenching algorithm and the time step of 0.1 fs. Af-ter structure optimization the substrate was equilibratedusing the Langevin thermostat with a damping time of0.1 ps to the target temperatures of 77 K and 300 Ksuch that obtained atomic velocities corresponded to theMaxwell-Boltzmann distribution. O site, QM/MM [24] Mg site, QM/MM [24] O site, this work Mg site, this work E ad s ( e V ) r Na-i (¯)
FIG. 1. Adsorption energy of a Na atom on the MgO sur-face as a function of the distance between Na and O or Mgsites. Present MD empirical force field calculations (lines) arecompared to QM/MM results (symbols) [24].
B. Sodium atom–substrate interaction
The interaction between Na and O or Mg atoms wasdescribed by means of the following pairwise potential: U ij ( r ) = D ij h e − β ij ( r − r ,ij ) − e − β ij ( r − r ,ij ) i − C ij r + q i q j ε r , (1)where r is the distance between atoms i and j . The firstterm on the r.h.s. describes a Morse-type interaction,whereas the second and third terms describe the long-range polarization and electrostatic interactions, respec-tively. Parameters D ij , r ,ij and β ij of the Morse poten-tial represent the depth of the potential well, equilibriuminteratomic distance and steepness of the potential foreach pair of atoms. q i and q j are atomic partial charges, ε is the effective charge screening factor (set equal to1 in present simulations) and C i is an empirical param-eter determining the intensity of the polarization forcesfor the interaction of Na with a particular atom i . Thelong-range Coulomb interaction was calculated by meansof the Ewald algorithm [40]. The partial charge q Na onthe Na atom was treated as an additional free parameterto account for additional attracting forces between theatom and the surface due to its strong polarizability inthe field of the ionic crystal.The parameters for the potential (1) were obtained bya trial and error procedure in an iterative manner until areasonable agreement with reference data from QM/MMsimulations [24, 25] has been reached. In particular, theadsorption energy and dynamics of a single Na atom ontop of Mg and O sites on the surface of MgO were eval-uated and compared with the reference ab initio calcu- TABLE I. Parameters of the force field for Na–O and Na–Mginteractions, Eq. (1), used in the simulations.Atom pair D ij (eV) β ij (˚A − ) r ,ij (˚A) C i (eV)Na–O 0.099 1.5 2.94 0.5Na–Mg 0.001 1.32 4.85 0.5 lations [24, 25]. The best agreement in terms of bothenergetics and dynamics of a Na atom atop MgO wasachieved with the parameters summarized in Table I. Theresulting adsorption energy curves for a Na atom on topof Mg and O sites are shown in Fig. 1. A partial charge q Na = +0 . | e | was assigned to the Na atom in this case.Note that, when the Na–O and Na–Mg interactions weredescribed only with the Morse potential or when the par-tial charge on a sodium atom was set equal to zero, simu-lated MD trajectories deviated from the earlier QM/MMresults, which can be attributed to underestimation ofthe attractive forces. The calculated adsorption energiesreproduce the shape of the QM/MM results, presentingonly slightly deeper potential wells. It should be notedthat variation on the order of a few tenths of eV canarise by considering different exchange-correlation func-tionals and basis sets in quantum-chemistry calculations.We have ensured that the derived empirical potential,Eq. (1), provides a reasonable agreement with the dy-namics of a single Na atom on MgO and gives resultsconsistent with QM/MM simulations of Na and Na cluster deposition [24]. The constructed potential is thusdeemed suitable for the simulation of the deposition of ananometer-sized sodium nanoparticle. C. Sodium nanoparticle
A spherical sodium nanoparticle of 4 nm diameter, con-taining 1067 atoms, was cut out from the correspondingbulk crystal by means of MBN Studio [43]. Geometryof the nanoparticle was first optimized using the velocityquenching algorithm with a 0.1 fs time step. The many-body Gupta potential was used to describe the inter-atomic interactions with parameters taken from Ref. [45].After initial energy minimization the nanoparticle wasannealed to create a more energetically favorable start-ing geometry for the deposition simulations. Several an-nealing cycles were simulated following the proceduresreported in Ref. [46]. The first cycle consisted of heatingfrom 0 K to 400 K at a rate of 0.08 K/ps, followed by aconstant temperature simulation at 400 K for 2 ns, andcooling down to 0 K at a rate of 0.08 K/ps. The follow-upcycles were similar but the nanoparticle was heated upto 200 K, a temperature slightly below the cluster melt-ing temperature, to allow surface reorganization with-out complete melting of the nanoparticle. The nanopar-ticle melting temperature of 260 K was determined bysimulating heating of the annealed structure. The eval-uated melting temperature of the nanoparticle is closeto the experimentally determined values of 280 −
290 K
FIG. 2. MD snapshots of the Na nanoparticle depositedon a MgO substrate at several deposition energies. Panel (a)shows the initial geometry of the system. Panels (b-f) illus-trate three different deposition regimes, see text for details:(b)-(c) nanoparticle structure after 50 and 200 ps, respec-tively, for E dep = 0 . E dep = 0 .
136 eV/atom; (f) snapshot of the fragmentednanoparticle deposited at E dep = 1 .
36 eV/atom after 5 ps. [41]. Cohesive energy of the nanoparticle converged to1.0194 eV/atom after three annealing cycles. A nearlyicosahedral shape was obtained (see Fig. 2(a)), compris-ing a mix of face centered cubic and hexagonal compactstructures, in accordance with the structures known forsodium clusters [2]. The nanoparticle was equilibratedto the target temperatures of 77 K and 300 K by meansof the Langevin thermostat and the resulting velocitiescorresponded to the Maxwell-Boltzmann distribution.
D. Deposition simulations
MD simulations of the nanoparticle deposition on thesubstrate were performed in the
N V E microcanonicalensemble, thus ensuring conservation of the total energyof the system. The nanoparticle was placed in the centerof the simulation box at a 10 ˚A distance from the surface,such that initial nanoparticle–substrate interactions werenegligible (see Fig. 1). Initial atomic velocities were takenfrom pre-equilibration simulations at 77 K and 300 K inorder to simulate the deposition process at low tempera-ture at which the nanoparticle is solid, as well as at roomtemperature at which the nanoparticle has the shapeof a liquid droplet. The partial charge of +0 . | e | , de-rived in the fitting procedure described in Sect. II B, wasequally distributed among all the atoms in the nanopar-ticle. Additional velocities in the direction normal to thesurface were given to every atom of the nanoparticle suchthat the nanoparticle was deposited with kinetic energies E dep = 0 . . × .
21 nm substrate was used for most ofthe simulations. However, a larger substrate of 24 . × .
43 nm was employed for deposition energies largerthan 0.34 eV/atom to avoid interaction of the heavilydeformed or fragmented nanoparticle with its periodicimages. Two bottom MgO layers were fixed to avoid thedisplacement of the substrate upon nanoparticle impact.Simulations were performed using the leapfrog algorithmwith a 1 fs time step, which ensured that variation of thetotal energy did not exceed 0.1%. As described in Sec-tion III, most of the phenomena arising during depositionhave been observed within the first 50 ps of the simula-tions. Nonetheless, longer simulations up to 500 ps wereconducted in some cases to analyze the dynamics of thesystem on the longer time scale. III. RESULTS AND DISCUSSION
In subsection III A we briefly discuss the time evolutionof the shape and structure of the nanoparticle at differentdeposition energies. Then, in subsection III B we analyzehow fast the nanoparticle has reached thermal equilib-rium with the substrate after the deposition. In subsec-tion III C we quantitatively characterize the final shapeacquired by the nanoparticle at different deposition en-ergies and evaluate the contact angle with the substrate.Finally, the longer-term changes of the internal structureof the nanoparticle are analyzed in subsection III D.
A. Shape and structure of the depositednanoparticle
Figure 2 presents several MD snapshots of the Na nanoparticle (pre-equilibrated at 77 K) deposited onMgO at different deposition energies E dep . Three dis-tinct deposition regimes have been observed as shown inFig. 2(b-f). In the “soft” collision regime (e.g. at E dep =0 . E dep = 0 .
136 eV/atom,see panels (d) and (e)) the nanoparticle undergoes acollision-induced melting phase transition followed by itssubsequent re-crystallization. Finally, the “hard” colli-sion regime (e.g. at E dep = 1 .
36 eV/atom, panel (f))leads to rapid multifragmentation of the nanoparticle. Adetailed analysis of the nanoparticle shape is presentedbelow in Section III C.Figures 2(b-e) show snapshots from the simulationstaken at time instances of 50 ps and 200 ps. The shapeof the nanoparticle is stabilized by about 50 ps and doesnot change significantly at larger simulation times. Thenanoparticle deposited at E dep = 0 . and the formation of singlesodium atoms and small clusters was observed within thefirst 5 ps of the simulation. B. Thermalization of the nanoparticle on thesurface
Figure 3 illustrates how instantaneous temperaturesof the nanoparticle and the substrate evolved in thecourse of simulations. Results shown in panels (a) and(b) were obtained for the nanoparticle pre-equilibratedat 77 K and deposited at E dep = 0 . -2 -1 (b) I n s t an t aneou s t e m pe r a t u r e ( K ) (a) Na melting temp. Na nanoparticle MgO substrate
Time (ps)
Na nanoparticle MgO substrate Na nanoparticle MgO substrate Exponential fits (c) M a x i m u m t e m pe r a t u r e ( K ) Deposition energy (eV/atom) Na melting temp. FIG. 3. Instantaneous temperature of the Na nanopar-ticle and the MgO substrate as a function of simulation timefor deposition energies of (a) 0.0068 eV/atom and (b) 0.068eV/atom. Panel (c) shows maximum instantaneous temper-ature of the nanoparticle and the substrate as a function ofdeposition energy.
The calculated dependencies of the maximal tempera-tures T max on E dep were fitted with an exponential func-tion: T max = T + T e ω E dep , (2)where T , T and ω are the fitting parameters listed inTable II.Dashed horizontal lines in Fig. 3(b-c) illustrate themelting temperature of Na obtained from simula-tions, T m ≈
260 K. Temperature of the nanoparticle
TABLE II. Parameters for the fitting functions, Eqs. (2), (7)and (8), describing, respectively, the maximum instantaneoustemperatures of the nanoparticle and the substrate after de-position, the nanoparticle aspect ratio and the contact anglebetween the nanoparticle and the substrate at the end of sim-ulations. Na
MgO T (K) − . − . T (K) 16754.7 40842.4 ω (eV − /atom) 0.28 6 . × − δ δ − . α (eV − /atom) 16.62 θ (deg.) 68.59 θ (deg.) 23.61 γ (eV − /atom) 33.92 deposited at E dep = 0 . E dep the instantaneous temperature of the nanopar-ticle exceeded the melting temperature so that thenanoparticle experienced the melting phase transitionfollowed by re-crystallization. Deposition at energiesabove 0.68 eV/atom led to fission or multifragmentationof the nanoparticle. C. Contact angle and aspect ratio of thenanoparticle
Coordinates of sodium atoms, extracted from each sim-ulated MD trajectory, were used to parameterize theshape of the deposited nanoparticle.As follows from the simulated trajectories, the sodiumnanoparticle is, to a good approximation, radially sym-metric with respect to the main axis z . In order to char-acterize the shape of the nanoparticle, cylindrical coor-dinates ( ρ, z ) were introduced, where: ρ = q ( x − x CM ) + ( y − y CM ) , (3)with x CM and y CM being x - and y -projections of the cen-ter of mass of the nanoparticle. The ρ -axis lies in the z c oo r d i na t e ( ¯ ) (¯) dep (eV/atom): FIG. 4. Radial profile of the Na nanoparticle depositedon MgO. Symbols show the distribution of sodium atoms ofthe nanoparticle pre-equilibrated at 77 K and deposited at E dep = 0 . E dep , obtained by means of Eq. (4). MgO surface plane, whereas z -axis is perpendicular tothe surface.Figure 4 shows by symbols ( z, ρ ) projections of atomsin the Na nanoparticle deposited at 77 K with an en-ergy of 0.0034 eV/atom. The shown distribution of atomswas evaluated at the 50 ps time instance by which theshape of the deposited nanoparticle has been stabilized.The horizontal dotted line of ordinate z ≈ ρ ( z ) = q a ( z − z ) + b ( z − z ) + c , z ≥ z , (4)where a , b and c are fitting parameters. Solid coloredcurves in Fig. 4 show the inverse dependence z ( ρ ) for theNa nanoparticle equilibrated at 77 K and depositedat different energies as indicated by labels. From thisdependence the nanoparticle contact angle θ was deter-mined as [47]: θ = arctan dzdρ (cid:12)(cid:12)(cid:12)(cid:12) ρ = ρ ′ ! , (5)where ρ ′ is the point of intersection between the fittingcurve and the line z = z representing the bottom-mostNa atomic layer parallel to the substrate surface. N anopa r t i c l e a s pe c t r a t i o Deposition energy (eV/atom): 0.0034 0.0068 0.0136 0.034 0.068 0.102 0.136 0.34 (b)(a)
Time (ps)
FIG. 5. Nanoparticle aspect ratio δ , Eq. (6), as a function ofsimulation time for deposition at (a) 77 K and (b) 300 K. Solid curves shown in Fig. 4 illustrate a gradual changeof the nanoparticle shape from a truncated prolatespheroid ( z max > ρ max ) to a semi-spheroid ( z max = ρ max )to a truncated oblate spheroid ( z max < ρ max ) as thedeposition energy increases. The dashed orange curvein Fig. 4 shows the z ( ρ ) dependence for the nanopar-ticle equilibrated at 300 K and deposited at E dep =0 . nanoparticle equi-librated at 77 K and deposited at 0.34 eV/atom (lightgreen curve) is very similar to that of the nanoparticleequilibrated at 300 K and deposited at much lower energyof 0.0034 eV/atom (dashed orange curve). In the for-mer case the nanoparticle experienced a collision-inducedmelting phase transition while in the latter case it wasdirectly deposited as a liquid droplet. The equilibriumshape of the nanoparticle deposited at 300 K dependsvery little on deposition energy. As a result, the shapes ofthe nanoparticles deposited at higher values of E dep (be-low the threshold energy for nanoparticle fragmentation)are similar to the one shown for E dep = 0 . E dep is discussed in greater detail below.We have also analyzed the ratio of the diameter of the bottom-most atomic layer of the nanoparticle (see reddotted line in Fig. 4) to the nanoparticle height, i.e. thenanoparticle aspect ratio [16]. It is defined as: δ = 2 ρ ( z ) z ( ρ = 0) . (6)The temporal evolution of δ for the nanoparticle de-posited at different energies is illustrated in Fig. 5. Pan-els (a) and (b) illustrate the deposition at 77 K and300 K, respectively. We consider here the values of E dep up to 0.34 eV/atom for 77 K and up to 0.136 eV/atom for300K, energies at which the nanoparticle have not frag-mented upon collision with the substrate. Values of δ be-tween 1 and 2 indicate that the nanoparticle acquires theshape of a truncated prolate spheroid, whereas δ = 2 de-scribes a perfect semi-spheroid. A parallel can be drawnwith the snapshots of nanoparticle deposition at 77 K,shown in Fig. 2(b-e). As discussed above, the overallshape of the nanoparticle practically does not change af-ter 50 ps of simulation. This observation is confirmedby constant values of δ observed after 50 ps. As shownin Fig. 5 this trend has been observed for all depositionenergies within the range considered.Two distinct deformation regimes have been revealedin the simulations. At high deposition energies ( E dep ≥ .
068 eV/atom), the nanoparticle aspect ratio increasesover the first few picoseconds and eventually convergesto the value of δ ≈ E dep . In contrast, a gradual increase of the final value of δ from about 1 to 2 is typical for low-energy depositionsat energy E dep ≤ .
034 eV/atom.As discussed above, the nanoparticle deposited in the“soft” regime remains in the solid phase in the course ofdeposition; this regime corresponds to a small increaseof δ over time. When temperature of the nanoparticleexceeds its melting temperature, the nanoparticle expe-riences the melting phase transition and thus it becomesmore susceptible to stress-induced deformation. This oc-curs at more energetic collisions when the nanoparticlewets the substrate and acquires the shape of an oblatespheroid, see Fig. 2(d-e). On this basis it is straight-forward to explain the dependencies shown in Fig. 5(b),which describe the aspect ratio for the nanoparticle pre-equilibrated at 300 K. A liquid sodium drop depositedon MgO experiences strong deformation upon contactingthe substrate; this behavior is characterized by a rapidincrease of δ (i.e. flattening of the nanoparticle) withinthe first 5 −
10 ps of the simulations. At larger time in-stances δ saturates at a constant value of about 3 whenthe nanoparticle reaches thermal equilibrium with thesubstrate. Note that the nanoparticles deposited at dif-ferent values of E dep have a similar shape as the aspectratio asymptotically approaches a constant value δ ≈
77 K 300 K Exponential fit (b) N anopa r t i c l e a s pe c t r a t i o (a)
77 K 300 K Exponential fit C on t a c t ang l e ( deg r ee s ) Deposition energy (eV/atom)
FIG. 6. The final aspect ratio (a) and contact angle (b) ofthe nanoparticle as functions of the deposition energy. Sym-bols denote time averages, while error bars represent standarddeviations. Data for the deposition at 77 K were fitted withexponential functions, Eqs. (7) and (8). in a broad deposition energy range below the nanopar-ticle fragmentation threshold. For each E dep , averagevalues and standard deviations in time were calculatedfor the parts of the trajectories at which δ and θ oscillatearound a constant value. Symbols represent the meanvalues while error bars correspond to a standard devia-tion.Two different trends are clearly seen for different de-position temperatures. At 300 K (black squares) boththe aspect ratio and the contact angle are practicallyconstant within the studied energy range. Both param-eters fluctuate around some characteristic values, δ ≈ θ ≈ ◦ , which are indicated by horizontal dottedlines. For deposition at 77 K (blue circles), aspect ratiogrows exponentially with E dep while the contact angledecreases. The dependencies of δ and θ on E dep werefitted with the following functions: δ = δ + δ e α E dep , (7) θ = θ + θ e − γ E dep . (8)Parameters of this fit are listed in Table II. At deposi-tion energies larger than 0.068 eV/atom both aspect ratioand the contact angle for the nanoparticle deposited at 77 K approach the values obtained at 300 K as a conse-quence of the collision-induced melting phase transitionabove this energy. These results agree with the conclu-sions made above that the resulting shape of a liquiddroplet deposited at 300 K depends very little on the de-position energy. In contrast, structure and contact anglefor the nanoparticle deposited at the lower temperaturestrongly depend on E dep . D. Collision-induced structural and phasetransformations
Further insights into the change of internal structure ofthe deposited nanoparticle can be drawn from the anal-ysis of the radial distribution function (RDF). Figure 7shows RDFs for the nanoparticle equilibrated at 77 Kand deposited at E dep = 0 . E dep = 0 .
068 eV/atom (panel (b)). The plottedRDFs were averaged over different periods of time, whichare marked by dashed vertical lines in Fig. 3(a,b). TheRDFs for the deposited Na nanoparticle are com-pared with the experimentally determined distributionfor liquid sodium [49].As shown in Fig. 7(a) main peaks in the RDFs remainfrom the initial stage of deposition (within the first 3 psof simulation) up to the last stage when the nanoparti-cle has reached thermal equilibrium with the substrate( >
420 ps). Moreover, at all stages of the simulated tra-jectory the calculated RDFs differ from the one for liquidsodium (see open symbols). The melting phase transitionand subsequent re-crystallization of the nanoparticle isclearly seen by the variation of RDFs in Fig. 7(b). Tem-perature of the nanoparticle deposited at 0.068 eV/atomreaches the maximum value at about 2.5 ps and it re-mains higher than the melting temperature up to about40 ps, see Fig. 3(b). In between these time instancesthe nanoparticle transforms into a liquid droplet as con-firmed by the close similarity of the calculated RDF withthe one for liquid sodium [49]. At larger time instances(the region of 40 −
300 ps) the nanoparticle re-crystallizesand eventually reaches thermal equilibrium with the sub-strate at about 400 ps. Re-crystallization manifests itselfin reappearance of peaks in the RDF, which resemblesthe RDF in the very beginning of the simulation.
IV. CONCLUSION
In this paper deposition of a nanometer-sized sodiumnanoparticle, containing 1067 atoms, on a MgO substratewas studied by means of molecular dynamics simulationsusing the MBN Explorer and MBN Studio software pack-ages. We focused on the broad deposition energy rangeof 0 . − .
36 eV/atom, which covers both the soft-landing and nanoparticle fragmentation regimes. Sim-ulations were performed at two different temperatures,77 K and 300 K, at which the nanoparticle is either in (b) R ad i a l d i s t r i bu t i on f un c t i on (a) Interatomic distance (¯)
FIG. 7. Radial distribution function (RDF) for the Na nanoparticle deposited at (a) 0.0068 eV/atom and (b) 0.068eV/atom. Different curves describe the RDF evaluated overdifferent time periods as indicated. These time periods aremarked by dash-dotted vertical lines in Fig. 3(a,b). the solid state or it forms a liquid droplet.A force field describing the interaction of sodium atomswith the MgO surface was developed and used in thesimulations. The force field was validated via calculat-ing adsorption energies of a single Na atom on differentMgO sites. The presented results agreed with the re-sults obtained previously in ab initio
QM/MM simula- tions [24, 25, 27].The process of nanoparticle deposition and subsequentrelaxation on the surface has been studied in detail as afunction of deposition energy. In particular, variation ofthe nanoparticle shape as a function of simulation time,its wetting properties, as well as the energy transfer be-tween the nanoparticle and the substrate were analyzedon the timescale of up to several hundreds picoseconds.This study provides detailed insights into the dynam-ics of sodium nanoparticle deposition on MgO substrateswhich complement the information already gathered forsmall clusters from QM/MM simulations and which maybe useful for experimental studies. A similar analysisfor other metallic aggregates deposited onto experimen-tally relevant oxide surfaces might reveal atomistic-levelinsights into the structure and shape of the depositedmetal systems, which might be useful for technologicalapplications.
ACKNOWLEDGEMENTS
This work was supported in part by DeutscheForschungsgemeinschaft (Project no. 415716638); by theEuropean Union’s Horizon 2020 research and innovationprogramme – the Radio-NP project (GA 794733) withinthe H2020-MSCA-IF-2017 call and the RADON project(GA 872494) within the H2020-MSCA-RISE-2019 call;by the Spanish Ministerio de Ciencia e Innovaci´on andthe European Regional Development Fund (Project no.PGC2018-096788-B-I00); by the Fundaci´on S´eneca –Agencia de Ciencia y Tecnolog´ıa de la Regi´on de Murcia(Project No. 19907/GERM/15); and by the Conselleriad’Educaci´o, Investigaci´o, Cultura i Esport de la Gener-alitat Valenciana (Project no. AICO/2019/070). PdVgratefully acknowledges the Alexander von HumboldtFoundation/Stiftung and the Spanish Ministerio de Cien-cia e Innovaci´on for their financial support by means of,respectively, Humboldt (1197139) and Juan de la Cierva(FJCI-2017-32233) postdoctoral fellowships. The possi-bility to perform computer simulations at Goethe-HLRcluster of the Frankfurt Center for Scientific Computingand the Scientific Computing Service of the University ofMurcia is gratefully acknowledged. [1] W. A. de Heer, Rev. Mod. Phys. , 611 (1993).[2] H. Haberland, ed., Clusters of Atoms and Molecules.Theory, Experiment, and Clusters of Atoms (Springer-Verlag, Berlin Heidelberg, 1994).[3] U. Kreibig and M. Vollmer,
Optical Properties of MetalClusters (Springer-Verlag, Berlin Heidelberg, 1995).[4] K.-H. Meiwes-Broer, ed.,
Metal Clusters at Surfaces.Structure, Quantum Properties, Physical Chemistry (Springer-Verlag, Berlin Heidelberg, 2000).[5] J.-P. Connerade and A. Solov’yov, eds.,
Latest Advancesin Atomic Cluster Collisions: Structure and Dynamics from the Nuclear to the Biological Scale (Imperial CollegePress, London, 2008).[6] P. M. Dinh, P.-G. Reinhard, and E. Suraud,
An Intro-duction to Cluster Science (Wiley-VCH, Berlin, 2013).[7] F. Calvo, Phys. Chem. Chem. Phys. , 27922 (2015).[8] S. T. Bromley and S. M. Woodley, eds., ComputationalModelling of Nanoparticles, Frontiers of Nanoscience vol.12 (Elsevier, 2018).[9] D. B. Janes, M. Batistuta, S. Datta, M. R. Mel-loch, R. P. Andres, J. Liu, N. P. Chen, T. Lee,R. Reifenberger, E. H. Chen, and J. M. Woodall, Superlattices Microstruct. , 555 (2000).[10] L. Sun, H. Lin, K. L. Kohlstedt, G. C. Schatz, and C. A.Mirkin, Proc. Natl. Acad. Sci. U.S.A. , 7242 (2018).[11] U. Heiz and U. Landman, eds., Nanocatalysis (Springer-Verlag, Berlin Heidelberg, 2007).[12] X. Han, K. Xu, O. Taratula, and K. Farsad, Nanoscale , 799 (2019).[13] K. Haume, S. Rosa, S. Grellet, M. A. ´Smia lek, K. T.Butterworth, A. V. Solov’yov, K. M. Prise, J. Golding,and N. J. Mason, Cancer Nanotechnol. , 8 (2016).[14] A. V. Verkhovtsev, A. V. Korol, and A. V. Solov’yov,Phys. Rev. Lett. , 063401 (2015).[15] Y. Liu, P. Zhang, F. Li, X. Jin, J. Li, W. Chen, andQ. Li, Theranostics , 1824 (2018).[16] V. N. Popok, I. Barke, E. E. Campbell, and K. H.Meiwes-Broer, Surf. Sci. Rep. , 347 (2011).[17] P. Milani and M. Sowwan, eds., Cluster Beam Depositionof Functional Nanomaterials and Devices, Frontiers ofNanoscience vol. 15 (Elsevier, 2020).[18] D. N. Poenaru, R. A. Gherghescu, A. V. Solov’yov, andW. Greiner, EPL (Europhys Lett.) , 63001 (2007).[19] D. N. Poenaru, R. A. Gherghescu, I. H. Plonski, A. V.Solov’yov, and W. Greiner, Eur. Phys. J. D , 379(2008).[20] H. H¨akkinen, Adv. Phys. X , 467 (2016).[21] H. H¨ovel, B. Grimm, M. Pollmann, and B. Reihl, Phys.Rev. Lett. , 4608 (1998).[22] A. V. Verkhovtsev, Y. Erofeev, and A. V. Solov’yov,Eur. Phys. J. D , 205 (2020).[23] H. H¨akkinen and M. Manninen,J. Chem. Phys. , 10565 (1996).[24] M. B¨ar, Non-Linear Dynamics of Metal Clusters onInsulating Substrates , PhD Thesis, Friedrich-Alexander-Universit¨at Erlangen-N¨urnberg (2008).[25] P. M. Dinh, P. G. Reinhard, and E. Suraud, Phys. Rep. , 43 (2010).[26] M. B¨ar, L. V. Moskaleva, M. Winkler, P. G. Reinhard,N. R¨osch, and E. Suraud, Eur. Phys. J. D , 507 (2007).[27] M. B¨ar, P. M. Dinh, L. V. Moskaleva, P. G. Reinhard,N. R¨osch, and E. Suraud, Phys. Rev. B , 195404(2009).[28] C. R. Henry, Surf. Sci. Rep. , 231 (1998).[29] G. Barcaro, A. Fortunelli, F. Nita, and R. Ferrando,Phys. Rev. Lett. , 246103 (2005).[30] G. Barcaro and A. Fortunelli, J. Chem. Theory Comput. , 972 (2005). [31] D. Paradiso and J. Z. Larese, J. Phys. Chem. C ,14564 (2020).[32] F. Buend´ıa, A. T. Anzaldo, C. Vital, and M. R. Beltr´an,J. Chem. Phys. , 024303 (2020).[33] T. T. Magkoev, G. S. Grigorkina, V. B. Zaalishvili, O. G.Burdzieva, E. N. Kozyrev, G. E. Tuaev, and K. Fuku-tani, Rus. J. Phys. Chem. A , 401 (2020).[34] S. Nigam and C. Majumder, Appl. Surf. Sci. , 144963(2020).[35] R. S. Averback and M. Ghaly,Nucl. Instrum. Meth. B , 191 (1994).[36] S. J. Carroll, P. D. Nellist, R. E. Palmer, S. Hobday, andR. Smith, Phys. Rev. Lett. , 2654 (2000).[37] T. J. Colla, R. Aderjan, R. Kissel, and H. M. Urbassek,Phys. Rev. B , 8487 (2000).[38] P. Me, A. Perez, Q. Hou, M. Hou, L. Bardotti, andB. Pre, Phys. Rev. B , 2835 (2000).[39] M. Moseler, O. Rattunde, J. Nordiek, and H. Haberland,Nucl. Instrum. Meth. B , 522 (2000).[40] I. A. Solov’yov, A. V. Korol, and A. V. Solov’yov, Multi-scale Modelling of Complex Molecular Structure and Dy-namics with MBN Explorer (Springer International Pub-lishing AG, 2017).[41] T. P. Martin, U. N¨aher, H. Schaber, and U. Zimmer-mann, J. Chem. Phys. , 2322 (1994).[42] I. A. Solov’yov, A. V. Yakubovich, P. V.Nikolaev, I. Volkovets, and A. V. Solov’yov,J. Comput. Chem. , 2412 (2012).[43] G. B. Sushko, I. A. Solov’yov, and A. V. Solov’yov,J. Mol. Graph. Model. , 247 (2019).[44] A. M. Dongare, D. D. Hass, and L. V. Zhigilei, Clustersand Nano-Assemblies: Physical and Biological Systems:Richmond, Virginia, U.S.A., 10-13 November, 2003 , 329(2005).[45] Y. Li, E. Blaisten-Barojas, and D. A. Papaconstantopou-los, Chem. Phys. Lett. , 331 (1997).[46] T. Ellaby, J. Aarons, A. Varambhia, L. Jones, P. Nel-list, D. Ozkaya, M. Sarwar, D. Thompsett, and C.-K.Skylaris, J. Phys.: Condens. Matter , 155301 (2018).[47] N. Giovambattista, P. G. Debenedetti, and P. J. Rossky,J. Phys. Chem. B , 9581 (2007).[48] J. ˇSkv´ara, J. ˇSkvor, and I. Nezbeda, Molec. Simul. ,190 (2018).[49] R. D. Murphy and M. L. Klein,Phys. Rev. A8