Development of a coarse-grained model for surface-functionalized gold nanoparticles: towards an accurate description of their aggregation behavior

Understanding the dispersion stability and aggregation propensity of self-assembled monolayer gold NPs at a molecular level is crucial to guide their rational design and to inform about the optimal surface functionalization for specific applications. To reach this goal, in silico modeling via coarse-grained (CG) molecular dynamics (MD) simulations is a fundamental tool to complement the information acquired from experimental studies since CG modeling allows to get a deep knowledge of the molecular interactions that take place at the nanoscale in this kind of systems. Unfortunately, current CG models of monolayer-protected AuNPs present several drawbacks that limit their accuracy in certain scenarios. We here develop a CG model that is fully compatible and extends the SPICA/SDK (Shinoda–DeVane–Klein) force field. Our model allows reproducing the behavior of AuNPs functionalized with hydrophobic as well as charged and more hydrophilic ligands. This model improves upon results obtained with previously derived CG force fields and successfully describes NPs aggregation and self-assembly in aqueous solution.


Supplementary Results
The dynamics of multiple NPs aggregation depends on their shell composition To understand the dynamics of the aggregation between multiple NPs, we introduced the concept of geometrical cluster. We defined a geometrical cluster as a group of neighboring NPs which are at a given distance.
A geometrical cluster, in general, does not describe a physical aggregate but it is merely used to investigate how close the NPs are in the simulation box. We evaluated geometrical clusters at increasing distances (from 2 to 13 nm). Figure S3 (panels A, C, E, G) shows, as heatmaps, the number of geometrical clusters that can be found at a given cut-off distance and as a function of time (rows) for each of the three replicas (columns). Among all possible distances, we further selected a specific subset that corresponds to the geometrical cluster at a 3.5 nm cut-off. This subset, in particular, was used to define physical aggregate. We then analyzed the number of NPs within each geometrical cluster Figure S3 (panels B, D, F, H) for this particular subset. The maps displayed in these panels allow to describe the aggregation process qualitatively as a function of time.
As shown in Figure S3, the overall aggregation dynamics depends on the NP type and CG forcefield.
For 100%OT NPs, the heatmaps obtained with both force fields (SPICA and MARTINI) are qualitatively similar: the number of aggregates reaches its maximum in the early stage of the simulations (∼50 ns) ( Figures   S3A and S3C), mainly generating dimers and trimers until complete aggregation, which was found on average at ∼700 ns ( Figures S3B and S3D). In addition, in this kind of NPs, geometric clusters can join together with high probability if their distance is between 4 and 6 nm. This is graphically represented in Figures S3A and S3C. There, each region of continuous color constitutes a distance interval characterized by the same number of geometrical clusters and, therefore, the highlighted circles at the right side of these regions serve as a qualitative assessment of the aforementioned regrouping.
Regarding the 50%OT NPs, the heatmaps corresponding to the geometrical clusters show an entirely different profile for SPICA and MARTINI. In the SPICA model, the maximum number of geometrical clusters (∼10) is reached in ∼ 1µs ( Figure S3E), and dimers and trimers persist until the end of simulation time ( Figure   S3F). In contrast, in the simulations with MARTINI, the maximum number of groups (∼8) is reached in ∼80 ns. Interestingly, the heatmaps corresponding to the geometrical clusters are similar to those obtained for the 100%OT NPs, but the NPs do not aggregate into a single group (Figures S3E and S3F) in any of the three replicas.
In summary, these analyses highlight that cluster formation is generally slower for 100%OT NPs, and that 50%OT MARTINI and 50%OT SPICA NPs exhibited disparate behaviors. Specifically, 50%OT SPICA NPs remained stable in the solution, whereas 50%OT MARTINI NPs rapidly formed aggregates.