Machine learning of atomic dynamics and statistical surface identities in gold nanoparticles

Metal nanoparticles exhibit size-, shape-, and surface-dependent properties distinct from bulk materials, but probing their atomistic structure and dynamics is challenging. Gold nanoparticles (Au NPs) are particularly interesting due to catalytic activity (e.g., CO oxidation), surface plasmon resonance, and size-dependent melting behavior, which implies significant surface atomic mobility even at relatively low temperatures. The research question is how to characterize, at atomistic resolution and under relevant thermal conditions, the stability, lifetimes, and interconversion dynamics of the different atomic environments (AEs) on Au NP surfaces and interiors. The purpose is to develop and apply a robust, data-driven machine learning pipeline to decode NP atomic dynamics from MD simulations, quantify native vs non-native AEs, and estimate a statistical surface identity and kinetic exchange among AEs across temperatures. This is important for understanding and predicting functional properties such as catalysis, where transient surface sites may govern activity.

Experimental advances (e.g., single-particle and single-atom resolution HAADF-STEM) have directly observed atomic motion in supported Au NPs at finite temperatures, yet ensemble averaging and lack of atom identity tracking limit quantitative dynamical insight. Prior theoretical/computational work has successfully used molecular simulations to study metals and Au clusters/NWPs, and ML methods have been applied to analyze MD trajectories. SOAP descriptors with unsupervised clustering have recently enabled reconstruction of structural/dynamical complexity and data-driven materials classification. However, comprehensive quantification of AE lifetimes, transition rates, and native vs non-native AE formation across NP morphologies and temperatures remained limited, motivating the present combined bottom-up and top-down SOAP-based approach.

- Systems: Ideal Au nanoparticles of multiple morphologies and sizes, including Ih309 (icosahedral, 309 atoms), Dh348 (decahedral, 348 atoms), and TO309 (truncated octahedral/cuboctahedral, 309 atoms). Additional larger Ih, Dh, and To NPs (Ih923; Dh1086, Dh1734; To807; To976) were used to build a comprehensive SOAP dictionary of environments. - Simulations: Classical MD with the second-moment tight-binding (SMATB) potential (Gupta/Rosato–Guillope–Legrand/Cyrot-Lackmann formalism) in LAMMPS. Initial energy minimization, short thermalization (20,000 steps, 1 fs), then production runs at 300 K, 400 K, and 500 K for 2 μs in the canonical ensemble using a Langevin thermostat (timestep 5 fs, damping 100 ps). Steady state was reached, and analyses were performed on 1000 frames sampled every 1 ns from the last 1 μs. - Descriptors: SOAP power spectra (dscribe) computed for each atom in each frame with rcut = 4.48 Å, lmax = 8, nmax = 8. Each SOAP vector has 576 components (324 unique). SOAP chosen for its non-discrete, robust capture of local environments and surface details. - Bottom-up analysis: Dimensionality reduction via PCA (scikit-learn) to first 3 PCs (~99% variance). Unsupervised clustering via HDBSCAN* (min_cluster_size=125, selection_method="eom") trained on Ih309 at 300 K; predictor applied to Ih309 at 400 K and 500 K. Noise points assigned to nearest exemplar cluster. Transition matrices were computed by counting i→j and i→i events between consecutive sampled frames (Δt = 1 ns) and row-normalized to probabilities. - Top-down dictionary approach: Built a dictionary of 47 SOAP AEs from minimized ideal NPs of Ih, Dh, and To families. Pairwise SOAP distances computed; hierarchical clustering (scipy, complete linkage) produced a global dendrogram. A cut at dSOAP = 0.08 coarse-grained the 47 AEs into 10 macro-clusters (b, ss, ss', c, c', s, e, e', v, v'). During analysis, each observed AE is assigned to the nearest dictionary entry (by SOAP distance) and then to its macro-cluster. Transition matrices and chord diagrams were constructed similarly for these top-down AEs to quantify residence and exchange. - Exemplary kinetic analysis: Transition rates k_ij ≈ p_ij/Δt, lifetimes τ_i estimated from dominant exit pathways; counts of observed transitions and AE populations are provided in Supplementary figures.

- Bottom-up identification of AEs in Ih309 at 300 K: Eight distinct SOAP clusters corresponding to structural AEs were detected, including Ico (central atom), Bulk, SubSurf, 5foldedSS (under vertexes), surface Faces (close-packed FCC(111)), Edges, Vertexes, and an emergent Concave AE (rosette centers) not present in the ideal NP. - Free-energy/density landscape and PCA show compact (less mobile) bulk/subsurface zones and diffuse (more mobile) surface AEs at 300 K; Concave AEs bridge interior and exterior, indicating formation of local point-defects (rosettes). - Transition matrices (Δt=1 ns) reveal high stability in core AEs (p_ii ≈ 1) and significant exchange among surface AEs even at 300 K. Example: Faces→Edges transition probability p ≈ 0.12 (×10−2 units matrix corresponds to 0.12), giving k ≈ 0.12 ns−1 and τ ≈ 8.3 ns; observed ~6700 events in the last 1 μs. - Temperature effects (Ih309): At 400 K and 500 K, surface AE minima merge on the FES; communication between Faces, Edges, and Vertexes increases. At 500 K, surface residence probabilities drop to near or below 50%, indicating surface pre-melting and liquid-like dynamics, while overall icosahedral shape persists. Rare central-void-like events cause transient similarity of deep bulk SOAP spectra to surface AEs at 500 K. - Direct atomic diffusion: Tracking a vertex atom at 300 K shows genuine surface diffusion and subsurface penetration, visiting multiple surface and subsurface AEs, including transient conversion into a concave (rosette center) AE around ~820 ns. - Top-down dictionary classification: Coarse-grained 10 AE groups (b, ss, ss', c, c', s, e, e', v, v') enable native vs non-native assignment. In Ih309, non-native surface AEs increase with temperature; concave c and c' and certain vertex/edge variants emerge at finite T. From transition matrices, c'→s has p ≈ 0.4 at 300 K (k ≈ 0.4 ns−1; τ ≈ 2.5 ns; ~4700 events in 1 μs). - Morphology comparison: Dh348 is more stable than Ih309 at all temperatures, with interior (b, ss) largely decoupled from surface even at 500 K and solid-like dynamics on the surface. TO309 is more dynamic than Dh348; at 300 K it shows ~24% of surface atoms in non-native AEs and is less stable than Ih309, but at 500 K TO309 is more stable than Ih309 with more discrete (solid-like) surface dynamics.

The study resolves how different atomic environments populate Au nanoparticles and interconvert over time, addressing the challenge of quantifying site-specific stability and dynamics. The bottom-up SOAP+PCA+HDBSCAN* approach uncovers the set of metastable AEs and their connectivity, while the top-down dictionary maps observed AEs to native vs non-native families across morphologies. Together, these reveal that surface dynamics arise from genuine atomic diffusion (not just local oscillations), with increasing exchange among surface AEs as temperature rises and evidence of pre-melting at 500 K for Ih309. Concave (rosette) environments act as bridges between interior and surface and emerge from both subsurface and surface states. The quantified lifetimes and transition rates provide a kinetic map that can be linked to potential catalytic behavior: only AEs with lifetimes comparable to or exceeding reaction characteristic times are likely to contribute significantly to reactivity. The framework thus connects intrinsic structural dynamics to functional properties and enables comparative assessment across NP shapes, showing higher stability for Dh348 and complex temperature-dependent stability for TO309 vs Ih309.

This work introduces a combined bottom-up and top-down machine learning pipeline to decode atomistic dynamics in Au nanoparticles from MD trajectories. Using SOAP descriptors, dimensionality reduction, and density-based clustering, the approach identifies and tracks native and emergent (non-native) atomic environments, quantifies their populations, interconversion probabilities, lifetimes, and temperature dependence. A general SOAP dictionary spanning icosahedral, decahedral, and truncated-octahedral AEs enables cross-morphology classification and estimation of a statistical surface identity. Key contributions include: (i) unbiased detection of surface and subsurface AEs and their kinetic interconnections; (ii) evidence of genuine atomic diffusion even at 300 K; (iii) temperature-driven merging of surface states and pre-melting behavior at 500 K for Ih309; and (iv) morphology-dependent stability trends (Dh348 most stable; TO309 highly dynamic at 300 K but more stable than Ih309 at 500 K). Future directions include expanding the AE dictionary to additional geometries and metals, integrating reactive environments under adsorbates, refining time resolution and force-field fidelity, and coupling AE kinetics with microkinetic models to predict catalytic performance under operando conditions.

- Force-field dependence: Although different force fields (e.g., EAM, ML potentials) yielded similar qualitative dynamics, exact lifetimes and probabilities may vary with the chosen potential. - Temporal resolution: Analyses are based on frames spaced by Δt = 1 ns; processes faster than this are not resolved and low-probability events (p < 0.01) are treated qualitatively. - Top-down vs bottom-up AE definitions: The two analyses use different AE partitions; direct kinetic comparisons are non-trivial and depend on the dendrogram cut (dSOAP = 0.08) in the dictionary. - System scope: Simulations consider ideal, unsupported NPs in vacuum with a thermostat; real supports, environments, and chemisorption could alter AE stability/dynamics and may introduce new AEs. - Generalizability of quantitative values: Reported probabilities and timescales pertain to the simulated sizes and conditions; extrapolation to other sizes or conditions should be done cautiously.