Full Atomistic Simulations of Nanocolloidal Solutions: Formation of Clusters, Aggregates, and Gels

In colloidal chemistry, attention is often paid to the interactions among colloidal units regardless of the details of structures of solvents and colloids. One of the popular approaches is based on the mean-field type approximation, in which the solution is treated as a homogeneous medium. However, the results of the recent MD simulations in atomistic-level have revealed that structural and dynamical changes of both solvent and salts play a role in spontaneous structure formations of aggregates and gels [1]. This situation will be discussed in this chapter.

Spontaneous Formation of Aggregates and Gels from the Nanocolloidal Silica in NaCl Solutions

In Ref. [1], coagulation of nanocolloidal silica in NaCl solutions is described based on molecular dynamics simulations using a

Molecular Dynamics of Nanostructures and Nanoionics: Simulations in Complex Systems Junko Habasaki

Copyright © 2021 Jenny Stanford Publishing Pte. Ltd.

ISBN 978-981-4800-77-8 (Hardcover), 978-1-003-04490-1 (eBook) www.jennystanford.com fully atomistic model. This approach enables us to examine the formation of clusters, aggregates, and gels without any assumption of mechanisms of formation and therefore new insights in colloidal chemistry can be brought.

For examining a colloidal system, a coarse-grained model is often used, and details of the solvent tend to be neglected, as already mentioned. A main reason for such treatments seems to be a large size difference between colloidal units and solvent molecules, which makes it difficult to treat the system in fully atomistic levels. Furthermore, dynamics during the formation of gels are slow, as observed in the glass transition of glass-forming materials, where heterogeneous dynamics are often found. In the case of the treatment of heterogeneous systems having several length scales, a large system and a long-time simulation are required. Therefore, molecular dynamics simulations of the colloidal system in fully atomistic models are challenging, although obviously, they are useful for a better understanding of coagulation processes in nanoscopic levels.

The effect of salt on the coagulation process is experimentally well known and tends to be explained by the DLVO (B. Derjaguin, L. Landau, E. Verwey, and T. Overbeek) model [2, 3] in many cases; however, many problems remain and are a subject of debate.

For example, the effect of electrolytes on dispersions depends on the kind of ions. Such a specific ion effect [4] is not fully explained by the theory. Furthermore, multivalent salt is known to be more efficient in reducing the effective charge of the colloidal solution. This multivalent ion effect [5, 6] is another example not to be well explained by the DLVO model.

In the present chapter, the origins of forces, the role of the solvent affected by salt for coagulation of nanocolloidal-silica— water-NaCl system will be discussed based on structural and dynamical changes.

Especially, the following changes in water structure cannot be clarified without atomistic and/or molecular level simulations. With the increased concentration of the NaCl, a pair correlation function of the OW-OW (OW means an oxygen atom in water) pair shows a sub-peak characteristic to the formation of high-density water (HDL) related to the so-called polyamorphism [7] of water.

Parallel with the structural changes of water, the network structure of the NaCl part also changes and these structures are affected by the presence of the silica part in the system as shown later. Therefore, the picture based on the interaction among colloidal units with a mean field approximation of solvent seems to need revision for considering coagulation processes.

Molecular Dynamics Simulation of Colloidal Systems

In the atomistic simulations shown in Chapters 8 and 9, potential parameters used in [1] are chosen to maintain consistency of quality and combinations as much as possible, although further improvement of these parameters might be necessary, depending on the purpose. For the colloidal silica and resultant silica part of aggregations and gels, the potential model developed by Tsuneyuki et al. [8] (the TTAM model) was used. It is worth to mention that the model of silica was not necessarily designed for the simulation of colloids and/or gels. However, it allows the breaks and restructures of the Si-0 bond because it is determined based on the ab initio MO calculation with the deformation of

4-

an Si0₄ + 4e⁺ cluster including changes in an Si-0 distance. For water, an SPC/E [9] model is used, and it is treated using the shake algorithm with a constraint. This model water has been well tested (for example, see [10]) and parameters for interactions with NaCl salt are available [11]. The interaction between silica and Na⁺ or Cl" was obtained from the CLAYFF Force Field [12].

Further details of the potential functions and parameters used are described in the Appendix and Ref. [1].

Ensembles used in each run will be shown in the next subsection. The DL_POLY (version 4.01) program [13] was used for the MD simulations.

Modeling of Colloidal Systems

One of the problems is how to prepare the complex structures formed by several components and or units used in MD simulations. In this work, equilibrated colloidal units surrounded by solvent molecules are prepared and then they are combined to form larger systems. This procedure will be described below.

Preparation of Colloidal Units

First, an initial colloidal unit composed of Si0₂ was prepared as an isolated partial system.

There are six atomic species in the system (Si, 0, Na, Cl, OW, HW). If we would like to distinguish oxygen and hydrogen in water from those belonging to other species, notations OW and HW are used, respectively. Among combinations of pairs, some repulsive interactions can be excluded from further treatments when direct contact is rare.

Here we consider the colloidal units formed by silica networks.

Figure 8.1 An example of an initial structure of a colloidal unit represented using a ball-and-stick model.

First, a colloidal silica unit consisting of 108 atoms was equilibrated in the NVE ensemble with a periodic boundary condition (PBC), where the volume was fixed to the value obtained from the NPT ensemble (using a Nose-Hoover thermostat and a barostat [14-16]) under atmospheric pressure. After equilibration at 300 К under constant pressure (0.1 MPa), the PBC was removed.

In the system with PBC, the Coulombic force is calculated by the particle-mesh method [17], while in the system without PBC, the Coulombic force was calculated as a direct sum. During this procedure, the temperature increased to allow the melting of the system and a cluster with a rounded shape was gradually formed. In Fig. 8.1, an example of a structure of a colloidal unit is illustrated using a ball-and-stick model. The characterization of the colloidal unit thus prepared can be found in Ref. [1].

The colloidal unit was solvated with SPC/E water within a cubic basic box. Then Na⁺ and CD ions were added at random positions.

Preparations of the Complex System with Colloidal Units, Water and Salt

A solvated system containing one colloidal unit and salt was equilibrated at 300 К in the NPT condition with PBC. After equilibration, an isotropic basic cell was prepared for the large-scale MD simulations by creating 64 (4x4x4) repetitions of the solvated colloidal unit cell.

Six systems shown in Table 8.1 with PBCs were examined in Ref. [1]. Some descriptions of the systems, including the number of colloidal units and atoms, the number of Na⁺-Cl" pairs, and the concentration of salt in water, with the status of the system, are provided in this table. These results and some new results will be discussed in this chapter. (Here the concentration is given as the ratio to water and not to solution. This choice seems to be better for comparison among systems because the status of after the system also depends on the amount of colloidal silica).

For each system, the MD simulation was first performed in the NPT ensemble, and the run in the NVE ensemble then followed the equilibration of the system.

The I-A, I-B, and I-C and II-A, II-B, and II-C systems examined are characterized as follows. In series I, with smaller water content, we have observed the systematic formation of clusters, aggregates, and a gel. In series II with larger water content, each system shows no large aggregates, including the reference system II-A without salt.

Table 8.1 Characteristics and status of systems examined in Ref. [1], Chapters 8 and 9

Each system except I-C was equilibrated at 400 К until the energy of the system became nearly constant, and the temperature was then decreased to 300 K. The size of the basic cell for the MD simulations, L, of the complex system addressed the present work was 65-90 A. To examine g[r) at a long length scale, an additional run was performed in each system with a long cutoff length of ~L/2 (for both the repulsive force in the real space and the real part of the electrostatic force), whereas a cut-off length of 12 A was used in all other cases. Each time step was set to 1 fs. After the construction of the initial configurations, the systems were equilibrated for 1-3 ns. When the energies and structures of one 100,000-step run (1 step=l fs) were essentially maintained in the subsequent 100,000-step run, the system was regarded as sufficiently equilibrated one, and the simulations were continued on a longer time scale. After the equilibration of the system, the ensemble was changed to NVE, using the volume obtained from the NPT run. The results for the time range of 400 ps to 2 ns (400,000-2,000,000 steps) at 300 К after the equilibration run were used in analyses.

The formation of a gel was observed in the following condition in I-C system. Beginning from the I-B configuration, the temperature was again increased to 400 К to accelerate coagulation. We observed the developments of infinitive networks in this system, in which the atoms became immobile as shown by the mean-squared displacement of the system (see Section 8.6). This system is regarded as a gel in the present model. The gel thus obtained was quasi-equilibrated one because of the slow dynamics.

Structural Changes with Coagulation Caused by Salt

The structural characteristics of each system in Table 8.1 will be discussed further in this subsection.

As naturally expected, the structural changes by the development of aggregates of colloidal silica were observed with increasing content of salt.

In the case of examining the structural changes on the colloidal surface or reaction with the coagulation process itself, the works by ab initio MD or classical MD using a reactive force field can be a useful and better approach; however, our concern here is equilibrated systems near the potential minimum after the coagulation processes. In such cases, potential parameters of silica, which allows the rearrangement of Si-0 bonds, are good enough, even if they were not designed to represent the reaction process or chemical states on the surface. Both bonds breaking and reconnection allow structures after the long relaxation time to be meaningful. In fact, changes of the (quasi)- equilibrated structures are systematic as shown in the following subsections.

In MD simulations, initial configurations are not necessarily equilibrated ones, so that time development of the structures depends on the choice of them. Therefore, for a kinetic argument for the formation of aggregations or gels, different approaches might be required.

8.2.4.1 Changes in network structures of silica and water parts

To learn the characteristics of networks, structures of the silica part were visualized, and then their fractal dimensions were examined. The calculated dimension for silica part was comparable to the experimental observations for the aggregation of silica colloids, and at longer length scales, super-aggregation was evident in the gelation process.

Clear changes in both the structure and mobility of the water were observed as the NaCl concentration increased, suggesting the importance of the solvent structures in these processes. The structural change in the water that was observed in the pair correlation function is interpreted as the change from a lower- density liquid to a high-density liquid (HDL) as shown later.