Wettability of organically coated tridymite surface – molecular dynamics study

Introduction

Alkylsiloxanes represent a group of self-assembled monolayer (SAM) substances. They are widely used in combination with oxide surfaces such as SiO₂ or glass [1–3]. Alkylsiloxanes easily react with hydroxylated oxide surfaces forming a self-assembled monolayer covalently bound to the surface via Si–O–Si bridge [4]. This leads to superior stability of the surface and allows further modification without deterioration of the monolayer. The first report of SAM preparation was by Sagiv [5]. After that followed studies investigating a mechanism of the SAM formation, e.g., [2]. Sagiv proposed that SAMs are formed by hydrolysis of the chlorosilanes followed by reaction with the hydroxyl groups on the surface and other silane molecules to obtain a polymerized network. Other groups proposed a mechanism which included the reaction of self-assembling molecules with the adsorbed water layer on the substrate [1]. This leads to the formation of a two-dimensional, cross-linked network of Si–O–Si bonding with only a few bonds to the substrate. Preparation of such surfaces is depending on several parameters like water content, temperature, the utilized solvent or age of the solution. For example, preparation of high quality SAMs of alkyltrichlorosilane derivatives are not simple to produce mainly because of the need to carefully control the amount of water in solution [4]. It was shown that in the absence of water incomplete monolayers are formed, whereas excess of water results in polymerization in solution and polysiloxane deposition of the surface. Another influence on forming the SAMs is the use of precursors. For example, in the work by Bierbaum et al., it was shown that the reaction mechanism of trichlorosilane and trimethoxysilane is different [6]. The methoxy groups of trimethoxysilanes are less reactive than chlorines in trichlorosilanes what leads to suppressed hydrolysis and oligomerization in solution. Instead of oligomers, mainly monomers react with the surface, and the two-dimensional polymerization is sterically hindered compared to trichlorosilanes. This leads to lower density of alkyl chains on the surface. Different density of the alkyl chains can affect directly the wettability of the surface. For example, Bachmann et al., showed by the contact angle (CA) measurement that a higher percentual coating of dimethylsiloxane on silt soil modifies the hydrophilic surface to almost hydrophobic [7].

Except experimental measurements of wettability, molecular simulation methods such as Monte Carlo (MC) or molecular dynamics (MD) can be used to characterize hydrophobicity/hydrophilicity of surfaces. It was shown that the theoretical approach is suitable for modeling of the contact angle and obtained results can be directly compared to experimental measurements [8]. So far, several theoretical studies of the surface wetting phenomena of silanized surfaces were performed. For example, a different degree of surface hydroxylation/silanization of amorphous silica was investigated by molecular dynamics simulation [9]. It was observed that the hydrophobicity becomes enhanced with an increase in the degree of surface silanization. In another work, MD simulations were performed to study the impact of different headgroups of SAMs bound to the (111) silicon surface [10].

In this paper we present a theoretical study using classical molecular dynamics simulations on models of silanized tridymite surfaces. The work is aimed to explain the observed difference between experimentally measured CA obtained for the same surface functionality on silica-glass surfaces, dimethylsiloxane (DMS), but prepared from two different silanes – dichlorodimethylsilane (DCDMS) or dimethyldiethoxysilane (DMDES). For this purpose, several models with a different coating by DMS were prepared to see an impact on the simulated CA. Particularly, two models with 50 and 100 % coverage by the monomeric DMS units bound to the surface through one Si–O–Si bridge (monodentate binding) were constructed. These models were completed with a model with 100 % coverage by DMS units having bidentate binding to the surface (two Si–O–Si bridges). Finally, a model mimicking a possible dimerization between DMS units was suggested as well.

Experimental details

Computational details

All models with a different degree of coating were derived from the (001) surface of tridymite (Fig. 1a). The SiO₂ polymorph tridymite was used instead of the more complicated structure of amorphous silica glass for the model preparation. A slab model with the (001) surface was cut from the monoclinic structure of tridymite (C1 symmetry, a= 10.0140 Å, b= 8.6004 Å, c= 8.2169 Å, α = 90.000°, β = 91.512°, γ = 90.000°). This slab with a thickness of about 9 Å consists from two linked planes of the SiO₄ tetrahedra and on its both sides top oxygen atoms are saturated by hydrogen atoms forming surface hydroxyl groups (Fig. 1a). These OH groups are substituted by siloxane functional groups during the silanization process. The area per chain of the active sites on the (001) tridymite surface is ∼21.5 Ų, which corresponds well to a typical full coverage of silica surface by SAM (area per chain of 22–25 Ų [13]). The extended model of the surface was created by a multiplication of the basic cell vectors in a and b directions (101.140 × 103.2048 Å) and imposing a vacuum in the c direction with a thickness of about 200 Å.

Fig. 1:

Structure of used models: (a) tridymite, (b) DMS_50, (c) DMS_100, (d) DMS_bi, (e) DMS_2. Yellow balls represent Si atoms, red O atoms, cyan C atoms, and white H atoms.

The models with the silanized surfaces were prepared by replacing surface OH groups with DMS units that have two Si–O– functionalities for binding. Considering also a distribution of the surface OH groups (distance between neighboring sites is about 5 Å) several models were prepared. In the first model with 100 % coating (DMS_100; Fig. 1c), each surface OH group was replaced by the DMS unit with a monodentate Si–O–Si binding. The second Si–O– site of DMS was saturated by a hydrogen atom forming a hydroxyl group. Thus, in this model, the surface is formed from two hydrophobic –CH₃ groups and one hydrophilic –OH group, respectively. As incomplete coating can be one of possible explanations for different CA of DMS-coated silica surfaces prepared from two different precursors, the model with only 50 % coating was constructed, too. In this model only each second surface OH group was replaced by a DMS unit (DMS_50 model in Fig. 1b). In the model DMS_bi (Fig. 1d), each DMS unit is bound to the surface through two Si–O–Si bridges forming bidentate two-site binding at a full coverage. Some experimental studies suppose also an additional lateral (poly-) oligomerization between SAM units. Therefore, the following model reflects this scenario. However, each DMS unit bound to the surface has only one –OH group available for binding, only a dimerization between the neighboring sites is possible with forming a lateral Si–O–Si bridge in a monolayer coating. This model (DMS_2) with a full coverage is displayed in Fig. 1e.

The simulation of the wetting properties of all five models was performed by using a water droplet consisting of 1000 water molecules placed in the vicinity of the surface (the size of the water droplet was shown to be sufficient in the study of the wetting properties of the kaolinite surfaces [14]). All MD simulations were performed using the LAMMPS force-field (FF) package [15]. The interactions between atoms in the tridymite slab were described by means of CLAYFF force field [16]. The CLAYFF field was developed for a study of hydrated crystalline compounds and their interfaces with fluid phases. The FF parameters for the siloxane units were of the OPLS-AA type [17–19]. Water molecules were expressed by a point charge FF model (SPC/E) [20]. The SHAKE [21] algorithm was used to keep the water molecules rigid during the simulations. Further, also Si and inner O atoms of the tridymite slab were fixed and only the rest of the slab was allowed to move (surface OH groups and coating moieties). The van der Waals interaction parameters between different atomic species i and j were calculated using the arithmetic combination rule.

In MD, the Newton’s equations of motions were numerically integrated by using the velocity Verlet algorithm [22] with a time step of 1.0 fs. The van der Waals interactions were truncated using a 12 Å cutoff. The Coulombic electrostatic interactions were computed using the Particle-Particle Particle-Mesh (PPPM) algorithm [23]. Cutoff for these interactions was always approximately half of the computational cell avoiding long range interactions between water droplet images in the neighboring cells. All simulations were performed at temperature T = 300 K controlled by rescaling atomic velocities during the MD run. The total simulation lengths were 4 ns for each model. First part of MD was an equilibration phase and last 2 ns were taken for an analysis (10 000 snapshots in the collected data).

The microscopic contact angle was calculated from a fit of a circle equation to the contour of the water nanodroplet. The droplet boundaries were obtained using the method of cylindrical bins. This procedure was described in detail in our previous paper [14]. Briefly, the water droplet was discretized to cylindrical bins using steps in a lateral dimension of ΔR≈ 1.0–1.5 Å and in a perpendicular dimension of Δz ≈ 0.5 Å, respectively. For each bin the average water density was calculated. The droplet boundary was considered when the water density of a bin was half of the bulk water density (0.5 g/cm³). From obtained points the contour curve was constructed and extrapolated to circle. Then the contact angle was calculated at the cross-section between boundary droplet curve and surface line.

Results

Tridymite surface

The Fig. 2a shows the front and top view of the last snapshot collected after a 4 ns MD run of a water droplet on the fully hydroxylated tridymite surface. Owing to polar hydroxyl groups the water nanodroplets almost diminished and the water molecules were spilled over the surface forming a thin film meaning that the pure (001) tridymite surface is strongly hydrophilic. It was observed formation of hydrogen bonds between water molecules and the hydroxyl groups from the surface and also between water molecules. No visible pattern of the orientation of water molecules was observed as in case of the octahedral surface of kaolinite [14]. Similar strong hydrophilic behavior was also observed for a fully hydroxylated amorphous silica surface [9, 24]. The one-dimensional z-density profile of water molecules (black dashed curve in Fig. 3; silicon atoms represent a reference plane) shows that molecules are distributed on the surface not only as an incomplete monomolecular layer but part of the water molecules also form a second water layer (see shoulder at ∼6 Å in Fig. 3). From the calculated average density of the water droplet it was only possible to estimate that the contact angle is below 5°. Both methods of the CA measurement of the soda-lime glass slides provide very similar CA of 18.4° (SDM) and 19.0° (WPM), respectively (Table 1). These values, although higher than the estimated CA from the MD simulations, show that the glass surfaces are hydrophilic due to surface OH groups. However, it has to be noted that the CA of the pure glass surfaces can strongly vary in time due to an aging of the surface and adsorption of contaminants from air [25]. Further, a difference between the measured and calculated CA can be explained by a lower density of OH groups on the soda-lime glass surfaces than on the regularly fully covered tridymite surface (4.6 OH groups per nm²). Leed and Pantano estimated the maximal density of the OH groups on the amorphous silica glass surface can be about 3.9 OH/nm² [26].

The MD simulation on the water droplet interacting with the partial coating (DMS_50 model) finished with a very similar result as achieved for the pure tridymite surface. The last snapshot collected after 4 ns MD run is displayed in Fig. 2b. The water molecules are spilled over the surface forming a thin incomplete water film. Thus, this surface is hydrophilic because water has access to two types of hydroxyl groups – surface –OH groups and hydroxyl groups from the DMS units that are not bound to the surface via the Si–O–Si bridge. It means that the overall surface density is the same as in case of the pure tridymite surface. It was again practically impossible to calculate the contact angle from the averaged water density profile and only an estimated value is given in Table 1. Figure 3 shows the z-density profile of the water molecules for the DMS_50 model (red dashed curve). The shape of the curve is different compared to the z-density profile of the pure tridymite surface. Although both curves have very similar width, three peaks are evident on the z-density profile of the DMS_50 model. First two peaks represent water molecules penetrating to the space between siloxane units forming multiple hydrogen bonds with both types of the polar hydroxyl groups. The third, largest peak represents the rest of water molecules that only weakly interact with the nonpolar methyl groups.

Fig. 2:

Front and top view of snapshots after 2 ns from MD simulation of water droplet (Nw = 1000) on (a) tridymite, (b) DMS_50, (c) DMS_100, (d) DMS_bi, and (e) DMS_2.

Fig. 3:

The z-density profiles ρ(z) of water oxygen atoms for water droplet from last 2 ns of MD simulation for all five models.

Table 1

Experimental and calculated contact angles (θ) in °. Experimental contact angles are given with the standard deviation.

	Glass slide		DMDESa				DCDMSa
Exp. SDM	18.4 ± 0.9		45.1 ± 5.1				93.9 ± 3.2
Exp. WPM	19.0 ± 1.8		31.6 ± 2.7				99.4 ± 3.8
Calculated	Tridymite	DMS_50		DMS_100	DMS_bi	DMS_2
	<5	<5		55.7	62.2	78.0

^aDMDES, prepared from dimethyldiethoxysilane; DCDMS prepared from dichlorodimethylsilane.

The situation completely changed for the models with the full coverage of the surface (DMS_100, DMS_bi, and DMS_2). The corresponding snapshots from the MD simulations are shown in Figs. 2c–e. It is clearly seen that the full surface coverage has a significant influence on the shape of the water droplet. In case of the DMS_100 and DMS_bi models the water droplets form a deformed sphere while in case of the DMS_2 model the water droplet keeps almost hemispherical shape. Figure 4 displays the averaged density profiles for all three models with the full coverage. It is evident that the profiles are similar what is reflected also in a similarity between z-density curves displayed in Fig. 3 (green (DMS_100), blue (DMS_bi), and red (DMS_2) in online version). The water density profiles are not homogeneous and show certain fluctuations near the surface. This is also confirmed by the z-density profiles with the peaks that correspond to domains with the higher water densities. The larger fluctuation observed for the DMS_100 model is a consequence of the presence of the hydroxyl groups from the DMS units. Although they are partially screened by the nonpolar methyl groups, water molecules from the droplet still feel them. Therefore, the calculated contact angle (55.7°, Table 1) indicates that this surface is partially hydrophilic. Very similar CA (62.2°, Table 1) was also obtained for the model with the bidentate binding of the DMS units to the surface (DMS_bi). However, in this case no hydroxyl groups, responsible for wettability, are available on the surface to form hydrogen bonding with water molecules. In this model, the DMS units are bound to the surface in a bidentate manner and these units become more rigid. Moreover, owing to the arrangement of the DMS units on the surface, a larger portion of the free silica surface then in case of two other models stays accessible for water molecules. Thus, a part of water molecules from the droplet penetrates close to the silica surface as it is evidenced by a sharp peak below 5 Å in the z-density profile (blue line in Fig. 3). Consequently, the surface of the DMS_bi model has similar wettability as the DMS_100 model.

Fig. 4:

Visualization of average water density profiles (in g/cm3) for water droplet from last 2 ns of MD simulation: (a) DMS_100, (b) DMS_bi, and (c) DMS_2.

Since the surface in the DMS_2 model is formed only from the dimerized DMS units, the water density fluctuations are less pronounced than in case of the DMS_100 model. The fluctuation in the water density nearby the surface comes mainly from a polarization effect of the tridymite slab. It still has a certain impact on a preferential orientation of the water molecules in the first water layer on the surface because the thickness of the layer of the nonpolar methyl groups is small and their screening effect is not complete. The calculated contact angle for the DMS_2 model of 78.0° shows that the surface fully covered with the dimerized DMS units is distinctly water repellent.

Table 1 collects also the CA results obtained experimentally for the silanized glass surfaces. Experiments provide two distinctly different values for two different precursors with the same (CH₃)₂Si(O–)₂ functionality. Both CAs (SDM and WPM) for the glass coated by using DMDES are relatively small (<50°) indicating an only slightly water repellent surface. In opposite, if the DCDMS precursor is used, the surface becomes hydrophobic (CA > 90°). Our computed CAs lay in between the measured CAs. The calculated CA for the DMS_100 and DMS_bi models are closer to the DMDES-coated surface, whereas the CA for the DMS_2 model is closer to the DCDMS-coated surface. Obtained CA values indicate that the application of CLAYFF force field is appropriate for the simulation of the CAs on siloxanes. This is in agreement with previous works presenting the CLAYFF as the best choice from the common force-fields for the simulation of the mineral–water interfaces (including SiO₂) [27–29]. It seems that two important factors have a crucial impact on the surface wettability – the percentage of the coating and lateral oligomerization, whereas a form of binding (i.e., mono- or bidentate) is less important. Therefore, for these two real chemical processes reaction kinetics will be a crucial factor. Generally, chlorosilanes are much more reactive than alkoxysilanes, so the hydrolysis in a first stage of coating will be a driving force for the different types of the surface coating. For example, Bierbaum et al. [6] showed that the preparation of the coated silicon surface using trichlorosilane or trimethoxysilane leads to a different surface arrangement. In case of the two precursors used in our work it is also not possible to exclude also more complicated reaction scenarios of coating, e.g., a branching and forming of –[O–Si(CH₃)₂–O–Si(CH₃)₂–O]– chains, especially in case of the more reactive DCDMS.

Conclusions

In this work we conducted molecular dynamics simulations aimed to explain a different wettability of silica-glass surfaces coated by the same organosiloxane moiety but prepared from two different precursors. For this purpose, the interactions of the water nanodroplets with several models of the tridymite surface coated by dimethylsiloxane units were simulated. The simulations showed that the hydroxylated tridymite (001) surface is fully hydrophilic and the same was observed for a partial coating by DMS. Only fully covered surfaces showed a water repellency character. It was shown that the surface hydrophillicity is directly related to the concentration of the hydroxyl groups and their accessibility for the water molecules. The decreasing amount of the hydroxyl groups leads to increasing hydrophobicity and, consequently, larger contact angle.

Based on the findings from the simulations it can be concluded that the difference in the experimental CAs for two glass surfaces prepared from two different precursors with the same functionality is mainly due to a different arrangement of the coating moieties. Using the more reactive dichlorodimethylsilane leads to a much more compact surface coverage than in case of using the dimethyldiethoxysilane precursor. Beside the replacement of the surface OH groups also a certain lateral oligomerization is supposed.