(n = 16–24) Clusters: A Case Study - ACS Publications - American

Article pubs.acs.org/JPCC

Multifaceted Thermodynamics of Pbn (n = 16−24) Clusters: A Case Study Anju Susan,†,‡ Vaibhav Kaware,†,§ and Kavita Joshi*,†,‡ †

Physical and Materials Chemistry Division, CSIR-National Chemical Laboratory, Pune 411008, India Academy of Scientific and Innovative Research (AcSIR), Anusandhan Bhawan, 2, Rafi Marg, New Delhi 110 001, India § Department of Physics, Savitribai Phule Pune University, Ganeshkhind, Pune 411007, India ‡

ABSTRACT: Thermodynamic response of small clusters is a challenging area of exploration, both experimentally and theoretically. In this article, we study the thermodynamic behavior of small Pb clusters (size 16−24) using Born− Oppenheimer molecular dynamics. A new ground state structure is reported for Pb20. Except for Pb21, all clusters fragment at temperatures above Tm[bulk] and show no signs of melting. Characteristic behavior like restricted diffusion and solid−solid transition is discussed in detail. Variation in the isomerization temperature of these clusters is explained using the bond length analysis. Root mean square bond length fluctuations (δrms) along with distribution of atoms about center of mass of the cluster as a function of time and distance−energy (DE) plots are used to bring out the essential features of Pb cluster thermodynamics. Analysis carried out using these parameters, and their interpretation regarding state of the system, are discussed in detail. We highlight that it is not possible to define “liquid state” for these small clusters, in the conventional frame of understanding.

■

INTRODUCTION Thermodynamics of clusters is an area that has few laid out rules and sorted out theories.1−20 Although much work has been done, newer explanations are needed as newer facts come into light with more work done progressively. Basic concepts like “melting” present problems of interpretation when applied to clusters directly. A comparison of melting of clusters with that of solids brings forth a few prominent differences. First, although almost all solids, with the exception of those that sublimate, exhibit melting, the same cannot be said about clusters. Although melting is a first-order transition for solids, the finite size of small clusters renders a finite range of temperatures over which the clusters melt. Though it takes conscious and considerable effort to maintain a solid at the boundary of solid and liquid phase, there exist extended ranges of temperatures over which clusters show a dynamical solid− liquid phase existence.6,7 Further, some of the small sized clusters of tin and gallium melt at higher temperatures than their bulk counterpart.18,19,21,22 A few clusters of Ga are known to melt without any latent heat, which is interpreted as a second-order transition.23 Although some solids exhibit the unusual phenomenon like the solid−solid phase transition,24 many clusters, along with solid−solid phase transition exhibit other unusual thermodynamic processes like partial melting, phase coexistence, fragmentation, and half-solidity.3,11−17 Atomic clusters can be divided into two classes on the basis of their size. The class of comparatively larger clusters, containing more than few hundred atoms follows the 1/r rule with their melting temperature decreasing with size.1,5,8 © 2015 American Chemical Society

However, clusters that consist of only a couple of hundred, or less number of atoms, show an erratic melting behavior, with no such simple rule as 1/r. Experimental studies for small clusters of Na, Ga, Al, and Sn have proven this fact.9,12,18−20 When all atoms in the cluster are on the surface, bonding among individual atoms is the most dominant contributor to the cluster’s stability and its finite temperature behavior is highly influenced by even small changes in positions of atoms or in the number of electrons.15,25−28 In slightly larger clusters, it is possible to make a distinction of internal and surface atoms. In these clusters, the response of surface and internal atoms toward external temperature is different.3,12,14,29 For still larger clusters, changes in positions of individual atoms lead to minor alterations in their properties, and variation in properties subsides in amplitude henceforth.20 Thus, we note that small clusters are a complex system that cannot be modeled without taking into account above considerations. One cannot ignore the trends in global minimum energy structures (GS) of small clusters while explaining their thermodynamics, because their GS structures sit at the center stage, when it comes to explaining their finite temperature behavior.30−33 It has been demonstrated for gallium as well as aluminum clusters that variation in the structural motif of the GS and melting temperatures are correlated.30,33 It must be noted, however, that although it is not possible to explain the Received: June 2, 2015 Revised: September 22, 2015 Published: September 22, 2015 23698

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707

Article

The Journal of Physical Chemistry C Table 1. Details about the Simulation Dataa cluster size no. of temperatures total simulation time (ns) Tfrag (K)

16 17 5.37 900

17 16 4.26 850

18 15 4.37 750

19 16 4.38 650

20 15 3.45 750

21 15 3.76 550

22 14 5.17 700

23 15 5.03 750

24 15 5.22 750

The first row is the size of the system. Total number of temperatures for each system with corresponding total simulation time (in ns) are in the second and third rows, respectively. The last row contains their fragmentation temperatures (Tfrag). a

simultaneously (i) heated the system at the rate of 50 K/15 ps to the next temperature and (ii) maintained the system at that temperature for at least 210 ps. Because the time needed for thermalization is more in the transition region, the system was maintained for 300 ps or more in this region. We observed that all the clusters undergo fragmentation at higher temperatures, instead of melting. Hence, the simulations were performed for a number of temperatures ranging from 14 to 17, in the temperature range 100 K, to the fragmentation temperature (Tfrag). The resulting total simulation time varied from 3.4 to 5.37 ns for different sizes. Details about the number of temperatures and the total simulation time for each system, along with their fragmentation temperatures, is presented in Table 1. The finite temperature data were analyzed using rootmean-square bond length fluctuations (δrms), distribution of atoms with respect to their distance from center of mass (DCOM) of the cluster, mean square displacements (MSD), distance−energy (DE) plots, and heat capacity curves computed using multiple histogram technique (not discussed here). Local optimizations were carried out using the geometries picked from the finite temperature molecular dynamics (MD) data. These initial geometries were selected in such a way that they are unbiased and equally spaced (each geometry after a period of 3 ps), which resulted into at least 1000 initial geometries per system. An extensive isomer analysis was carried out using these optimized structures. For local optimization, the convergence criteria for force on each ion was taken as 0.005 eV/Å. Vibrational frequencies were calculated for GS of all the sizes, to confirm that the obtained geometries are indeed a local minima. Root Mean Square Bond Length Fluctuations. Root mean square bond length fluctuations (δrms) is one of the prominent theoretical tools used to characterize thermodynamic behavior of solids as well as clusters.4 It is defined as

liquid-like state of clusters from their GS structures, their behavior until the point of isomerization and before melting is closely related to their GS for clusters with few tens of atoms.32 Another example could be elements Si and Sn, belonging to group IV-A, which exhibit fragmentation of their clusters with increasing temperature.15,16,19,22 A possible reason for such behavior could be the covalent type bonding among atoms of these elements, which leads to formation of locally bonded, more stable smaller subunits at finite temperatures. These subunits are bonded with less strong intersubunit bonding, which could lead to fragmentation.16,34,35 Previous investigations about Pb clusters in this size range concern their geometries, bonding, electric dipole moments, and electronic structures.35−43 However, with the exception of work done by Pushpa et al., there exist no first principle studies that discuss the thermodynamics of small clusters of Pb.44 In their work, Pushpa et al. relate the stiffness of bonds in small (j

(⟨rij2⟩t − ⟨rij⟩t2 )1/2 ⟨rij⟩t

where N is the total number of atoms in the system, rij is the distance between particle i and j, and ⟨ ⟩t denotes time average over the entire trajectory. δrms is a measure of fluctuations in bond lengths averaged over time t that occur in a system at a specific temperature. A schematic depicting typical behavior of δrms for atomic clusters, is shown in Figure 1. As shown in the figure, the regions of solid (δrms < 0.05), liquid (δrms > 0.2), and the jump from ∼0.05 to ∼0.2 can be distinctly identified for a cluster that exhibits a well-defined melting behavior. The width of the transition region defines if the cluster is a “magic” melter, or broad melter. Because an increase in temperature expands the system, one expects an increase in average bond lengths as the temperature increases. A cluster in its solid state basically retains its GS structure, with atoms oscillating about their mean positions. A change in the slope of δrms occurs due to breaking 23699

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707

Article

The Journal of Physical Chemistry C

note that the behaviors of δrms are almost identical whether we remove the initial 30 ps/temp (total simulation time 4.86 ns) or 60 ps/temp (total simulation time 4.35 ns) of data to calculate it. This implies that the conclusions drawn from the analysis of these data will not change substantially upon extending simulations further. We also note that, the simulation time required for converged results varies from 200 to 500 ps depending upon the temperature. The convergence for each temperature along with error bars is shown in Figure 2b. The error bars have been computed by dividing the complete data in two sets, and calculating δrms for each of the set separately. It is worth noting that, for temperatures associated with structural transitions or fragmentation, the error is relatively larger.

■

Figure 1. A schematic of the simple δrms variation as seen in atomic cluster that exhibits melting. δrms has a constant slope in the solid-like region and in the liquid-like region. A steep increase in the value of δrms occurs during the transition region, indicating the isomerization. The finite width of temperature over which the transition occurs implies that the clusters melt over a range of temperatures gradually, and not at single melting temperature like in solids. The dotted red line indicates the schema when the melting transition is narrow, and the bold black line is the schema when the cluster melts over a broad range of temperatures.

RESULTS AND DISCUSSION Ground State Structures. We begin our discussion by presenting the putative global minima for Pb clusters in the size range 16−24, shown in Figure 3. We began our search for the GS from the previously published data.35,38−41 Nearly 3−5 ns of BOMD was performed for each size, the details of which are noted in Table 1. About a thousand initial positions were optimized from the MD for each size, before concluding the search for their GS. The revised, or otherwise, GS structures from this search were used as the starting point of the thermodynamics runs. With such a large data set, the GS structures used to initiate the thermodynamics become most reliable to explain the structure dependent finite temperature properties of these small clusters. In this section, we discuss the salient features of these structures that are interesting by themselves and will be used later to explain their thermodynamics. A new GS structure was obtained for size 20 in our search that has not been reported in the literature. Although the earlier reported structure of Pb20 is composed of a single added atom at the center of Pb19 GS (shown in Figure 3), the newly found structure is without an internal atom, and its energy is 0.17 eV less than that of the previously reported GS. This makes Pb20 the only other structure that is without an internal atom, besides Pb16. Pb20 is also the only other GS that is not prolate, along with Pb17. The other nonprolate geometry, Pb17, has the most symmetric GS among all other sizes. All the geometries except Pb17 and Pb20 are prolate, and formed out of stacking of rings of Pb atoms. Number of atoms in the rings vary as the size changes, and the three stacks of rings grow to four rings and an apex atom as the cluster grows in size from 16 to 24 atoms. The more stable Pb6 and Pb7 units are seen in almost all of these structures.

of bonds, when the oscillations of atoms extend beyond their own mean positions, and the cluster begins isomerizing. In the liquid-like state it regains the linear behavior similar to that of the solid-like state but with a different slope. Convergence. Before discussing the convergence of these simulations we note that these are ab initio MD simulations and time scales of these simulations could not be compared with that of classical MD simulations that are not computationally demanding. The accuracy of any statistical data depends critically on its convergence. Hence, any MD should be tested for its convergence as well. As a first measure in this regard, we always discarded first 30 ps (for each temperature) of MD for thermalization. To demonstrate the degree of convergence in our calculations, we have plotted δrms computed over three different subsets of our data, as shown in Figure 2a. Black points in the graph represent the “nonconverged” subset of the data, which includes δrms calculated over the first 30 ps/temp (in addition to leaving initial 30 ps for thermalization) of MD. The other two graphs, red and blue circles, respectively, plot δrms calculated over data subsets with all data points except the initial 30 ps/temp (in addition to 30 ps/temp left for thermalization), and all data points except the initial 60 ps/ temp (in addition to 30 ps/temp left for thermalization). We

Figure 2. (a) δrms calculated over subsets of the entire data. One of the graphs shows the nonconverged results (black circles in the graph) to contrast it with other two larger subsets of data that show converged values of δrms. (b) Graph of δrms for the entire data, in comparison with δrms of the initial half of the data (lower point of error bar), and that of the later half of the data (upper point of error bar). 23700

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707

Article

The Journal of Physical Chemistry C

Jellium trend. We find that MOs formed out of atomic s orbitals are the only ones that follow the Jellium model (in terms of degeneracy as well as shape of MOs). The rest of the MOs mostly form the p-complexes in accordance with the symmetry of the cluster. Because there are two valence s electrons per Pb atom, only n orbitals can be occupied by these 2n s electrons in the cluster Pbn. Table 2 gives the detailed statistics on the extent to which the MOs of Pb clusters follow Jellium orbitals. We find that only for Pb17, a maximum of all 17 MOs (formed out of atomic s orbitals) followed the Jellium orbitals. For the rest of them, less than n (between 50% to 90%) MOs (formed out of atomic s orbitals) followed Jellium, in an n atom cluster. It is interesting to note that such a partial/ characteristic similarity of atomic clusters with Jellium has been previously seen in clusters of gallium and of aluminum.32,55 Thus, similarity of “MOs formed out of atomic s” with those of the Jellium orbitals is not just Pb specific but seems to be a trait common to other small clusters in general. Next, we plot the fraction of atoms with their coordination number against the size of the cluster in Figure 4. The plots are

Figure 3. Ground state structures of Pbn (n = 16−24) clusters. Sizes of clusters are mentioned at the top of each structure, along with their symmetry groups in the parentheses. The internal atoms are shown in pink. In the case of Pb20, the four atoms shown in cyan color are additions to the Pb16 motif.

Figure 4. Numbers at the top are the coordination numbers plotted as corresponding colored vertical bars against the fraction of atoms with that coordination for all the clusters.

constructed with a cutoff bond length of 3.8 Å. The cutoff is chosen on the basis of the distribution of bond lengths, which shows a gap at value of 3.8 Å, indicating the bifurcation of first nearest neighbors from the further nearest neighbors. The number at the top of each histogram indicates the coordination and on the y-axis we have plotted fraction of atoms with a specific coordination. For example, in the case of Pb17, the same fraction of atoms (8 out of 17) has 5- and 6-coordination. Further, we also observe that the coordination of these clusters shows an increasing trend in general (Figure 4). A lower coordination of 4 occurs only in Pb16, after which the fraction of atoms with coordination 5 and 6 goes on increasing with the cluster size. Because metals have a higher coordination, it can be concluded that Pb clusters become relatively “more” metallic with increasing size. Finite Temperature Behavior. Having discussed the important static properties of the Pb clusters, we now proceed to discuss the finite temperature behavior of these clusters. What follows in the rest of the results and discussion is a detailed account of thermodynamic behavior of these clusters, discerned using different theoretical quantities like the δrms,

Thus, in summary, all the studied clusters have a prolate structural motif, except Pb17 and Pb20, with Pb17 being the most symmetric structure among all. Pb16 and Pb20 are the only ones without any explicit central atom. This specific feature is the cause of the solid−solid transition observed in these two sizes during their thermodynamics, which will be discussed later. Bonding. Bonding in small Pb clusters has been an issue of debate in the past35,38,39 and is not yet resolved. Toward this end, we have analyzed molecular orbitals (MOs), coordination of atoms within the cluster, and electron localization function (ELF) to gain insight into bonding of Pb clusters. First, we bring out some interesting facts about molecular orbitals of Pb clusters. We note that two factors are considered while concluding whether a cluster follows the Jellium model or not. The first one is the observed degeneracy in the eigenvalue spectra and its resemblance with that of Jellium spectra, and the second is the shape of the MOs and their resemblance with the characteristic shapes of Jellium MOs. From our analysis, we conclude that the clusters follow the Jellium model, only if the degeneracy pattern, as well as the shape of MOs, follow the

Table 2. Molecular Orbitals of Pb Clusters Compared with Jellium MOsa size s MOs† p MOs‡ a

16 12/16 (75%) 0

17 17/17 (100%) 3

18 13/18 (72%) 1

19 12/19 (63%) 0

20 15/20 (75%) 3

21 16/21 (76%) 0

22 14/22 (63%) 2

23 13/23 (57%) 1

24 10/24 (42%) 3

s MOs†: MOs formed out of atomic s, and follow Jellium orbitals. p MOs‡: MOs formed out of atomic p, and follow Jellium orbitals. 23701

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707

Article

The Journal of Physical Chemistry C

have the least fragmentation energy. Thus, Pb16 has two possible channels of fragmentation, (7,9) and (8,8). However, (7,9) is the most probable channel observed in MD simulations at the fragmentation temperature and indeed has much less fragmentation energy compared to the other possible fragment (8,8). However, when the difference in the fragmentation energy is less than 0.2 eV, the corresponding fragments are observed with similar probability. For example, in the case of Pb23, fragmentation channels (8,15), (11,12), and (9,14) are observed with same probability during the MD. Further, we note that fragments of size 7, 9, or 10 are ubiquitously present for all the clusters. It has been found experimentally that Pb clusters of sizes 7, 9, and 10 are the most abundant in their mass spectra, which implies that they are the most stable sizes of Pb clusters.56 Root Mean Square Bond Length Fluctuations (δrms). δrms is a popular tool to describe “melting” of bulk solids, which has been adapted for finite sized systems as well. Figure 7 plots the δrms for Pbn (n = 16−24) clusters, as a function of temperature. We observe that δrms plots for these clusters look quite atypical when compared with the schematic shown in Figure 1. We also note that the solid state of clusters Pb17, Pb18, and Pb20 exists over an extended range of temperatures. Other clusters either do not exhibit any pronounced solid state, like Pb19, Pb23, and Pb24, or have a very brief one, like Pb16, Pb21, and Pb22. A oneto-one correspondence of δrms of Pb clusters with the schematic behavior (Figure 1) of δrms ends here, because the response of all of these clusters to increasing temperature is entirely different thereafter. The thermodynamic behavior of these clusters can be said to follow the general trend as

radial distribution function about center of mass of the cluster, detailed analysis of MD trajectories, and isomer analysis using distance energy plots. Fragmentation. All Pb clusters between sizes 16 and 24 undergo fragmentation at elevated temperatures. Figure 5

Figure 5. Fragmentation temperatures (Tfrag) of Pbn (n = 16−24) clusters. Tm[bulk] is the melting temperature of bulk Pb. Fragmentation temperatures for all sizes except Pb21 are above that of melting temperature of bulk Pb (600 K).

shows the fragmentation temperatures of these clusters, as a function of their size. We note that all these Pb clusters, except Pb21, fragment above the melting temperature of bulk Pb (Tm[bulk]). This agrees well with the prediction from previous work that Pb clusters will be stable above Tm[bulk].44 The fragmentation temperatures vary in the range of about 350 K, from 550 K for Pb21, to 900 K for Pb16. We have calculated the fragmentation energy (Efrag) for each of the clusters for different fragmentation channels observed during MD. For a cluster of size n that fragments into two smaller clusters of sizes m1 and m2, the fragmentation energy is defined as

GS geometry/solid state → marginally modified GS 2−4 atoms displace → major rearrangement of atoms ∼50% or more atoms displace → fragmentation

Efrag = (EmGS1 + EmGS2 ) − EnGS

(1)

Each of these transitions (→) is signified by a jump in δrms. Except for Pb16, Pb18, and Pb20, all clusters follow this trend faithfully. Although it is desirable to quantify terms like “major” and “minor” rearrangements, one can only get as close to their quantification as describing them more discretely. In this sense, by “marginally modified” structures, we mean a class of isomers with a structural motif similar to that of GS, but with minor

GS where EGS m1 and Em2 are the energies of the GS structures of the fragments of the cluster, EGS n is the energy of the GS of cluster of size n, and n = m1 + m2. A plot of this quantity for all Pb clusters between sizes 16 and 24 is shown in Figure 6. Multiple data bars for the cluster implies that there is more than one possible way in which the particular cluster has fragmented during the MD. The most probable fragmentation channels

Figure 6. Fragmentation energies of Pbn (n = 16−24) clusters. Numbers in the brackets correspond to the sizes of the fragments. Multiple data bars imply that cluster has fragmented in more than one ways during the MD. 23702

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707

Article

The Journal of Physical Chemistry C

Figure 7. Variation in δrms of Pb clusters as a function of their temperature. Red circles, marked Tiso, indicate “high” isomerization temperatures, whereas golden circles indicate the “low” isomerization temperatures. Magenta rectangles (sizes 16 and 18) highlight the temperature ranges in which “restricted diffusion” occurs. Figures highlighting the motion of atoms during restricted diffusion are shown alongside. Solid−solid transition in Pb16 and Pb20 are marked with an arrow.

Figure 8. (a) Bond lengths within a cutoff distance of 3.8 Å for sizes Pb20 and Pb23, representing the classes of high and low isomerization temperature clusters, respectively. (b) Frequency distribution of data from (a), indicating the majority (nearly 55%) of atoms in Pb20 possess shorter bond length of 3.25 Å, in comparison to Pb23 that has more evenly spread out distribution of its shortest bonds.

sizes occur at much lower temperatures and are indicated by golden circles in the figure. (2) Within the “high isomerization temperature” class there occur two interesting cases of “restricted diffusion”, for size 16 and 18. The temperature range over which restricted diffusion occurs, is highlighted by magenta rectangles in the figure. (3) Solid−solid transition is observed for sizes 16 and 20, and is labeled accordingly in the figure. (4) δrms of Pb23 is more or less featureless. In what follows, we will explain all these observations in detail and bring out the fact that for such small clusters, it is not “always” possible to define liquid-like state of the cluster unambiguously. High vs Low Isomerization Temperatures. To understand the first observation of low versus high isomerization temperatures, we have carried out a detailed bond length analysis of these clusters. We perform the bond length analysis for clusters Pb20, representing the “high isomerization temperature” class, and Pb23, representing the “low isomerization temperature” class. The difference between the two classes begins to emerge, when we plot the shortest few bond lengths, within the cutoff of 3.8 Å, of these two representative clusters. It is evident from Figure 8a that bonds of Pb23 have a continuous distribution of bond lengths, whereas shortest bonds of Pb20 occur in discrete groups. This indicates that

rearrangements like hopping of a cap atom, or rearrangements of atoms due to processes like restricted diffusion (discussed later in the text). In contrast to this, by “major rearrangements”, we mean structural transitions like a solid−solid transition (discussed later in the text), which is observed in Pb16 and Pb20. In such transitions, the cluster is “trapped” into another structural motif, like hollow motif → core−shell motif, where although most of the atoms diffuse, the diffusion is limited to specific parts of cluster, such as its surface, whereas the cluster preserves its specific structural identity. However, major rearrangements do not imply that the cluster is in liquid-like state, where “all” atoms diffuse throughout the cluster, and its shape, or the structural identity is lost altogether. We note the interesting variations to the trend in eq 1, which show up in the δrms plots as a change of slope but need additional analysis, like MD trajectory analysis, to be singled out. We note the most interesting observations as follows: (1) These clusters can be divided into two classes on the basis of the temperature at which their first distinct isomer occurs. Sizes 16, 17, 18, and 20 belong to the class of clusters with “high” isomerization temperature. Distinct isomers of these clusters appear at relatively higher temperatures and are indicated by red circles in Figure 7, whereas first isomers for the remaining 23703

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707

Article

The Journal of Physical Chemistry C atoms in Pb20 are more organized, with “groups” of many atoms experiencing identical environments in terms of bonding with their neighboring atoms. The shortest bonds of Pb23, on the contrary, exhibit a lack of such grouping, with its shortest bonds forming a more or less continuous distribution of bond lengths. The same graph also shows that nearly all the short bonds of Pb20 are shorter in comparison with the corresponding short bonds of Pb23. Plotting a frequency distribution of these short bonds clears the picture further. Figure 8b shows the frequency distribution of the bond lengths plotted in Figure 8a. We note that the maxima of these two distributions lie at 3.25 Å and 3.35 Å, for Pb20 and Pb23, respectively. Also, nearly 55% of the short bonds of Pb20 are 3.25 Å in length, whereas only about 30% of total short bonds are 3.35 Å in length for Pb23. This makes it clear that the majority of short bonds of Pb20 are identical in length (3.25 Å) and shorter than those of Pb23, whereas bond lengths of Pb23 are more evenly spread over the first nearest neighbor (INN) range of 3.8 Å. Thus, majority of atoms at smaller separation with identical surroundings impart Pb20 greater stability than Pb23. This makes Pb20 isomerize at higher temperature than that of Pb23. It must be noted that although we have exemplified the comparison with representative cluster sizes here, the analysis was performed for the rest of the clusters as well, for which similar conclusions were arrived at. Restricted Diffusion. The second interesting observation marked in the δrms plot is the “restricted diffusion”, which is observed in Pb16 and Pb18. The temperatures over which this phenomenon occurs, are indicated by magenta boxes in Figure 7. In this phenomenon, the atoms in one of the two 5-atom rings of Pb16, and the three-atom ring in Pb18, swap positions within the ring during the MD. Meanwhile, the rest of the atoms in the respective clusters either oscillate about their mean positions or displace only marginally (see Figure 3 for details of the geometries). In both these cases, their geometries are left unaltered once this transition takes place, because swapping atoms within the ring does not alter the shape of the ring or the respective clusters. However, the process of swapping atoms does give rise to added fluctuations in the bond lengths of the clusters, which are reflected as the first jump in δrms, within the magenta boxes of Figure 7. Although “restricted diffusion” can be compared with other similar processes like “automerization”, it must be noted that restricted diffusion is a more specific process. Although automerization is defined as “any rearrangement that yields a degenerate form of the starting material”, restricted diffusion occurs among a subset of atoms (in this case the ring atoms only). Interestingly, there are examples of automerization in clusters as well.26 For instance, the GS−meta stable state−GS transition observed in the case of Sn10 cluster. At finite temperatures, Sn10 undergoes isomerization in such a way that at any instance of GS to GS transition four atoms swap their positions. This process leads to diffusion of “all” the atoms throughout the cluster (after sufficient time of MD) without changing the structural motif of the cluster. In contrast to this, diffusion of “all” atoms throughout the cluster does not occur in restricted diffusion, because the atoms that take part in it are a definite subset (pentagonal ring atoms in the case of Pb16) of the cluster. The reason for occurrence of the restricted diffusion is just as interesting as the phenomenon itself. Restricted diffusion is an attempt by the ring atoms to increase their coordination dynamically. This can be proved by investigating the coordination of atoms that take part in the process, as a

function of simulation time. Figure 9 plots the coordination number of two ring atoms of Pb16 during the process of

Figure 9. Coordination of two different atoms in the five-atom ring of Pb16 over an MD trajectory at 250 K, where restricted diffusion occurs. (a) Coordination of the ring-atom that has coordination 4 at the beginning, which later becomes 5 as the restricted diffusion proceeds. (b) Coordination of ring-atom that has coordination 5 initially, which becomes 4 during restrictive diffusion at the same instant at which that other atom’s, (a)’s, coordination increases from 4 to 5.

restricted diffusion at 250 K. We observe that their coordinations swap values from the initial 4 and 5, to 5 and 4, after about 70 ps of MD. Thus, at any instance during the transition, two of the atoms from the five-membered ring change their coordination, simultaneously. A similar process occurs for ring atoms of Pb18 also. Two of the three atoms in the three-atom ring of Pb18 start out with coordination 7 (not shown in the figure), whereas the third one with coordination 6. A similar graph for the two sets of atoms shows a similar swapping along their MD, where the 6-coordinated atom increases in coordination to 7 dynamically, and vice versa. It is noteworthy that swapping of atoms occurs in atoms of both the rings of Pb16, although not simultaneously. The process occurs multiple times during the course of MD, with about a 150 ps interval, and its frequency of occurrence goes on increasing with temperature. It was observed that the motion of atoms in “restricted diffusion” cannot be explained by a single normal mode of vibration but required 3−4 different normal modes to mimic it. One of the normal modes was found to be the torsion of the five-atom ring, but the frequency corresponding to this mode was only second highest among all normal modes, and not a lower one, thereby suggesting that the bonds of ring atoms are stiff modes, and not soft ones. The activation barrier involved in the process of restricted diffusion of Pb16 was also estimated using the nudged elastic band (NEB) method.57−59 The energy profile of the NEB is given in Figure 10. Geometric configurations corresponding to the intermediate images are displayed overlaid on the graph, where the ring atoms are highlighted in different colors for the visual aid. The activation barrier for the rotation of ring atoms was found to be about 0.16 eV. On a final note, although restricted diffusion may seem to be specific to Pb16 and Pb18, the reason it happens indicates otherwise. Because restricted diffusion is a dynamics attempt of low coordinated atoms to increase their coordination, it is expected to occur in clusters with slight imbalance of coordination in geometrically equivalent atoms. Moreover, it is more likely to take place in small clusters. Solid−Solid Transition. The third observation noted in the δrms plots is the solid−solid transition that is observed only for 23704

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707

Article

The Journal of Physical Chemistry C

the atom at the COM is colored magenta to aid the eye. Comparing the distribution function at all temperatures indicates that a single broad peak, which is the signature of liquid state, is clearly absent at all these temperatures up to fragmentation of the cluster. Thus, it becomes clear that the observed transition is a solid−solid transition, where the structure switches from the hollow GS motif, to the one with an internal atom, without becoming liquid at any temperatures, up to its fragmentation at 850 K. The new structural motif (with internal atom) remains the same up to 850 K. Isomerization occurs before as well as after the solid−solid transition. However, all isomers after the transition, have at least one internal atom, whereas all isomers before solid−solid transition have no internal atom. Diffusion of surface atoms, is also seen until the fragmentation temperature. It is this isomerization process, and not the liquification, that results in increased δrms values before and after the transition. The same can be said about the other cluster Pb20, whose geometries corresponding to this solid−solid transition are shown in Figure 12. A

Figure 10. NEB energy profile of the restrictive diffusion of Pb16. Intermediate image geometries, and their corresponding image index, are shown overlaid on the energy profile for representative changes in the structure. Rotation of atoms throughout the ring becomes clear by these 4 configurations, in which the rotating atoms are highlighted in various colors for the visual aid.

Pb16 and Pb20. During this transition the hollow GS structures of these clusters transform into structures with internal atoms. To differentiate this solid−solid transition from liquification, we have plotted the distribution of atoms about the center of mass (COM) of cluster, computed over 240 ps at four different temperatures for Pb16 and contrasted it with that of Pb17, in which no such transition occurs. Figure 11a shows the

Figure 12. Solid−solid transition in Pb20. Pb20 GS, which does not have any enclosed atoms initially, transits into a new motif with two internal atoms, via the solid−solid transition. Blue spheres represent atoms in the GS that are connected by the shortest bond but are not entirely internal to the structure. Magenta spheres, on the contrary, represent the atoms that are completely enclosed inside the transformed structure. Golden spheres represent the atoms on the surface of the structures. Figure 11. Graph of probability of finding an atom at a given distance from the center of mass, for Pb16 and Pb17. Temperatures for Pb16 are chosen, where major rearrangements of atoms occur. GS motif is preserved at 100 K, where actual geometry is shown alongside its graph. At 550 K, only a solid−solid structural transition is seen, and the transformed geometry is shown alongside the 550 K graph. The cluster does not melt at any temperature above 550 K, until it fragments at 900 K. Data are plotted for 240 ps of MD. Similar graphs for Pb17 show no such solid−solid transition at any temperature up to its fragmentation.

contrasting case of Pb17 is shown in Figure 11b, in which no such transition occurs at any temperature until its fragmentation. The motif of the Pb17 structure remains the same (one with an internal atom) throughout the range of temperatures from 100 K to its fragmentation at 850 K. Thus, to summarize, the solid−solid transition corresponds to the transition of structural motif of Pb16 and Pb20 from one without any internal atom, to a structure containing internal atoms. Because this solid−solid transition is an artifact of a specific GS structural motif, it can therefore occur in clusters that have the GS structure as cages. We also note that once this transition occurs, the cluster visits the previous motif only once, just before its fragmentation. The possible reason could be that the probability of occurrence of the hollow motif reduces with increasing temperature. Isomer Analysis. The fourth and the last observation in the δrms plots is the almost featureless increase in δrms values of Pb23. Although we do not have an exact understanding on why δrms of Pb23 shows a featureless increase, an empirical observation about it in the DE plots (shown in Figure 13) can be used to explain it. In Figure 13, the distance matrix in the Kohn−Sham eigenspace, as defined in refs 60 and 61, is plotted against the energy of isomers with respect to the GS energy, for all sizes. Only the isomers that occur before

distribution for Pb16 at temperatures 100, 500, 550 (solid− solid transition), and 850 K (near fragmentation). At 100 K, the distribution of atoms resembles that of the Pb16 GS, with temperature specific broadening that accounts for oscillation of atoms about their mean positions. The corresponding geometry is a hollow structure without any internal atom and is shown alongside the graph of 100 K in the figure. The width of the peaks goes on increasing until 500 K, which is consistent with the GS type isomers that we find until this temperature. At 550 K, however, there appears a peak at the COM of the cluster. This is a sign that the geometric structure at 550 K (and thereafter) contains an atom at the COM of cluster, which was absent until 500 K. The transformed geometry with an internal atom is shown alongside the graph of 550 K in the figure and 23705

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707

Article

The Journal of Physical Chemistry C

clusters. DE plots are used to explain the continuous thermodynamic transition of Pb23. Restricted diffusion and solid−solid transition can be observed in other clusters as well, depending upon the appropriateness of the cluster. However, the fragmentation observed in the case of Pb clusters is characteristic of clusters of group IV elements.

■

AUTHOR INFORMATION

Corresponding Author

*Phone: +91 20 25902476. E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The authors are grateful to CSIR-4PI and CDAC for availing their computational facility. V.K. thanks MSM (CSC-0129) for the financial support. A.S. acknowledges the funding agency CSIR, New Delhi.

Figure 13. DE plots for Pbn (n = 16−24). The distance in Kohn− Sham (KS) eigenspace is plotted against the energy of isomers (energy of GS = 0) for respective sizes. All sizes exhibit continuity along x-axis, and a jump from GS to the first isomer along the y-axis, except Pb23. Pb23 shows continuity in both directions, in continuum with the GS.

■

REFERENCES

(1) Pawlow, P. The dependency of the melting point on the surface energy of a solid body. Z. Phys. Chem. 1909, 65, 545−548. (2) Berry, R. S.; Smirnov, B. M. Configurational transitions in processes involving metal clusters. Phys. Rep. 2013, 527, 205−250. (3) Proykova, A.; Berry, R. Insights into phase transitions from phase changes of clusters. J. Phys. B: At., Mol. Opt. Phys. 2006, 39, R167− R202. (4) Kanhere, D. G.; Vichare, A.; Blundell, S. A. Reviews of Modern Quantum Chemistry; World Scientific: Singapore, 2002; pp 1568− 1605. (5) Aguado, A.; Jarrold, M. F. Melting and freezing of metal clusters. Annu. Rev. Phys. Chem. 2011, 62, 151−172. (6) Kunz, R. E.; Berry, R. S. Coexistence of multiple phases in finite systems. Phys. Rev. Lett. 1993, 71, 3987−3990. (7) Wales, D. J.; Berry, R. S. Coexistence in finite systems. Phys. Rev. Lett. 1994, 73, 2875−2878. (8) Baletto, F.; Ferrando, R. Structural properties of nanoclusters: Energetic, thermodynamic, and kinetic effects. Rev. Mod. Phys. 2005, 77, 371−423. (9) Schmidt, M.; Kusche, R.; von Issendorff, B.; Haberland, H. Irregular variations in the melting point of size-selected atomic clusters. Nature 1998, 393, 238−240. (10) Bulgac, A.; Kusnezov, D. Thermal properties of Na 8 microclusters. Phys. Rev. Lett. 1992, 68, 1335−1338. (11) Starace, A. K.; Cao, B.; Judd, O. H.; Bhattacharyya, I.; Jarrold, M. F. Melting of size-selected aluminum nanoclusters with 84−128 atoms. J. Chem. Phys. 2010, 132, 034302. (12) Neal, C. M.; Starace, A. K.; Jarrold, M. F. Melting transitions in aluminum clusters: The role of partially melted intermediates. Phys. Rev. B: Condens. Matter Mater. Phys. 2007, 76, 054113. (13) Cao, B.; Starace, A. K.; Judd, O. H.; Jarrold, M. F. Phase coexistence in melting aluminum clusters. J. Chem. Phys. 2009, 130, 204303. (14) Hock, C.; Bartels, C.; Straßburg, S.; Schmidt, M.; Haberland, H.; von Issendorff, B.; Aguado, A. Premelting and postmelting in clusters. Phys. Rev. Lett. 2009, 102, 043401. (15) Krishnamurty, S.; Joshi, K.; Kanhere, D.; Blundell, S. Finitetemperature behavior of small silicon and tin clusters: An ab initio molecular dynamics study. Phys. Rev. B: Condens. Matter Mater. Phys. 2006, 73, 045419. (16) Joshi, K.; Kanhere, D.; Blundell, S. Thermodynamics of tin clusters. Phys. Rev. B: Condens. Matter Mater. Phys. 2003, 67, 235413. (17) Kang, J.; Kim, Y.-H. Half-solidity of tetrahedral-like Al55 clusters. ACS Nano 2010, 4, 1092−1098. (18) Breaux, G. A.; Hillman, D. A.; Neal, C. M.; Benirschke, R. C.; Jarrold, M. F. Gallium cluster “magic melters”. J. Am. Chem. Soc. 2004, 126, 8628−8629.

fragmentation are considered for the purpose. We observe that Pb23 shows isomers that are located toward the bottom left corner of the graph, near the origin, whereas other sizes show a lack of it. Isomers of all the sizes show a continuous distribution of energies (x-direction in the DE plots) but exhibit a jump in eigenvalues (y-direction in the DE plots) of their first isomer. Isomers of Pb23 exhibit continuity in both directions. Because eigenvalues reflect the changes in structure, we believe that there exists a continuity of isomer structures with respect to isomer energies in Pb23 that is not present in other sizes. This continuity of Pb23 in eigenvalues, as well as total energies, implies that the isomers of Pb23 have a more continuous range of structure as well as energy that is not separated from its GS. This avails Pb23, a seamless transitions from its GS structure to its isomers, making its δrms featureless with respect to these transitions.

■

CONCLUSIONS The main focus of this paper is to discuss the thermodynamics of small Pb clusters in detail and compare it with the general features of thermodynamics of small clusters. Putative global minima of Pb clusters between size 16 and 24 were researched in this work, with a new GS structure found for Pb20. The revised GS structures were used to carry out the thermodynamics of these clusters. All GS structures of these Pb clusters are compact cages with internal atoms, except the GS of Pb16 and Pb20. All Pb clusters studied here undergo fragmentation at elevated temperatures and are stable (with an exception of Pb21) above Tm[bulk]. Melting, as we know it conventionally, does not occur in these clusters. Thermodynamics of these small Pb clusters usually follows the path from solid state GS structure, to its isomerization, followed by major rearrangements of atoms, and then fragmentation. Solid−solid transition occurs during thermodynamics of Pb16 and Pb20. The origin of this transition is in their GS structures with all atoms on the surface and it is expected to occur for other clusters with similar GS structural motifs. Restricted diffusion is observed in clusters with equivalent atoms having a slight imbalance in their coordination, which is a generic feature of small clusters with unequally coordinated but geometrically equivalent atoms. In this process, a specific subset of atoms swaps places during constant temperature MD, preserving the overall shape of 23706

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707

Article

The Journal of Physical Chemistry C

(40) Li, X.-P.; Lu, W.-C.; Zang, Q.-J.; Chen, G.-J.; Wang, C.; Ho, K. Structures and stabilities of Pbn (n ≤ 20) clusters. J. Phys. Chem. A 2009, 113, 6217−6221. (41) Li, X.-P.; Lu, W.-C.; Wang, C.; Ho, K. Structures of Pbn (n= 21−30) clusters from first-principles calculations. J. Phys.: Condens. Matter 2010, 22, 465501. (42) Kelting, R.; Otterstätter, R.; Weis, P.; Drebov, N.; Ahlrichs, R.; Kappes, M. M. Structures and energetics of small lead cluster ions. J. Chem. Phys. 2011, 134, 024311. (43) Schafer, S.; Heiles, S.; Becker, J. A.; Schafer, R. Electric deflection studies on lead clusters. J. Chem. Phys. 2008, 129, 44304. (44) Pushpa, R.; Waghmare, U.; Narasimhan, S. Bond stiffening in small nanoclusters and its consequences for mechanical and thermal properties. Phys. Rev. B: Condens. Matter Mater. Phys. 2008, 77, 045427. (45) Kresse, G.; Hafner, J. Ab initio molecular-dynamics simulation of the liquid-metal-amorphous-semiconductor transition in germanium. Phys. Rev. B: Condens. Matter Mater. Phys. 1994, 49, 14251. (46) Kresse, G.; Furthmüller, J. Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set. Phys. Rev. B: Condens. Matter Mater. Phys. 1996, 54, 11169. (47) Kresse, G.; Furthmüller, J. Efficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis set. Comput. Mater. Sci. 1996, 6, 15. (48) Blöchl, P. E. Projector augmented-wave method. Phys. Rev. B: Condens. Matter Mater. Phys. 1994, 50, 17953. (49) Kresse, G.; Joubert, D. From ultrasoft pseudopotentials to the projector augmented-wave method. Phys. Rev. B: Condens. Matter Mater. Phys. 1999, 59, 1758. (50) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized gradient approximation made simple. Phys. Rev. Lett. 1996, 77, 3865. (51) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized gradient approximation made simple [Phys. Rev. Lett. 77, 3865 (1996)]. Phys. Rev. Lett. 1997, 78, 1396. (52) Nosé, S. A unified formulation of the constant temperature molecular dynamics methods. J. Chem. Phys. 1984, 81, 511. (53) Shuichi, N. Constant temperature molecular dynamics methods. Prog. Theor. Phys. Suppl. 1991, 103, 1. (54) Bylander, D.; Kleinman, L. Energy fluctuations induced by the Nosé thermostat. Phys. Rev. B: Condens. Matter Mater. Phys. 1992, 46, 13756. (55) Kaware, V.; Joshi, K. Scaling up the shape: A novel growth pattern of gallium clusters. J. Chem. Phys. 2014, 141, 054308. (56) Mühlbach, J.; Pfau, P.; Sattler, K.; Recknagel, E. Inert gas condensation of metal microclusters. Z. Phys. B: Condens. Matter 1982, 47, 233−237. (57) Sheppard, D.; Terrell, R.; Henkelman, G. Optimization methods for finding minimum energy paths. J. Chem. Phys. 2008, 128, 134106. (58) Henkelman, G.; Jónsson, H. Improved tangent estimate in the nudged elastic band method for finding minimum energy paths and saddle points. J. Chem. Phys. 2000, 113, 9978. (59) Henkelman, G.; Uberuaga, B. P.; Jónsson, H. A climbing image nudged elastic band method for finding saddle points and minimum energy paths. J. Chem. Phys. 2000, 113, 9901. (60) Sadeghi, A.; Ghasemi, S. A.; Schaefer, B.; Mohr, S.; Lill, M. A.; Goedecker, S. Metrics for measuring distances in configuration spaces. J. Chem. Phys. 2013, 139, 184118. (61) De, S.; Schaefer, B.; Sadeghi, A.; Sicher, M.; Kanhere, D. G.; Goedecker, S. Relation between the dynamics of glassy clusters and characteristic features of their energy landscape. Phys. Rev. Lett. 2014, 112, 083401.

(19) Breaux, G. A.; Neal, C. M.; Cao, B.; Jarrold, M. F. Tin clusters that do not melt: Calorimetry measurements up to 650 K. Phys. Rev. B: Condens. Matter Mater. Phys. 2005, 71, 073410. (20) Pyfer, K. L.; Kafader, J. O.; Yalamanchali, A.; Jarrold, M. F. Melting of size-selected gallium clusters with 60−183 atoms. J. Phys. Chem. A 2014, 118, 4900−4906. (21) Breaux, G. A.; Benirschke, R. C.; Sugai, T.; Kinnear, B. S.; Jarrold, M. F. Hot and solid gallium clusters: too small to melt. Phys. Rev. Lett. 2003, 91, 215508. (22) Shvartsburg, A. A.; Jarrold, M. F. Solid clusters above the bulk melting point. Phys. Rev. Lett. 2000, 85, 2530−2532. (23) Breaux, G. A.; Cao, B.; Jarrold, M. F. Second-order phase transitions in amorphous gallium clusters. J. Phys. Chem. B 2005, 109, 16575−16578. (24) Pavone, P.; Baroni, S.; de Gironcoli, S. α ↔β phase transition in tin: A theoretical study based on density-functional perturbation theory. Phys. Rev. B: Condens. Matter Mater. Phys. 1998, 57, 10421− 10423. (25) Chacko, S.; Joshi, K.; Kanhere, D. G.; Blundell, S. A. Why do gallium clusters have a higher melting point than the bulk? Phys. Rev. Lett. 2004, 92, 135506. (26) Joshi, K.; Kanhere, D.; Blundell, S. Abnormally high melting temperature of the Sn10 cluster. Phys. Rev. B: Condens. Matter Mater. Phys. 2002, 66, 155329. (27) Ghazi, S. M.; Lee, M.-S.; Kanhere, D. The effects of electronic structure and charged state on thermodynamic properties: An ab initio molecular dynamics investigations on neutral and charged clusters ofNa39, Na40, and Na41. J. Chem. Phys. 2008, 128, 104701. (28) Lee, M.-S.; Chacko, S.; Kanhere, D. First-principles investigation of finite-temperature behavior in small sodium clusters. J. Chem. Phys. 2005, 123, 164310. (29) Cheng, H.-P.; Berry, R. S. Surface melting of clusters and implications for bulk matter. Phys. Rev. A: At., Mol., Opt. Phys. 1992, 45, 7969−7980. (30) Susan, A.; Kibey, A.; Kaware, V.; Joshi, K. Correlation between the variation in observed melting temperatures and structural motifs of the global minima of gallium clusters: An ab initio study. J. Chem. Phys. 2013, 138, 014303. (31) Joshi, K.; Krishnamurty, S.; Kanhere, D. “Magic melters” have geometrical origin. Phys. Rev. Lett. 2006, 96, 135703. (32) Susan, A.; Joshi, K. Rationalizing the role of structural motif and underlying electronic structure in the finite temperature behavior of atomic clusters. J. Chem. Phys. 2014, 140, 154307. (33) Starace, A. K.; Neal, C. M.; Cao, B.; Jarrold, M. F.; Aguado, A.; López, J. M. Electronic effects on melting: Comparison of aluminum cluster anions and cations. J. Chem. Phys. 2009, 131, 044307. (34) Li, H.; Chen, W.; Wang, F.; Sun, Q.; Guo, Z.; Jia, Y. Tin clusters formed by fundamental units: a potential way to assemble tin nanowires. Phys. Chem. Chem. Phys. 2013, 15, 1831−1836. (35) Li, H.; Ji, Y.; Wang, F.; Li, S.; Sun, Q.; Jia, Y. Ab initio study of larger Pbn clusters stabilized by Pb7 units possessing significant covalent bonding. Phys. Rev. B: Condens. Matter Mater. Phys. 2011, 83, 075429. (36) Sattler, K.; Mühlbach, J.; Recknagel, E. Generation of metal clusters containing from 2 to 500 atoms. Phys. Rev. Lett. 1980, 45, 821−824. (37) Rajesh, C.; Majumder, C.; Rajan, M.; Kulshreshtha, S. Isomers of small Pbn clusters (n= 2−15): Geometric and electronic structures based on ab initio molecular dynamics simulations. Phys. Rev. B: Condens. Matter Mater. Phys. 2005, 72, 235411. (38) Senz, V.; Fischer, T.; Oelßner, P.; Tiggesbäumker, J.; Stanzel, J.; Bostedt, C.; Thomas, H.; Schöffler, M.; Foucar, L.; Martins, M.; et al. Core-hole screening as a probe for a metal-to-nonmetal transition in lead clusters. Phys. Rev. Lett. 2009, 102, 138303. (39) Wang, B.; Zhao, J.; Chen, X.; Shi, D.; Wang, G. Atomic structures and covalent-to-metallic transition of lead clusters Pbn (n = 2−22). Phys. Rev. A: At., Mol., Opt. Phys. 2005, 71, 033201. 23707

DOI: 10.1021/acs.jpcc.5b05250 J. Phys. Chem. C 2015, 119, 23698−23707