Anisotropic Hydrolysis Susceptibility in Deformed Polydimethylsiloxanes

compounds with steric shielding by the two methyl groups that also protect the backside. .... Either the chain remained intact with a free nearby wate...
0 downloads 0 Views 10MB Size
Subscriber access provided by Macquarie University

B: Fluid Interfaces, Colloids, Polymers, Soft Matter, Surfactants, and Glassy Materials

Anisotropic Hydrolysis Susceptibility in Deformed Polydimethylsiloxanes Matthew P. Kroonblawd, Nir Goldman, and James P. Lewicki J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.9b07159 • Publication Date (Web): 27 Aug 2019 Downloaded from pubs.acs.org on August 30, 2019

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Anisotropic Hydrolysis Susceptibility in Deformed Polydimethylsiloxanes Matthew P. Kroonblawd,∗,† Nir Goldman,†,‡ and James P. Lewicki† †Physical and Life Sciences Directorate, Lawrence Livermore National Laboratory, Livermore, CA 94550, United States ‡Department of Chemical Engineering, University of California, Davis, California 95616, United States E-mail: [email protected]

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract Chemical reactions involving the polydimethylsiloxane (PDMS) backbone can induce significant network rearrangements and ultimately degrade macro-scale mechanical properties of silicone components. Using two levels of quantum chemical theory, we identify a possible electronic driver for chemical susceptibility in strained PDMS chains and explore the complicated interplay between hydrolytic chain scissioning reactions, mechanical deformations of the backbone, water attack vector, and chain mobility. Density functional theory (DFT) calculations reveal that susceptibility to hydrolysis varies significantly with the vector for water attacks on silicon backbone atoms, which matches strain-induced anisotropic changes in the backbone electronic structure. Efficient semiempirical density functional tight binding (DFTB) calculations are shown to reproduce DFT predictions for select reaction pathways and facilitate more exhaustive explorations of configuration space. We show that concerted strains of the backbone must occur over at least few monomer units to significantly increase hydrolysis susceptibility. In addition, we observe that sustaining tension across multiple monomer lengths by constraining molecular degrees of freedom further enhances hydrolysis susceptibility, leading to barrierless scission reactions for less substantial backbone deformations than otherwise. We then compute chain scission probabilities as functions of the backbone degrees of freedom, revealing complicated configurational inter-dependencies that impact the likelihood for hydrolytic degradation. The trends identified in our study suggest simple physical descriptions for the synergistic coupling between local mechanical deformation and environmental moisture in hydrolytic degradation of silicones.

1

Introduction

Silicones are widely used to form rubber and foam components for technological applications due to their desirable shape-filling factors, tunable mechanical properties, and chemical inertness across a wide range of operational conditions. 1,2 The dominant constituent in most silicones is polydimethylsiloxane (PDMS), which consists of a Si-O backbone with two methyl 2

ACS Paragon Plus Environment

Page 2 of 32

Page 3 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

side groups per monomer. Specific macroscopic properties of a silicone are imparted through the formation of crosslinks, addition of fillers and/or dissoluble pore formers, copolymerization with other backbone units, controlling the curing processes, and more recently through imparting engineered microstructures with additive manufacturing. 3–8 Despite their high degree of chemical inertness, silicones are not impervious to chemical degradation processes. Competing chain scission and crosslinking reactions and larger scale chain motions can result in creep, permanent set, and mechanical failure, which can irreversibly degrade the performance of components and potentially render them inoperable. 9–16 Degradation effects are often exacerbated under the application of complicated (and often sustained) stress loads typical of some operating conditions. Understanding the drivers for these processes at an atomistic level can inform models for aging characteristics over the service lifetime and potentially guide remediation of degradative effects. 9,17–20 It is well-known that exposure to water, acids, and bases all promote chain scission reactions that can alter the larger polymer network in silicones and degrade macroscopic mechanical characteristics. 12,21–23 We focus here on hydrolysis reactions involving the siloxane backbone as both unfilled and filled PDMS polymers are permeable to moisture and can uptake non-negligible quantities of water at room temperature and modest relative humidities. 24,25 Filler materials such as silica greatly exacerbate moisture uptake in PDMS and are thought to serve as sites for pooling. 17,25 Gauging possible connections between mechanical loading, environmental moisture, and hydrolytic chain-scissioning reactions would provide key characteristics for understanding larger network rearrangements in a number of real-time and accelerated aging scenarios. Detailed reaction mechanisms and energetics for chemical degradation processes can be difficult to elucidate from experiments alone. Atomistic modeling based on electronic structure theory can predict these properties and provides a robust and accurate means to probe silicon chemistry. 20,26–32 The chemistry of polysiloxanes under tension has received some attention at various levels of quantum chemical theory. 20,26,28,29,31 Most of these studies have

3

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

focused on understanding the response of single or multiple PDMS chains in the absence of environmental gases. A consistent finding is that dry chains broken through purely mechanical means typically have highly reactive termini at the cleavage site, which can go on to induce further degradation. Lupton et al. 28 modeled an atomic force microscopy experiment 33 in which a siloxane was placed under uniform tension using static calculations and dynamics simulations. They noted that the critical strain required to induce a chain scission in PDMS was reduced in the presence of water, but did not uncover the physical underpinnings for this synergistic environmental-mechanochemical effect. Here we investigate the susceptibility of siloxanes to chain-scissioning hydrolysis reactions under a variety of mechanical loading scenarios, to elucidate its coupling to the presence of water and local polymer deformations. Two levels of quantum chemistry theory are used, including highly accurate and computationally expensive density functional theory 34 (DFT) and very efficient semiempirical density functional tight binding 35,36 (DFTB). The DFTB Hamiltonian is derived from DFT and offers significant gains in computational efficiency while retaining much of the accuracy of DFT through a combination of approximate quantum mechanics and empirical functions. Due to the large number of degrees of freedom that lead to many possible backbone configurations, assessing the chemical susceptibility of deformed siloxanes with a purely DFT-based approach can require hundreds of thousands of CPU hours. In contrast, DFTB calculations for these systems are roughly 400 times less computationally intensive and have substantially lower wall-clock times as they can be run in batches on individual threads. A distinct advantage of DFTB is that its balance between accuracy, computational efficiency, and low wall-clock requirements greatly reduce turnaround when developing approaches to exhaustively explore high-dimensional chemical spaces. Through strategic application of expensive DFT calculations and high throughput sampling with DFTB, we uncover a possible electronic driver for increased chemical susceptibility in strained PDMS chains and then develop a probabilistic model that links hydrolytic chainscissioning reactions and deformations of the backbone degrees of freedom. Susceptibility

4

ACS Paragon Plus Environment

Page 4 of 32

Page 5 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

PDMS Chain

DMPS (dodecamethylpentasiloxane)

Lz



Figure 1: Model systems of an infinite PDMS chain and gas-phase DMPS molecule. Atoms are colored yellow, red, cyan, and white for silicon, oxygen, carbon, and hydrogen. The periodic simulation cell for the infinite PDMS chain is drawn with green lines. Snapshots were rendered using VESTA. 37 to hydrolysis is found to vary significantly with the approach vector for water attacks on silicon backbone atoms, which matches strain-induced anisotropic changes in the backbone electronic structure. An off-the-shelf DFTB model for silicon-containing systems is shown to reproduce DFT-level predictions for dry PDMS chains under strain and for selected hydrolysis pathways. Exhaustive explorations of configuration space made practical with DFTB reveal that concerted strains of the backbone must occur over at least a few monomer units to significantly increase hydrolysis susceptibility. Chain scission probabilities are extracted as functions of the backbone degrees of freedom, revealing complicated configurational interdependencies that impact the likelihood for hydrolytic degradation.

2 2.1

Methods Model Siloxane Systems

We considered two model polysiloxane systems (see Figure 1), namely a periodically infinite PDMS chain with four monomers and the methyl-terminated molecule dodecamethylpentasiloxane (DMPS). These two “reductionist” models were chosen as they include enough monomer units to provide physicochemical insights on how chain-level deformations and

5

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

mobility influence specific local reactivity. While this approach was chosen in part to reduce the configuration space (and the subsequent computational cost) of exhaustively exploring chemical reactivity, it also allows us to isolate potential local, short-ranged order parameters that govern hydrolysis of these systems. The periodically infinite PDMS chain was used to study polysiloxanes in a state of uniform strain (equivalently uniform tension). The chain was oriented so that it was periodic in the z-direction and uniform strains were imposed by varying the cell length Lz along that direction. The transverse cell lengths were set to Lx = Ly = 8 Å to reduce non-bonded self-interactions so that the model system corresponds closely to an isolated polymer strand, consistent with our reductionist approach. Gas-phase calculations of the terminated molecule DMPS were used to investigate non-uniform deformations of the siloxane backbone that were otherwise difficult to access.

2.2

DFT Calculations

Density functional theory 34 (DFT) level optimizations of PDMS and DMPS configurations were performed with VASP 38 using the Perdew-Burke-Ernzerhof 39 (PBE) generalized gradient approximation functional with projector-augmented wave (PAW) potentials 40,41 and Grimme D3 dispersion corrections. 42 All DFT-level calculations were performed using 3D periodic simulation cells. The electronic structure was evaluated at the Γ-point only using a 500 eV plane wave cutoff and Fermi-Dirac thermal smearing 43 with the electron temperature set to 25.85 meV (i.e., 300 K). The self-consistent field accuracy threshold was set to 10−6 eV. Conjugate-gradient optimizations were peformed with a force-based threshold set to 2 ×10−2 eV Å−1 .

2.3

DFTB Calculations

High-throughput quantum-based optimizations of PDMS and DMPS configurations were performed using the self-consistent charge DFTB method 35,36 as implemented in the DFTB+ code. 44 The DFTB total energy is derived from an expansion of the Kohn-Sham energy 34 6

ACS Paragon Plus Environment

Page 6 of 32

Page 7 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

about a reference electronic density to second order in charge fluctuations and is evaluated as EDFTB = EBS + ECoul + ERep + EDisp .

(1)

Here EBS is the electronic band structure energy, ECoul captures electrostatic interactions between fluctuating partial charges, and ERep and EDisp are empirical repulsion and dispersion terms. We used the DFTB parameter set pbc-0-3 (available at http://www.dftb.org), a typical off-the-shelf parameter set for silicon-containing systems that we previously applied 20 to study dephenylation reactions in copolymers of polydiphenylsiloxane and polydimethylsiloxane. Universal force field dispersion terms 45 were used for EDisp . DFTB employs a tight-binding framework 46 to efficiently construct the electronic Hamiltonian and overlap integrals from tables in terms of a nonorthogonal, atomic orbital basis set. Calculations involving PDMS were performed using 3D periodic simulation cells and calculations of DMPS were performed in the gas phase. The electronic structure was evaluated without spin polarization at the Γ-point only using Fermi-Dirac thermal smearing 43 with the electronic temperature set equal to 25.85 meV (i.e., 300 K). The self-consistent field accuracy threshold was set to 10−6 eV. Conjugate-gradient optimizations were peformed with a force-based threshold set to 2 ×10−2 eV Å−1 .

2.4

Structure Analysis

To characterize how deformations of the siloxane backbone impact hydrolysis susceptibilities, we performed many thousands of optimizations that started from specific configurations of PDMS and DMPS with a nearby water molecule. The primary metric of interest was whether the siloxane backbone remained intact or underwent scission during optimization. The vast number of optimizations involved prevented manual identification of the resulting structures, so we automated the structural characterization process by using a graph-based analysis of covalent connectivity. Two atoms were considered bonded if their separation distance fell

7

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Table 1: Bonding criteria for automated structure recognition Interaction Type Bonding Cutoff (Å) Si-Si 2.1 Si-O 1.8 Si-C 2.1 Si-H 1.8 O-O 1.8 O-C 1.8 O-H 1.3 C-C 1.8 C-H 1.3 H-H 1.0 within prescribed cutoffs taken to be slightly larger than typical equilibrium bond lengths (see Table 1). For a given structure, a square symmetric adjacency matrix (Ai,j ) representing bonded connectivity was formed in which Ai,j = 1 if atoms i and j were bonded and Ai,j = 0 otherwise. The connected components of matrix (Ai,j ) correspond to discrete non-connected molecular entities (e.g., DMPS and H2 O), and can be used to numerically distinguish between different siloxane backbone compositions. We identified all connected components of (Ai,j ) using the opensource package NetworkX 47 and compared the components against a library of known structures to determine the composition of final post-optimization configurations.

3 3.1

Results and Discussion Electronic Structure of Strained Siloxanes

Initial characterization of a PDMS chain under uniform strain revealed that its electronic structure is sensitive to deformations of the Si-O backbone. Figure 2 shows DFT-level predictions for the energy gap between the highest-occupied and lowest-unoccupied molecular orbitals (i.e., the HOMO-LUMO gap, EGap ) of an infinite PDMS chain as a function of strain, molecular representations of the HOMO and LUMO states at particular strains, and the HOMO and LUMO states of a uniformly strained DMPS configuration. A series of

8

ACS Paragon Plus Environment

Page 8 of 32

Page 9 of 32

HOMO-LUMO Gap (eV)

(a)

Strain

HOMO

(b)

LUMO

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Strain 1.10

1.15

1.20

(c) HOMO

LUMO

Figure 2: DFT-level predictions for (a) the HOMO-LUMO gap of a strained infinite PDMS chain, (b) representative HOMO and LUMO states of the PDMS chain spanning the change in LUMO character with increasing strain, and (c) the HOMO and LUMO states of DMPS with similar backbone geometry to the strained infinite PDMS chain. Atoms are colored yellow, red, cyan, and white for silicon, oxygen, carbon, and hydrogen, the periodic simulation cells for the infinite PDMS chain in (b) are drawn to scale with green lines, and molecular orbitals are rendered as gray. Snapshots were rendered using VESTA. 37

9

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

unconstrained DFT-level optimizations were performed for the periodically infinite PDMS chain (see Figure 1). Strains were imposed by setting the cell length Lz along the backbone to fixed values in 8 Å ≤ Lz ≤ 14 Å, which brackets the equilibrium chain length and is not sufficiently large to induce a chain scission. A minimum in energy occurs near Lz,0 = 10 Å, which we use to define strains of the chain as  = Lz /Lz,0 . We also performed a constrained optimization of DMPS in which the Si-O backbone was set to have the same bond lengths and angles as the infinite PDMS chain at  = 1.34 (when EGap = 3.0 eV) in order to mimic uniform strain. The starting DMPS backbone configuration was obtained using a classical FF (discussed below) and the constrained DFT-level optimization was performed while holding the Si and O atoms fixed in a large cubic cell with a side length of 15 Å. Two features in the electronic structure of the infinite PDMS chain are immediately apparent. The first feature is that EGap decreases substantially as the PDMS chain is strained beyond 1.10. The second feature is that the character of the LUMO state changes under strain near the inflection point of the EGap () plot. At strains of 1.20 and larger, LUMO density localizes to the interior of the O-Si-O bond angles along the backbone. LUMO localization coupled with a substantial reduction in EGap indicates that the interior O-Si-O bond angles may be potential sites for chemistry such as a hydrolytic attack by environmental water or by roaming acids or bases. This is in contrast to an attack on the backside of the OSi-O angle, where there is practically no density of low-energy unoccupied states. A distinct possibility, which we explore below, is that strained PDMS chains have an electronically driven anisotropy in their susceptibility to hydrolysis. This electronic susceptibility likely compounds with steric shielding by the two methyl groups that also protect the backside. Similar LUMO localization occurs in strained DMPS, albeit with a slightly larger EGap (3.6 eV) compared to the infinite PDMS chain with the same backbone geometry (3.0 eV). These electronic similarities indicate that increased chemical susceptibility may be a highly localized (or localizable) phenomena and that DMPS serves as a reasonable gas-phase analog to PDMS.

10

ACS Paragon Plus Environment

Page 10 of 32

Page 11 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Me

Si O1

O2 θ

x

y

Figure 3: Spherical coordinate system for vectored attacks on silicon. The azimuthal plane for the coordinate system is shown as a gray transparent circle and the azimuthal angle θ is defined with respect to rSi−O1 k x. Atoms are colored yellow, red, cyan, and white for silicon, oxygen, carbon, and hydrogen. Snapshots were rendered using VESTA. 37

3.2

PDMS Chain Hydrolysis

Potential anisotropy in the susceptibility of PDMS chains to hydrolysis was investigated at both the DFT and DFTB levels of theory. Since strain-induced changes to the electronic structure suggest that the interior of O-Si-O bond angles are potential sites for chemistry, we defined a silicon-centered coordinate system to analyze vectored attacks by water. This coordinate system, shown in Figure 3, is defined as a set of spherical coordinates whose azimuthal plane is given by rSi−O1 × rSi−O2 . Here, rSi−O1 and rSi−O2 are the separation vectors between the silicon atom center Si and the two bonded oxygen atoms O1 and O2 . We focus on vectored attacks by water within the azimuthal plane at prescribed azimuthal angle θ. By definition, oxygen O1 is located at θ = 0o , oxygen O2 lies within the azimuthal plane at some non-zero θ < 180o , and the methyl groups bonded to the silicon center lie out of the plane with θ > 220o . We first considered two vectors for attack by water on a silicon atom in uniformly strained PDMS chains. One vector probed hydrolysis susceptibility within the pocket formed by the O-Si-O bond angle along θ = 40o and the other probed susceptibility to attack on the backside along θ = 240o (approximately between the methyl groups). Two sets of analogous calculations were performed, with one using DFT and the other DFTB. For a given set of

11

ACS Paragon Plus Environment

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy ( kcal mol-1 monomer-1)

The Journal of Physical Chemistry

Page 12 of 32

50 DFTB (θ = 40 Deg.) DFTB (θ = 240 Deg.) DFTB (Dry) DFT (θ = 40 Deg.) DFT (θ = 240 Deg.) DFT (Dry)

40 30 20 10 0 1.0

1.1

1.2

1.3

1.4

1.5

Strain

Figure 4: Comparison of DFT and DFTB energies for optimized uniformly strained PDMS chains. Curves for the dry chain are plotted up to the highest strain that does not induce scission. Results for vectored attacks by water are shown using circles with green lines and diamonds with orange lines for attack angles 40o and 240o , respectively. calculations, we first obtained optimized dry chain configurations for strains  starting at 1.00 in increments of 0.02. Both DFT and DFTB predict similar equilibrium Lz,0 , which we set to 10 Å. The dry chain scissions at very high strains during optimization, so the upper limit was taken to be the largest strain for which the chain remains intact following optimization. Starting from the optimized dry chain configurations, we spawned two separate series of calculations in which a single water molecule was placed at r = 2.0 Å and either θ = 40o or 240o in the coordinate system described above. The water molecule was oriented so that its oxygen atom was at r = 2.0 Å with the hydrogen atoms pointed away from the silicon center at positions that were orthogonal to the azimuthal plane. These hydrated configurations were optimized without constraints, which resulted in only two distinct final states across all of the calculations. Either the chain remained intact with a free nearby water molecule or hydrolysis occurred resulting in a scissioned chain with Si(OH)(CH3 )2 terminal ends. Figure 4 collects the combined energy results from the six series of optimizations just described. The reference energy for a given curve is set to be the minimum in energy Emin () for that curve over the sampled strain interval, which is either at or very near  = 1.00 for all cases. Focusing first on the results for the dry chain, it is clear that DFT and DFTB predict very similar strain energies (compare solid and dashed black lines). We take the critical strain

12

ACS Paragon Plus Environment

Page 13 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

threshold ∗ to be the largest strain for which a barrierless scission reaction does not occur, and the corresponding critical strain energy to be E(∗ ) − Emin (). The only difference of note is that DFTB predicts a slightly higher critical strain (1.46) required for a dry chain scission than DFT (1.40). The predicted strain energies per monomer at the breaking point are 41.2 and 48.8 kcal mol−1 monomer−1 for DFT and DFTB, respectively. (Note that per-bond quantities are half the per-monomer value.) These critical strain energies are similar in magnitude to results obtained from semiempirical PM3-based 48 quantum chemistry optimizations of gas-phase dry PDMS analogs under tension. 26 Dynamical effects would almost certainly reduce the apparent critical strains at finite temperatures. 28 Abrupt decreases in the curves for the hydrated configurations coincide with the critical strain for barrierless hydrolytic chain scission, with the values beyond that drop corresponding to the energy of a scissioned chain with Si(OH)(CH3 )2 terminal ends. Anisotropy in hydrolysis susceptibility is immediately apparent. Both DFT and DFTB predict that substantially lower strains are required for barrierless hydrolysis along θ = 40o than along θ = 240o . DFT predicts barrierless hydrolysis along θ = 40o for strains over 1.24, which corresponds to a strain energy of 14 kcal mol−1 monomer−1 . This is in contrast to the DFT result for θ = 240o , which exhibits a larger critical strain (1.36) and corresponding input strain energy (35 kcal mol−1 monomer−1 ). Analogous but larger results were obtained with DFTB, with critical strains (strain energies) of 1.32 (23 kcal mol−1 monomer−1 ) and 1.40 (40 kcal mol−1 monomer−1 ) for the 40o and 240o cases, respectively. Both levels of theory also predict very similar energies for scissioned chains with Si(OH)(CH3 )2 terminal ends. Based on these results, we anticipate that DFTB will predict similar hydrolysis trends as DFT for other geometries with critical strains and strain energies that are modest overestimates. Using DFTB, we were able to extend our anisotropic hydrolysis analysis to cover the full domain of the azimuthal angle. We performed a series of DFTB-level optimizations following the same procedure used to generate Figure 4, but with finer resolution in the strain  and azimuthal angle θ. The strain interval was set to 1.00 ≤  ≤ 1.46 in 0.01 increments and the

13

ACS Paragon Plus Environment

The Journal of Physical Chemistry

Energy (kcal mol-1 monomer-1)

(a) 0

10

20

30

40

1.45

50

*

1.40

**

1.35 1.30

Strain

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 32

1.25 1.20 1.15 1.10 1.05

O1

1.00 0

O2 45

90

135

Me 180

225

O1 270

315

360

Azimuthal Angle (Deg.) (b)

*

**

1.88 Å

Figure 5: (a) Optimized energy surface for a PDMS chain under uniform strain and vectored attack by water at prescribed azimuthal angles. Positions of selected atoms in the initial optimized dry chain configuration are plotted as magenta lines, including oxygen atoms O1 and O2 used to define the azimuthal plane and the projected positions of the methyl (Me) carbon atoms (see Figure 3). Regions with sterically excluded configurations are rendered as black. (b) Snapshots of specific hydrolysis products denoted by ? and ?? rendered using OVITO 49 with the same atom coloring scheme used previously. Regions in which products ? and ?? were found are denoted in (a). azimuthal angle interval was set to 0o ≤ θ ≤ 360o in 2o increments. Hydrated configurations with an interatomic separation distance less than 1.2 Å were excluded from consideration to avoid sterically unfavorable geometries. The results from this analysis, which involved 3834 separate optimizations, are shown in Figure 5(a). Two distinct products that form following a chain scission event are shown in Figure 5(b), which we denote by ? and ??. Focusing first on low θ values between atoms O1 and O2 , we predict that there is a single critical strain (1.32) beyond which attack by water leads to formation of a favorable hydrolysis product ?? that corresponds to a scissioned chain with Si(OH)(CH3 )2 terminal ends. Insensitivity to the specific angle of attack, provided it is within the O-Si-O bond angle, is consistent with hydrolysis susceptibility being driven by an abrupt strain-induced change 14

ACS Paragon Plus Environment

Page 15 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

to the local electronic structure in that region. Product ?? also manifests for other θ at larger strains including just outside the O-Si-O bond angle (near 180o ) and for vectors right between the methyl groups (near 240o ). Another unfavorable product also forms for very large strains in 270o ≤ θ ≤ 315o , which we denote by ?. One terminal end is an Si(CH3 )2 O group and the other a Si(CH3 )2 group that separated by 1.88 Å from the adsorbed attacking water molecule. While product ? is lower in energy than the strained chain, it is nearly 100 kcal mol−1 (or 25 kcal mol−1 monomer−1 ) higher in energy than product ??. One possible interpretation is that the presence of water can destabilize highly strained chains to induce scissions, even if a high-θ attack vector does not leave a viable downhill path to the hydrolysis product ??. A distinct possibility is that ? and ?? may both be intermediate species at finite temperatures, with the cyclic species octamethylcyclotetrasiloxane being a reasonable final degradation product obtained though a chain-backbiting-type reaction. 23 Our predicted reductions in the critical strain required for hydrolytic chain scission verses dry scission qualtitively match the results of Lupton et al. 28 In their study, quantum-based molecular dynamics simulations of a very rapidly strained terminated gas-phase PDMS chain revealed that scission occurs at a relative extension that is ≈9% less than the dry-chain maximum. Although the specific vector for water attack was not reported, their result would be consistent with practically all of the attack vectors considered here. The mechanically driven changes in PDMS electronic structure uncovered here, namely a simultaneous reduction in band gap and localization of the LUMO state with increasing tension, could explain the synergistic coupling between water and tension reported by Lupton et al.

3.3

DMPS Hydrolysis

Real silicones have highly complicated microstructures in which individual PDMS chains might adopt unusual geometries owing to their significant conformational flexibility. Increased susceptibility to hydrolysis of PDMS chains under uniform strain was established above, but it remains unclear whether other deformations to the siloxane backbone might 15

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

also facilitate hydrolysis chemistry. Under uniform strains, a PDMS chain at T = 0 K adopts specific configurations that span only a small subspace of the possible values for the three internal coordinate types that describe the backbone, namely the Si-O bond lengths, the O-Si-O angles, and the Si-O-Si angles. To better gauge hydrolysis susceptibility within a larger configuration subspace, we focus here on hydrolysis of DMPS. Working with a gas-phase molecule allows for substantially greater control of the backbone internal coordinates than is possible with a periodically infinite chain. In all calculations described below, we considered a water attack on the central silicon atom in DMPS along a vector oriented at θ = 40o with r = 2.0 Å in the azimuthal plane defined using the coordinate system introduced in Figure 3. This particular angle for water attack was chosen based on our results for the infinite PDMS chain that showed attacks within the O1 -Si-O2 bond angle led to hydrolysis for lower applied strains than for attacks on the methylated backside. Positioning of the water molecule next to pre-optimized dry DMPS configurations followed the same procedure used in Section 3.2 for PDMS. Based on the localized nature of the LUMO state in uniformly strained PDMS, a reasonable hypothesis is that only three internal coordinates strongly influence hydrolysis susceptibility, namely the bond lengths Si-O1 and Si-O2 and the bond angle O1 -Si-O2 for a silicon under attack by water. We assessed the susceptibility to hydrolysis reactions in DMPS with prescribed values for the above internal coordinates. The bond lengths Si-O1 and Si-O2 were kept equal and were set to values in the interval 1.60 Å ≤ rSi−O ≤ 2.40 Å with 0.05 Å increments. The bond angle O1 -Si-O2 was varied independently of rSi−O with values in the interval 90o ≤ θO−Si−O ≤ 170o with 1.0o increments. For reference, DFTB predicts equilibrium values of 1.64 Å and 113o for the Si-O bond length and O-Si-O angle in DMPS. A total of 1440 dry DMPS configurations were generated with the above-prescribed positions for the O1 , Si, and O2 atoms. Each configuration was optimized while holding those three atoms fixed, after which we positioned an attacking water molecule and then optimized without constraints. The final optimized structures were analyzed using the graph-based approach

16

ACS Paragon Plus Environment

Page 16 of 32

Page 17 of 32

Scission? No

Yes

2.4

Si-O Bond Length (Å)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

2.3 2.2 2.1

O2

2.0 1.9 1.8

O1

1.7 1.6 90

100

110

120

130

140

150

160

170

O-Si-O Angle (Deg.)

Figure 6: Probability map for barrierless hydrolytic scission of DMPS as a function of the two equal Si-O bond lengths Si-O1 and Si-O2 and the single O1 -Si-O2 angle. For comparison, a parametric plot of the O1 and O2 atom positions in the uniformly strained infinite PDMS chain are shown as purple and orange curves, respectively. Dashed and solid portions of those curves differentiate between values for which barrierless hydrolysis will and will not occur in uniformly strained PDMS. described in Section 2.4 to determine whether a hydrolytic scission of the siloxane backbone occurred. Results from this analysis are shown in Figure 6. Regions in Si-O bond and O-Si-O angle space where non-uniformly deformed DMPS exhibits barrierless scission are considerably different than those where the uniformly strained PDMS chain will scission. This could indicate a degree of non-locality and/or non-uniformity in connections between a structural order parameter and the potential for reactivity. Such a “reactive order parameter” would be useful to assess the formation of likely reaction sites in larger scale simulations of siloxane-based materials under complicated loading scenarios. At very high strains, there is a breakage of symmetry in the Si-O1 and Si-O2 bond lengths in uniformly strained PDMS, which is seen in the diverging dashed purple and orange curves. The minimum Si-O1 and Si-O2 bond lengths required for barrierless scission in DMPS are 1.95 Å, which occurs at the smallest value of θO−Si−O considered. This corresponds to a compressed bond angle with highly elongated bond lengths relative to the equilibrium configuration. The bond lengths required for scission steadily increase with increasing O-

17

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Si-O angle. A small region of stability is seen near the equilibrium value for θO−Si−O at rSi−O ≈ 2.2 Å. At least two factors contribute to the differences between PDMS and DMPS scission probabilities seen in Figure 6. The first difference is that the deformations of DMPS were non-uniform. In uniformly strained PDMS, all of the Si-O bond lengths and O-Si-O bond angles in the chain adopt deformed periodically equivalent geometries, not just those near a silicon atom under attack by water. The second difference is that the gas-phase DMPS molecule cannot support a tensile state, unlike the periodically infinite PDMS chain. We explored the separate influences of uniform deformation and tension using DMPS as a model system. An initial DMPS backbone configuration was prepared using an FF-based approach where configurations for particular optimization series were obtained using a classical united atom force field (FF) for siloxanes that we implemented in LAMMPS. 50 Using a FF offers a convenient route to prepare initial backbone configurations due to its extremely high computational efficiency as well as its definition in terms of the same internal coordinates targeted in our conformational scans. As a result, a desired conformation can be easily obtained by modifying FF parameters (e.g., coordinate equilibrium values) and performing a preliminary computationally inexpensive FF-level optimization, followed by subsequent higher-level DFTB optimization. The siloxane FF used here combined the covalent bond, angle, and dihedral terms from the hybrid/UA model of Frischknect and Curro 51 with the Lennard-Jones non-bonded potentials from Thol et al. 52 Two-center bond and three-center angle potentials were modeled as harmonic oscillators and particular backbone geometries were obtained by performing FF-level optimizations in which the FF equilibrium values for the Si-O bonds, Si-O-Si angles, and O-Si-O angles were set to prescribed values with artificially increased harmonic force constants. We used harmonic force constants 2000 kcal mol−1 Å−2 and 2000 kcal mol−1 deg−2 for the backbone bonds and angles, respectively. Following FF-level optimization, we hydrogenated the carbon atoms and performed constrained optimizations at higher levels of theory in which the Si and O atoms in the backbone were held fixed.

18

ACS Paragon Plus Environment

Page 18 of 32

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Free Optimization (No Scission)

Initial

Constrained Optimization (Scission)

1.5

H2O (θ = 40o) 0.0

Displacement Magnitude (Å)

Page 19 of 32

Figure 7: Snapshots of an initial uniformly deformed DMPS molecule under attack by water and a comparison of the final structures that result from an unconstrained optimization and a constrained optimization in which the terminal silicon atoms were held fixed. Atoms are colored yellow, red, cyan, and white for silicon, oxygen, carbon, and hydrogen. Color-coded displacement vectors highlight the movement of atoms away from the initial configuration after optimization. Snapshots were rendered using OVITO. 49 In this work, Si-O bonds, O-Si-O angles, and Si-O-Si angles were uniformly deformed to values near those at the critical strain in the infinite PDMS chain (that is, near the beginning of the dashed curves in Figure 6). These internal coordinates were respectively set to 1.9 Å, 141o , and 179o . The actual critical Si-O bond length in PDMS was 1.8 Å, which was slightly lower than that in the present example. After completion of an anhydrous DFTB-level optimization, we placed an attacking water molecule near the central silicon atom at θ = 40o . Two subsequent optimizations were performed, with one being unconstrained and the other holding the two terminal silicon atoms in fixed positions. Constraining the terminal silicon atoms allows for local relaxations under a supported strain state more typical of a real (condensed-phase) PDMS system. Real silicones exhibit complicated crosslinked structures that could provide “pinning sites” such as those defined here to support local uniform tensile states. Dynamical effects that reduce chain mobilities, such as chain entanglements, might also produce local tensile states that are supported on chemical reaction time scales. Treating the terminal silicon atoms as either unconstrained or fixed corresponds to two extreme bracketing cases for chain mobility. Figure 7 shows the initial and final structures resulting from the unconstrained and

19

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

constrained optimization of uniformly deformed DMPS under attack by water. Optimization without constraints resulted in free molecules of DMPS and water whereas the constrained optimization led to clean hydrolytic scission of DMPS into two molecules with electronically stable terminal ends. The influence of pinning points that maintain a tensile state is obvious. Displacement vectors reveal that during a free optimization, the entire DMPS molecule contracts to partially relieve localized strains near the central silicon atom under attack by water. In contrast, when the ends are fixed, there is an notable relaxation towards those ends that pulls at the central atom. While not shown, these separate types of concerted motions manifest within the first few dozen optimization steps and likely lock-in a particular fate for the chain. It is perhaps worth noting that both types of optimization lead to an unscissioned state when the initial Si-O bond lengths are set to 1.8 Å, the actual critical value determined for the infinite PDMS chain. This is possibly due to the fact that pinning the terminal silicon atoms is not strictly equivalent to supporting tension through a fixed periodic boundary. Larger gas-phase analogs to PDMS than the DMPS molecule considered here might also better converge to the infinite chain results. One interpretation of the above results is that hydrolysis susceptibility increases due to semi-localized uniform deformations of the siloxane backbone that occur on length scales of a few monomers or more. Deformations of polymer backbones in real silicones are not necessarily uniform. A chain under a complicated stress load could adopt any number of geometries, many of which might increase the susceptibility of the chain to degradative chemical reactions. As a final exploration of hydrolysis susceptibility in deformed DMPS, we extended the analysis presented in Figure 7 to encompass a wider range of possible values for the Si-O bonds, O-Si-O angles, and Si-O-Si angles. From this analysis, we extract probability distributions for hydrolytic scission to obtain a reactive order parameter for hydrolysis susceptibility that is a function of the backbone degrees of freedom. All of the Si-O bonds, O-Si-O angles, and Si-O-Si angles were set to the same value in a given configuration, which reduces the phase space dimensionality from the 15 individual internal coordinates to a subspace that captures the three

20

ACS Paragon Plus Environment

Page 20 of 32

Page 21 of 32

Free

Constrained 190o

190o

Si-O-Si Angle

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

2.5 Å

2.5 Å 1.5 Å

190o

80o 80o

1.5 Å

190o

80o 80o

Si-O Bond O-Si-O Angle

Figure 8: Volume representations for configurations of uniformly deformed DMPS under attack by water where barrierless scission reactions occur during free and constrained optimization. Coordinate values within the volumes correspond to a scission probability of one. The volumes were rendered with the surface interpolation tools packaged in OVITO 49 using an atomic-style representation of binary scission probabilities defined on a rectangular Cartesian grid mapped to (rSi−O , θO−Si−O , θSi−O−Si ) space. coordinate types. This three-dimensional internal coordinate space was scanned on the intervals 1.6 Å ≤ rSi−O ≤ 2.4 Å in 0.1 Å increments, 90o ≤ θO−Si−O ≤ 180o in 5o increments, and 90o ≤ θSi−O−Si ≤ 180o in 5o increments, which corresponds to 3249 unique configurations. Similar to before, we performed two sets of optimizations with one being unconstrained and the other fixing the positions of the two terminal silicon atoms. Post-optimization configurations were analyzed using the graph-based approach described in Section 2.4 to determine which DMPS chains underwent scission reactions. The full volumetric representations shown in Figure 8 reveal a complicated topology for the bounding surface that delineates between a scission probability of zero and one. It is clear that the scission probability cannot be reduced to a simple function, much less guessed a priori, and that a practically usable reactive order parameter based on these probabilities would be best expressed in a tabulated form. Cursory inspection of the volumes reveals a number of similarities between the two cases for chain mobility. Supporting a tensile state in the chain through pinning its ends leads to the same hydrolysis susceptibility as a freely mobile chain for the vast majority of configurations explored in our search. To facilitate quantitative comparisons between the two sets of calculations, we generated 2D heat maps and 1D line plots of scission probability by averaging over one or two of the 21

ACS Paragon Plus Environment

The Journal of Physical Chemistry

Constrained

150 60 135 40 120 20

105

(a) 90 165

40 120 20

105

(d)

100

Si-O Bond (Å)

80 2.2 2.1

60

2.0 40

1.9 1.8

20

(b)

1.6 120

135

150

165

2.1

60

2.0 40

1.9 1.8

20 1.7

(e)

1.6

60

2.0 40

1.9 1.8

20

(c)

1.6

0 150

165

Si-O-Si Angle (Deg.)

120

135

150

180

165

80 2.2 2.1

60

2.0 40

1.9 1.8

20 1.7

(f)

1.6

0 120

135

150

2.4

60

40

20

165

180

(h) 105

120

135

150

165

180

Si-O-Si Angle (Deg.) 100

105

2.2

Free Constrained

80

90

180

2.3

90

2.0

0

0 105

2.4

Si-O Bond (Å)

2.1

Scission Probability (%)

80 2.2

135

(g) 1.8

100

O-Si-O Angle (Deg.) 100

120

20

Si-O Bond (Å)

80

90

2.3

105

40

1.6

180

2.2

180

2.4

90

165

100

O-Si-O Angle (Deg.)

1.7

150

2.3

0 105

135

2.4

Si-O Bond (Å)

2.3

Scission Probability (%)

2.4

90

120

60

Si-O-Si Angle (Deg.)

Si-O-Si Angle (Deg.)

1.7

105

Free Constrained

80

0

0 90

180

Scission Probability (%)

150

60 135

100

100

Scission Probability (%)

135

150

Scission Probability (%)

120

80

Scission Probability (%)

105

165

90

0 90

100

Scission Probability (%)

80

180

O-Si-O Angle (Deg.)

165

Scission Probability (%)

O-Si-O Angle (Deg.)

100

Scission Probability (%)

Free 180

Si-O Bond (Å)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 32

Free Constrained

80

60

40

20

(i)

0 90

Si-O-Si Angle (Deg.)

105

120

135

150

165

180

O-Si-O Angle (Deg.)

Figure 9: Comparison of different averages taken over the 3D volumes presented in Figure 8. Averages were taken over all the internal coordinates (rSi−O , θO−Si−O , and/or θSi−O−Si ) that are not plotted in a given panel. Columns one and two show 2D heat maps of scission probability as a function of two coordinates that were averaged over the third coordinate for the free and constrained sets of calculations. Contours for the heat maps are plotted in 20% intervals. Column three shows a comparison of free and constrained scission probabilities as a function of one coordinate that were obtained by averaging over the other two coordinates. internal coordinates. These averaged results are captured in Figure 9. Focusing first on the results presented in the first row (panels (a), (d), and (g)) reveals three key features. The first is that there are some combinations of the two backbone angles, such as near (θO−Si−O = 165o , θSi−O−Si = 135o ), for which barrierless hydrolysis occurs even when the rSi−O is near its equilibrium value. That is, there are regions in panels

22

ACS Paragon Plus Environment

Page 23 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(a) and (d) where the scission probability is 100%. These regions are somewhat larger in the constrained case compared to the free optimization case. Approximately 30% of the possible combinations of the backbone angles lead to scissions at rSi−O = 1.6 Å, regardless of whether a tensile state is supported through pinning. Panel (g) shows that there is a strong correlation between increases in rSi−O with increasing likelihood for scission events. Pinning in the constrained optimizations leads to a modest increase in susceptibility that manifests most strongly when the Si-O bond lengths are close to equilibrium values. Thus, pinning may still play a considerable role in determining hydrolysis susceptibilities for realistic backbone configurations. Row two of Figure 9 reveals that the scission probability is largely insensitive to the backbone angle θSi−O−Si . It is perhaps surprising that the curves in panel (h) exhibit a subtle maximum in scission probability when θSi−O−Si is near its equilibrium value of 144o . However, there is not a strong correlation between scission probability and θSi−O−Si , especially compared to coordinate rSi−O . Pinning leads to only modest differences in susceptibility that are greatest for large θSi−O−Si , which is also captured by the significant similarities in surface topology shown in panels (b) and (e). The insensitivity to θSi−O−Si is perhaps unexpected as strain-induced changes to the electronic structure localize instead near the O-Si-O angles. Analyzing row three of Figure 9 shows modest sensitivity in hydrolysis susceptibility to the angle θO−Si−O , with hydrolysis being most likely when rSi−O is large. Unlike θSi−O−Si , there is a minimum in scission probability when when θO−Si−O is near its equilibrium value of 113o . The likelihood for scission is highest when θO−Si−O approaches 180o , which is consistent with the changes in configuration and electronic structure that occur when PDMS is uniformly strained.

23

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

4

Conclusions

Hydrolysis reactions inducing chain scission of the Si-O polydimethylsiloxane (PDMS) backbone were investigated using quantum chemistry calculations performed at two levels of theory. Density functional theory (DFT) calculations revealed that deformation of the backbone under tensile loading reduced the energy gap between the highest-occupied and lowestunoccupied molecular orbitals (the HOMO-LUMO gap) while LUMO density localized to silicon sites in the interior of the O-Si-O backbone angles. Water significantly reduced the critical strain required to induce chain scissions and the susceptibility of strained chains to hydrolysis was found to vary considerably with the water attack vector. The lowest critical strain required for barrierless hydrolysis was found to occur for an attack vector that coincided with the localized low-energy LUMO state that developed under tension. High-throughput semiempirical density functional tight binding (DFTB) calculations were coupled with an automated graph-based structure recognition protocol to more exhaustively explore links between hydrolysis susceptibility, water attack vector, and deformations of particular backbone degrees of freedom. Chain scissions in uniformly deformed backbones resulting from most water attack vectors corresponded to a hydrolysis product with Si(OH)(CH3 )2 terminal ends, although a narrow range of attack angles induced nonhydrolytic scissions with unstable termini. A gas-phase analog to PDMS, namely dodecamethylpentasiloxane (DMPS), was used as a reductionist model system to isolate effects of changes in specific backbone internal coordinates on chemical susceptibility. Using DMPS, we probed hydrolysis of siloxanes with non-uniform backbone deformations and assessed the importance of chain mobility in the scissioning process. Sustaining tension across multiple monomers through pinning sites was found to increase hydrolysis susceptibility, leading to barrierless scission reactions for less substantial backbone deformations than were required to scission chains without constraints. Probabilities for barrierless scission were obtained as functions of three backbone degrees of freedom, including the Si-O bond length and the O-Si-O and Si-O-Si angles, to generate a “reactive order parameter” that links the backbone 24

ACS Paragon Plus Environment

Page 24 of 32

Page 25 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

geometry to changes in chemical susceptibility. The likelihood for chain scissions were most clearly correlated to deformations of the Si-O bond length and O-Si-O angle, consistent with the predicted strain-induced changes in electronic structure. The trends identified in our study suggest a plausible physical origin for synergistic coupling between mechanical deformation and environmental moisture in hydrolytic degradation of silicones. Probability functions for hydrolysis susceptibility can be used as inputs to help analyze and refine coarse-grained or constitutive modeling approaches for polymer material degradation. Our calculations could have a significant impact on aging studies of siloxane polymers, where the interplay between environmental moisture and complicated and/or sustained stress loads drives network rearrangements resulting in irreversible changes in material strength and compressive properties.

Acknowledgement The authors thank Rebecca K. Lindsey for guidance in selecting a classical force field for siloxanes. This work was performed under the auspices of the U.S. Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC52-07NA27344.

References (1) Arkles, B. Look what you can make out of silicones (biomedical applications of silicones). 1983, 13, 542–555. (2) Brooke, M. A. Silicon in organic, organometallic, and polymer chemistry; Wiley: New York, USA, 1999. (3) Boonstra, B. Role of particulate fillers in elastomer reinforcement: A review. Polymer 1979, 20, 691 – 704. 25

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(4) Sun, C.-C.; Mark, J. Comparisons among the reinforcing effects provided by various silica-based fillers in a siloxane elastomer. Polymer 1989, 30, 104 – 106. (5) Lewicki, J. P.; Harley, S. J.; Finnie, J. A.; Ashmore, M.; Bell, C.; Maxwell, R. S. In Progress in Silicones and Silicone-Modified Materials; Clarson, S. J., Owen, M. J., Smith, S. D., Van Dyke, M., Brook, M., Mabry, J., Eds.; American Chemical Society, 2013; Chapter 10, pp 133–154. (6) Lewicki, J. P.; Maxwell, R. S.; Mayer, B. P.; Maiti, A.; Harley, S. J. In Concise Encyclopedia of High Performance Silicones; Tiwari, A., Soucek, M. D., Eds.; Wiley, 2014; Chapter 11, pp 151–176. (7) Maiti, A.; Small, W.; Lewicki, J. P.; Weisgraber, T. H.; Duoss, E. B.; Chinn, S. C.; Pearson, M. A.; Spadaccini, C. M.; Maxwell, R. S.; Wilson, T. S. 3D printed cellular solid outperforms traditional stochastic foam in long-term mechanical response. Sci. Rep. 2016, 6, 24871. (8) Fatona, A.; Moran-Mirabal, J.; Brook, M. A. Controlling silicone networks using dithioacetal crosslinks. Polym. Chem. 2019, 10, 219–227. (9) Rottach, D. R.; Curro, J. G.; Grest, G. S.; Thompson, A. P. Effect of strain history on stress and permanent set in cross-linking networks: A molecular dynamics study. Macromol. 2004, 37, 5468–5473. (10) Chinn, S.; DeTeresa, S.; Sawvel, A.; Shields, A.; Balazs, B.; Maxwell, R. S. Chemical origins of permanent set in a peroxide cured filled silicone elastomer - tensile and 1H NMR analysis. Polym. Degrad. Stab. 2006, 91, 555 – 564. (11) Maiti, A.; Gee, R.; Weisgraber, T.; Chinn, S.; Maxwell, R. Constitutive modeling of radiation effects on the permanent set in a silicone elastomer. Polym. Degrad. Stab. 2008, 93, 2226 – 2229.

26

ACS Paragon Plus Environment

Page 26 of 32

Page 27 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(12) Patel, M.; Morrell, P.; Cunningham, J.; Khan, N.; Maxwell, R. S.; Chinn, S. C. Complexities associated with moisture in foamed polysiloxane composites. Polym. Degrad. Stab. 2008, 93, 513 – 519. (13) Patel, M.; Chinn, S.; Maxwell, R. S.; Wilson, T. S.; Birdsell, S. A. Compression set in gas-blown condensation-cured polysiloxane elastomers. Polym. Degrad. Stab. 2010, 95, 2499 – 2507. (14) Maiti, A.; Weisgraber, T.; Dinh, L. N.; Gee, R. H.; Wilson, T.; Chinn, S.; Maxwell, R. S. Controlled manipulation of elastomers with radiation: Insights from multiquantum nuclear-magnetic-resonance data and mechanical measurements. Phys. Rev. E: Stat., Nonlin., Soft Matter Phys. 2011, 83, 031802. (15) Maiti, A.; Weisgraber, T. H.; Gee, R. H.; Small, W.; Alviso, C. T.; Chinn, S. C.; Maxwell, R. S. Radiation-induced mechanical property changes in filled rubber. Phys. Rev. E: Stat., Nonlin., Soft Matter Phys. 2011, 83, 062801. (16) Clough, J. M.; Creton, C.; Craig, S. L.; Sijbesma, R. P. Covalent bond scission in the mullins effect of a filled elastomer: Real-time visualization with mechanoluminescence. Adv. Funct. Mater. 2016, 26, 9063–9074. (17) Gee, R. H.; Maxwell, R. S.; Balazs, B. Molecular dynamics studies on the effects of water speciation on interfacial structure and dynamics in silica-filled PDMS composites. Polymer 2004, 45, 3885 – 3891. (18) Chenoweth, K.; Cheung, S.; van Duin, A. C. T.; Goddard, W. A.; Kober, E. M. Simulations on the thermal decomposition of a poly(dimethylsiloxane) polymer using the ReaxFF reactive force field. J. Am. Chem. Soc. 2005, 127, 7192–7202. (19) Lacevic, N. M.; Maxwell, R. S.; Saab, A.; Gee, R. H. Molecular dynamics simulations of ordering of poly(dimethylsiloxane) under uniaxial stress. J. Phys. Chem. B 2006, 110, 3588–3594. 27

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(20) Kroonblawd, M. P.; Goldman, N.; Lewicki, J. P. Chemical degradation pathways in siloxane polymers following phenyl excitations. J. Phys. Chem. B 2018, 122, 12201– 12210. (21) Grassie, N.;

Macfarlane, I. The thermal degradation of polysiloxanes - I.

Poly(dimethylsiloxane). Eur. Polym. J. 1978, 14, 875 – 884. (22) Lewicki, J. P.; Albo, R. L.; Alviso, C. T.; Maxwell, R. S. Pyrolysis-gas chromatography/mass spectrometry for the forensic fingerprinting of silicone engineering elastomers. J. Anal. Appl. Pyrolysis 2013, 99, 85 – 91. (23) Lewicki, J. P.; Maxwell, R. S. In Concise Encyclopedia of High Performance Silicones; Tiwari, A., Soucek, M. D., Eds.; Wiley, 2014; Chapter 13, pp 191–210. (24) Harley, S. J.; Glascoe, E. A.; Maxwell, R. S. Thermodynamic study on dynamic water vapor sorption in Sylgard-184. J. Phys. Chem. B 2012, 116, 14183–14190. (25) Sun, Y.; Harley, S. J.; Glascoe, E. A. Modeling and uncertainty quantification of vapor sorption and diffusion in heterogeneous polymers. ChemPhysChem 2015, 16, 3072– 3083. (26) Nikitina, E. A.; Khavryutchenko, V. D.; Sheka, E. F.; Barthel, H.; Weis, J. Deformation of poly(dimethylsiloxane) oligomers under uniaxial tension: Quantum chemical view. J. Phys. Chem. A 1999, 103, 11355–11365. (27) Gibbs, G. V.; Jayatilaka, D.; Spackman, M. A.; Cox, D. F.; Rosso, K. M. Si-O bonded interactions in silicate crystals and molecules: A comparison. J. Phys. Chem. A 2006, 110, 12678–12683. (28) Lupton, E. M.; Achenbach, F.; Weis, J.; Bräuchle, C.; Frank, I. Modified chemistry of siloxanes under tensile stress: Interaction with environment. J. Phys. Chem. B 2006, 110, 14557–14563. 28

ACS Paragon Plus Environment

Page 28 of 32

Page 29 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(29) Lupton, E. M.; Achenbach, F.; Weis, J.; Bräuchle, C.; Frank, I. Molecular origins of adhesive failure: Siloxane elastomers pulled from a silica surface. Phys. Rev. B: Condens. Matter Mater. Phys. 2007, 76, 125420. (30) Nangia, S.; Garrison, B. J. Reaction rates and dissolution mechanisms of quartz as a function of pH. J. Phys. Chem. A 2008, 112, 2027–2033. (31) Lupton, E. M.; Achenbach, F.; Weis, J.; Bräuchle, C.; Frank, I. Origins of material failure in siloxane elastomers from first principles. ChemPhysChem 2009, 10, 119–123. (32) Rimsza, J. M.; Yeon, J.; van Duin, A. C. T.; Du, J. Water interactions with nanoporous silica: Comparison of ReaxFF and ab initio based molecular dynamics simulations. J. Phys. Chem. C 2016, 120, 24803–24816. (33) Schwaderer, P.; Funk, E.; Achenbach, F.; Weis, J.; Bräuchle, C.; Michaelis, J. Singlemolecule measurement of the strength of a siloxane bond. Langmuir 2008, 24, 1343– 1349. (34) Kohn, W.; Sham, L. J. Self-consistent equations including exchange and correlation effects. Phys. Rev. 1965, 140, A1133 – A1138. (35) Elstner, M.; Porezag, D.; Jungnickel, G.; Elsner, J.; Haugk, M.; Frauenheim, T.; Suhai, S.; Seifert, G. Self-consistent-charge density-functional tight-binding method for simulations of complex materials properties. Phys. Rev. B: Condens. Matter Mater. Phys. 1998, 58, 7260 – 7268. (36) Koskinen, P.; Mäkinen, V. Density-functional tight-binding for beginners. Comput. Mater. Sci. 2009, 47, 237 – 253. (37) Momma, K.; Izumi, F. VESTA3 for three-dimensional visualization of crystal, volumetric and morphology data. J. Appl. Crystallogr. 2011, 44, 1272–1276, VESTA is available at http://jp-minerals.org/vesta/en (accessed Aug 27, 2019). 29

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(38) 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 – 50. (39) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized gradient approximation made simple. Phys. Rev. Lett. 1996, 77, 3865 – 3868. (40) Blöchl, P. E. Projector augmented-wave method. Phys. Rev. B: Condens. Matter Mater. Phys. 1994, 50, 17953–17979. (41) Kresse, G.; Joubert, D. From ultrasoft pseudopotentials to the projector augmentedwave method. Phys. Rev. B: Condens. Matter Mater. Phys. 1999, 59, 1758–1775. (42) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu. J. Chem. Phys. 2010, 132, 154104. (43) Mermin, N. D. Thermal properties of the inhomogeneous electron gas. Phys. Rev. 1965, 137, A1441 – A1443. (44) Aradi, B.; Hourahine, B.; Frauenheim, T. DFTB+, a sparse matrix-based implementation of the DFTB method. J. Phys. Chem. A 2007, 111, 5678 – 5684, DFTB+ is available at https://www.dftb-plus.info (accessed Aug 27, 2019). (45) Rappe, A. K.; Casewit, C. J.; Colwell, K. S.; Goddard III, W. A.; Skiff, W. M. UFF, a full periodic table force field for molecular mechanics and molecular dynamics simulations. J. Am. Chem. Soc. 1992, 114, 10024–10035. (46) Slater, J. C.; Koster, G. F. Simplified LCAO method for the periodic potential problem. Phys. Rev. 1954, 94, 1498 – 1524. (47) Hagberg, A. A.; Schult, D. A.; Swart, P. J. Exploring network structure, dynamics, and function using NetworkX. Proceedings of the 7th Python in Science 30

ACS Paragon Plus Environment

Page 30 of 32

Page 31 of 32 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Conference. Pasadena, CA USA, 2008; pp 11 – 15, NetworkX is available at https://networkx.github.io (accessed Aug 27, 2019). (48) Stewart, J. J. P. Optimization of parameters for semiempirical methods I. Method. J. Comp. Chem. 1989, 10, 209 – 220. (49) Stukowski, A. Visualization and analysis of atomistic simulation data with OVITO the open visualization tool. Model. Simul. Mater. Sci. Eng. 2010, 18, 015012, OVITO is available at https://www.ovito.org (accessed Aug 27, 2019). (50) Plimpton, S. Fast parallel algorithms for short-range molecular dynamics. J. Comput. Phys. 1995, 117, 1 – 19, LAMMPS is available at http://lammps.sandia.gov (accessed Aug 27, 2019). (51) Frischknecht, A. L.;

Curro, J. G. Improved united atom force field for

poly(dimethylsiloxane). Macromol. 2003, 36, 2122–2129. (52) Thol, M.; Dubberke, F.; Rutkai, G.; Windmann, T.; Köster, A.; Span, R.; Vrabec, J. Fundamental equation of state correlation for hexamethyldisiloxane based on experimental and molecular simulation data. Fluid Ph. Equilibria 2016, 418, 133 – 151.

31

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Graphical TOC Entry Oriented H2O Attack

Tension Scission Under Tension

32

ACS Paragon Plus Environment

Page 32 of 32