Modeling Molecular Relaxation Mechanisms in ... - ACS Publications

Article pubs.acs.org/Macromolecules

Modeling Molecular Relaxation Mechanisms in Amorphous Polymers: Application to Polyamides Anthony Bocahut,† Jean-Yves Delannoy,‡ Didier R. Long,‡ and Paul Sotta*,‡ †

Advanced Polymers and Materials Department, Solvay, UMR5268, and ‡Laboratoire Polymères et Matériaux Avancés, CNRS/Solvay, UMR5268, Axel’One, 87 avenue des Frères Perret, 69192 Saint Fons, Cedex, France S Supporting Information *

ABSTRACT: We present a new multiscale hierarchical numerical approach combining ab initio, molecular dynamics, and metadynamics to model molecular motions in polymers over widely different time scales. We demonstrate that molecular motions corresponding to γ and β secondary (subglassy) relaxations can be described precisely in terms of free enthalpy, activation energy, and time scales and may be identified with relaxations detected experimentally by spectroscopic and other techniques. This method is illustrated in the case of aliphatic and semiaromatic polyamides. The highlight of this innovative approach in polymers is the ability to simulating long time scale motions while retaining an all-atomistic description of the polymer material. A complete description of the free energy landscape, in terms of conformational (intramolecular) and interactional (intermolecular) contributions, is obtained (energy barriers and time scales). This method thus may provide a predictive tool to understand secondary relaxations in amorphous polymers.

1. INTRODUCTION In amorphous polymers or in the amorphous phase of semicrystalline polymers, the main relaxation is the glass transition, often (not always) denoted as the α relaxation, which corresponds to large-scale motions involving several monomers collectively. At the glass transition temperature Tg, the amorphous phase of polymers changes from the glassy (rigid) state to a rubbery or fluid macroscopic state. While below Tg the polymer cannot flow, some local motions may still be activated, corresponding to subglassy, or secondary, relaxation processes generally denoted as β and γ (in order of increasing frequencies or decreasing temperatures). It is known that secondary relaxations, specifically the β relaxation, affect mechanical properties of polymers, particularly deformation, yield, and fracture behavior, toughness, impact strength, and brittle−tough transition.1 One way to improve thermomechanical properties in polymers is therefore to modulate these relaxations. For this, it is of crucial importance to identify and characterize in details the molecular mechanisms which are involved in the relaxations of interest. Experimentally, relaxations in polymers can be characterized by dynamic mechanical analysis (DMA)2 and broadband dielectric spectroscopy (BDS).3,4 These techniques allow obtaining a global overview of the various transitions through the evolution of the relaxation times as a function of temperature. The α (glass transition) relaxation follows a Williams−Landel−Ferry (WLF) empirical time vs temperature law.5 β and γ subglassy relaxations generally follow apparent Arrhenius laws, from which associated free enthalpy barriers may be determined. In addition, information on the spatial © XXXX American Chemical Society

scale of the considered motions may be obtained by quasielastic neutron scattering (QENS),6 while solid-state nuclear magnetic resonance (NMR)7,8 may give a more detailed description of the type and amplitude of motions. Nevertheless, identifying the molecular mechanisms involved in those relaxations is always a challenge. A way to tackle this problem is to use molecular modeling simulations. However, in many polymers, the different relaxations occur at very different time scales (in normal conditions (300 K, 1 bar), the β relaxation occurs in the range of μs−ms and γ in the range of ns−μs, and they cannot be modeled with a single method. Actually, relaxation times relevant in mechanical experiments are generally too slow to be conveniently treated by molecular dynamics. Generally speaking, both accessing long relaxation times (relevant for the mechanical behavior) and providing a detailed molecular picture of relaxation mechanisms are two antinomic requirements to meet at the same time. We have developed an innovative multiscale approach in order to describe motions of different amplitudes and time scales in polymers. Our approach is inspired by the one which has been developed to study relatively slow kinetic mechanisms and/or molecular flexibility in biological systems (diffusion of small molecules, large conformational changes, etc.).9−11 The hierarchical multiscale scheme which is proposed combines ab initio → molecular dynamics → metadynamics. Each step provides input parameters to the upper level. Ab initio Received: September 7, 2015 Revised: February 12, 2016

A

DOI: 10.1021/acs.macromol.5b01963 Macromolecules XXXX, XXX, XXX−XXX

Article

Macromolecules calculations (QM) are first used to optimize the force field parameters, mainly dihedral angles and charge parameters, that are then fed into molecular dynamics simulations. In strongly interacting polymers, these parameters require to be determined for each particular case. Such calculations are heavily time-consuming and have to be performed on small molecules (homologous of the monomer) in gas phase. At atomic scale, microstructural aspects in polymers involve both chain conformations, described by conformational degrees of freedom and intermolecular interactions. Energy barriers for the various degrees of freedom, which contain only conformational and intramolecular contributions, are obtained at this step. In a second step, classical atomistic molecular dynamics (MD) are performed to determine fast relaxation processes, including interactions such as hydrogen bonds. Electrons are not explicitly described, and punctual atomic charges remain constant along the simulation. Thus, in this approximation, all conformations are treated with the same set of charges. The accessible time scale is of order 10−100 ns. This level of simulation does not allow to observe large-scale conformational changes and long time-scale relaxation processes involved in secondary relaxations. Such motions may be reached either by coarse-grained or biased simulations.12 In our approach, for this third step, we have chosen to keep the atomic-scale molecular description of polymer chains and to bias MD simulations by using metadynamics simulations (MTD).13,14 MTD is a biased dynamical method in which a history-dependent bias potential is used to simulate rare events and reconstruct the free energy surface along a set of degrees of freedom. Relaxation times associated with free energy (or strictly speaking, free enthalpy, which is equivalent here) barriers ΔG, which may be well beyond the time scale directly accessible in MD simulations, are then inferred through an Arrhenius thermal activation equation τ = τ0 exp(ΔG/kBT). The initial conformation of MTD is obtained from MD. Therefore, the approach combines two advantages: accessible time scales may be very long, while a precise description of the particular type of motion is provided. Note that energy barriers determined in our multiscale approachcombine both conformational and interaction aspects. In this paper, we apply the proposed multiscale method to polyamide systems (PAs). To illustrate the relevance of the approach, we discuss the influence of an aromatic group inserted in PA6,6 chains on the local dynamics, namely on γ and β relaxations. Thus, two different PAs are considered and compared: PA6,6 (aliphatic) and PA6I (semiaromatic). Both aliphatic and amide parts of PAs will be studied by MD and MTD, at ns and μs−ms time scales, respectively. Whereas PA6I is amorphous, PA6,6 is semicrystalline. These systems have been studied experimentally.15 In dry PAs, in addition to the α (glass transition) relaxation, two secondary (subglassy) relaxations, denoted β and γ, occur. It has been suggested that γ involves the aliphatic part of the polymer and β the amide groups.16−19 The dynamics in the amorphous phase only is considered in this study. We argue that it may be of high potential interest in the particular case of polyamides due to their relatively low crystallinity index. Indeed, the crystallinity of polyamides is often in the range of roughly 30%, and it is well recognized that the molecular mobility in the amorphous phase has generally a drastic influence on mechanical and physical properties, such as solvent barrier and impact strength.3,20 The presence of water, even in relatively small amount, has a significant influence on both the amplitudes and time scales of the relaxations in PAs.

In some cases, and specifically in the presence of water, it has been shown that several distinct relaxation processes may coexist in the region of the relaxation map associated with β relaxation.15 Thus, accessing this time scale region by molecular modeling is of high potential interest to elucidate these various and complex relaxation processes, which in turn have a strong impact on the properties of the materials.

2. METHODS The method is based on a multiscale approach inspired by the approach used to describe protein dynamics in the field of molecular biophysics. Simulations are performed at increasing space and time scales, with increasing levels of coarse-graining approximation, from ab initio quantum molecular dynamics (QM) calculations, up to molecular dynamics (MD), up to metadynamics (MTD) calculations. Each level of calculations is used hierarchically to feed up the next level with finely tuned values of simulation parameters. 2.1. Ab Initio Calculation Protocol. The purpose of ab initio calculations (QM) is to optimize the force field parameters, mainly dihedral angles and charge parameters, that are then used in MD simulations. In flexible hydrogen bonded polymers like PAs, these parameters are very important and require to be previously determined. QM calculations are time-consuming and therefore can only be performed on noninteracting small molecules, i.e., in gas phase. The smallest molecules (PA analogues) of PA6,6 and PA6I, Nhexylpentanamide and N-phenylpentanamide, respectively, are schematized in Figure 1. These molecules are built using the molecule

Figure 1. Schematic representation of the different degrees of freedom involved in polyamides: PA6,6 (A) and PA6I (B). Aliphatic degrees of freedom are described by α and β dihedral angles in red. Amide degrees of freedom are described by ϕ (blue), ψ (green), and ω (orange) dihedral angles. Because of the presence of aromatic ring, β angles are not present in PA6I. builder of MAPS.21 The geometry of the molecules is optimized with an optimization run of 100 steps using steepest descent algorithm and the Dreiding force field22 included in MAPS.21 The force convergence criterion is 0.01 kJ/mol. QM is done with Turbomole.23 Calculations are carried out at the Hartree−Fock(HF)/6-31G* level of theory, which offers a satisfactory compromise between precision and computation time and is generally used for strongly polarizable molecules.24,25 Possible organizations of polymer chains in bulk amorphous state can be determined at the scale of all-atom MD. QM calculations are used prior to MD in order to determine authorized local conformations, based on the quantitative determination of energy minima and energy barriers. For this, given pairs of degrees of freedom (dihedral angles) are chosen, while maintaining the other ones fixed at their most probable values. Each minima of energy may be determined and optimized and full energy maps (so-called Ramachandran maps) for all values of the pair of chosen angles can be drawn. Around local B


Article

Macromolecules energy minima, the potential wells are constructed by rotating the chosen angles in a window of typically ±30° to ±40° with an increment of typically 10°. Mulliken charges26 are calculated on the optimized PA analogue molecules by single point energy calculations of the electronic wave functions with HF method in split-valence polarization (SVP) bases using Turbomole software.23 The calculated Mulliken charges are reported in the Supporting Information. Charge determination by QM calculations was necessary for PA6I (using N-phenylpentanamide), in which no data have been reported to our knowledge. Then, to ensure consistent and robust comparison between PA6I and PA6,6, we have kept exactly the same approach for PA66, even though optimized force fields have been already proposed in this case.27 As a matter of fact, by contrast to what was proposed in ref 27, we do not impose a priori equal charges on all aliphatic carbons. It turns out indeed that the charge of a given aliphatic carbon depends significantly on its position with respect to the amide group (see Supporting Information). While aliphatic carbons in α position to amide groups have a charge comparable to those in proteins, aliphatic carbons located further away are less charged. 2.2. Molecular Dynamics Protocol. 2.2.1. Amorphous Box Generation. One chain of polyamide (PA6,6 or PA6I) containing 80 monomers (MN ≃ 20 000 g/mol) is constructed using the polymer builder included in MAPS. The geometry of the polymer chain is optimized using steepest descent algorithm and the Dreiding force field with the optimizer tools included in MAPS. The bulk polymer is constructed using an amorphous cell program included in MAPS, in which the amorphous phase of a glassy polymer is created using excluded volume interactions and nonbonded interactions to avoid overlaps between chains.28 The initial system configuration is composed of 20 polymer chains in coiled conformations built inside a box with periodic boundary conditions, at low density (0.5 g/cm3). The total number of atoms is about 60 000. An optimization run of maximum 10 steps with 0.01 kJ/mol of force convergence criterion using a steepest descent algorithm and the Dreiding force field is performed. 2.2.2. Details of Simulations. Molecular Dynamics (MD) simulations have been performed using the GROMACS software (GROningen MAchine for Chemical Simulations)29−32 on the Solvay computer cluster. We use the all-atoms OPLS (optimized potentials for liquid simulations) force field,33 in which all the amide, aliphatic, and phenyl groups can be described. Atomic charges determined previously by QM calculations are used as input in OPLS. Long-range electrostatic interactions are treated using the particle mesh Ewald (PME) method,34 with a grid spacing of 0.12 nm and a nonbond pair list cutoff of 1.0 nm. The pair list is updated every five steps. The time step is 1 fs, and the system is kept stable by constraining bond lengths involving H atoms with the LINCS algorithm (Linear Constraint Solver).35 50 ps of molecular dynamics is performed in the NVT ensemble at 600 K, followed by 100 ps in the NPT ensemble at 1 bar and 600 K to equilibrate the density. The temperature of the system is then dropped to 300 K using a simulated annealing process, with a cooling rate of 1 K/ps.36 The simulated annealing combines a decrease of temperature and Monte Carlo steps to equilibrate the density at each temperature. It ensures that the whole energy landscape is explored, with properly thermally equilibrated populations. It is necessary to obtain a system with an equilibrated density at 300 K. Finally, a molecular dynamics run of 10 ns in the NPT ensemble at 300 K and 1 bar is performed Only the last 2 ns is considered for analysis (production step). It was checked that the simulated PA6,6 does not crystallize. 2.2.3. Force Field Validation. In order to obtain representative simulated polymeric systems, degrees of polymerization (DP) which reproduce experimental data have to be determined. Furthermore, in our polar systems, dipolar interactions (hydrogen bonds) between chains modify the chain configurations. The ratio of intra- and intermolecular interactions will therefore be affected by the number of monomers in the chains. Thus, different DPs (4, 20, 40, 80, 100, and 200) have been simulated (note that one monomer contains two amide groups), and some local properties have been determined.

Obtaining the proper density in quenched systems is generally considered to be the primary criterium of reaching proper equilibration in molecular modeling of polymers. The density is found to be 1.06 g/cm3 for DP ≥ 20. An experimental density value of 1.07 g/cm3 for the amorphous phase of PA6,6 has been obtained by extrapolating the experimentally measured density as a function of the crystalline ratio to 0% crystallinity.20 Thus, the value found in the simulations corresponds to a relative error between experimental and simulated data not larger than 0.8%. In the case of PA6,I we find a density value of 1.14 g/cm3 at DP = 80, while the experimental density of amorphous PA6I is around 1.15 g/cm3,37 which corresponds to a relative error 0.9%. Thus, in both cases, simulated data are in excellent agreement with experimental ones. In the amorphous cell approach, it is important to check that satisfactory chain conformations are obtained. At a semilocal scale, chain conformations may be described for instance by the characteristic ratio, CM, which relates the chain end-to-end average distance to the monomer number M. CM has been calculated for PA chains having 4−100 monomers by considering that PA chains behave as Gaussian chains.38,39 It was found that CM decreases as M increases, in agreement with the Flory−Fox empirical law.40 Experimentally, Saunders found C∞ = 5.941 and Flory et al. C∞ = 6.1,42 and a similar, though a little larger, value C∞ = 6.18 is reported by Fetters et al. for PA6 at 543 K.43 Our data for PA6,6 are compatible with an extrapolated value C∞ ≃ 6 and are thus in agreement with experimental values (see Figure 4a in Supporting Information). Also, the computed squared radii of gyration for PA6 give a satisfactory linear variation with M (Figure 4b in Supporting Information). We therefore infer that the OPLS force field used in this work well reproduces local parameters of aliphatic and semiaromatic PAs, as reflected by the density and characteristic ratio. 2.3. Metadynamics Protocol. While MD simulations cannot describe relaxation processes with time scales longer than a few nanoseconds, either coarse-grained or biased simulations allow accessing such motions.12 We have chosen to retain an atomic-scale approach and to bias MD simulations by using metadynamics simulations (MTD).13,14 MTD is a biased dynamical method able to simulate rare events and reconstruct free energy surfaces along a set {Ψi} of degrees of freedom (or so-called collective variables, CVs), using a history-dependent potential. This quite recent method is described in detail in the literature.13,14 Along a MTD simulation, an external potential VG, which disfavors the conformations {Ψi} that have already been visited, is added recursively. This potential is constructed as a sum of Gaussian functions centered around the values of the {Ψi(t)} explored during the dynamics: VG({Ψi(t )}) = w

∑ t ′= τG ,2τG ,...; t ′< t

⎡ (Ψ − Ψ(t ′))2 ⎤ i i ⎥ 2δ 2 ⎣ ⎦

∑ exp⎢− i

(1)

where Ψi(t′) is the value of the CV Ψi at time t′. Three parameters are used to describe the potential VG: the Gaussian height w, the Gaussian width δ, and the deposit time interval of the Gaussian τG. At the end, the effective free energy surface constructed through the MTD is flat. To recover the “real” energy landscape, the potential VG must be eliminated and the free energy of the system is given by

lim VG({Ψi(t )}) ≃ − F({Ψi(t )})

t →∞

(2)

The initial conformation of MTD is obtained by running 10 ns MD simulations of a single chain of PAs polymer (DP = 80). The polymer chain is equilibrated by using the protocol described in section 2.2. The Plumed package,44 a portable plug-in included in Gromacs software, has been used. A monomer in the middle of the chain is chosen for conformational analysis and the CVs are identified on this monomer for each type of PA. The width and height of the Gaussian potential functions used in MTD calculations are conventional one for dihedrals angle variables,44 with parameters optimized to obtain a good compromise between computing time and accuracy on the obtained free energy map. A Gaussian width δ = 0.35 rad = 20° has been chosen. It gives the C


Article

Macromolecules angular resolution on the definition of potential wells in the map. The Gaussian height w = 0.25 kJ/mol = 0.1 kBT is constant throughout the MTD simulation. It gives the resolution (uncertainty) on the measured barrier height. Finally, a deposit time τG = 100 fs was chosen as a compromise between unbiased exploration of the map and computer time.

3. RESULTS 3.1. Results of ab Initio Calculations: Case of PA6,6. The conformational degrees of freedom in PA6,6 are represented by eight different dihedral (torsion) angles, as depicted in Figure 1A. αi (i = 1, 2, 3) and βj (j = 1, 2) correspond to the aliphatic segment of PA between amine groups (NH) and between carbonyl groups (CO), respectively. ω, ϕ, and ψ correspond to amide groups. ω describes the torsion angle O−C−N−H, ϕ the angle CH2−CH2−NH−CO, and ψ the angle CH2−CH2−C−O. The conformations of aliphatic carbons are well established in the literature. 45 Aliphatic chains can describe four conformations, denoted c (cis, α = 0°), g+ (gauche+, α = 60°), t (trans, α = 180°), and g− (gauche−, α = 300°). They mainly adopt the t (trans) conformations where X−CH2− CH2−X dihedral angles are equal to 180°, except β1 (CO− CH2−CH2−CH2) which adopts the three conformations t, g+, and g− in roughly the same proportions (∼30%).46 In order to determine the relative conformations of the three angles ω, ϕ, and ψ, QM calculations have been performed on N-hexylpentanamide (Figure 1A). 3.1.1. CONH Conformations. Amide bonds present a planar geometry due to orbital overlap between the π-orbitals of the sp2-hybridized atoms of the carbonyl group (CO) and the πorbitals of the sp2-hybridized amine group (NH). This delocalized π-orbitals system does not allow free rotation of the ω dihedral angle (see Figure 1). Ab initio calculations have shown that the trans (t) conformation (ω = 180°, with CO and NH pointing in opposite directions, as depicted in Figure 1) is more stable than the cis (c) conformation (ω = 0°, CO and NH pointing in the same direction) with an energy difference ΔE = 12 kJ/mol ≃ 5 kBT at room temperature (T = 300 K) and a high-energy barrier 125 kJ/mol ≃ 50 kBT. This result has also been established for short aliphatic segments (in N-methylacetamide47), suggesting that ω is not influenced by the size of the molecule, as well as in proteins, which can be seen as copolymers of PA, using various force fields such as OPLS, CHARMM, and AMBER.48−50 Thus, the ω angle must be almost exclusively in the trans position in PA6,6, with less than 1% of the bonds in cis position. This small fraction is not detected in practice. 3.1.2. ϕ−ψ Conformations. The authorized conformations of ϕ (dihedral angle on the NH side of amide group) and ψ (carbonyl side) in N-hexylpentanamide (Figure 1) have first been determined by QM calculations, αi (i = 1, 2, 3) and βj (j = 1, 2) being fixed at 180°. Each energy minimum has been determined and optimized. The full energy map (Ramachandran map) for all values of the ϕ and ψ angles is shown in Figure 2A. For both ϕ and ψ, conformations around the cis position (ϕ = 0 and ψ = 0) are excluded. A global basin of lower energy is observed in the range [60°:300°] for ϕ and in the range [90°:270°] for ψ. In this basin, there are three stable conformations for either ϕ (90°, 180°, and 270° or g+, t, and g−) or ψ (120°, 180°, and 240° or g+, t, and g−). Given that the map is centrosymetric with respect to (180°/180° or tt) (see

Figure 2. Enthalpy surface in kJ/mol as a function of the two dihedral angles ϕ and ψ at the HF/6-31G* level of theory of Nhexylpentanamide (A) and N-hexylphenylamide (B). The maps are centrosymmetric around (180°, 180°) modulo 360°. Additionally, in the map in (B), angles ψ = 0° and 180° are equivalent.

Figure 2A), this corresponds to five independent stable conformations of the (ϕ/ψ) couple: (90°/120° or g+g+), (90°/180° or g+t), (90°/240° or g+g−), (180°/120° or tg+), and (180°/180° or tt). The Newman representations of ϕ and ψ stable conformations are given in the Supporting Information. The activation energy between each ϕ local basin is ΔHϕ/ψ conf = 2 ± 1 kJ/mol ≃ 1 kBT. Thus, based on this conformational analysis, ϕ has little influence on ψ and may undergo free rotation at room temperature. 3.1.3. β1−ψ Conformations. As mentioned previously, β1 is the only aliphatic angle which could provide three isoenergetic conformations (g+, t, and g−). A possible way to investigate the influence of β1 on (ϕ, ψ) basins is to look at the (β1, ψ) 2D energy map while fixing ϕ. Indeed, the angle ϕ can rotate in the D


Article

Macromolecules [60°, 300°] interval with little restriction as compared to ψ (see Figure 2B) and almost not effect on ψ, and it may thus be fixed to any value within this interval. The potential energy map has been constructed by rotating β1 and ψ of ±180° from the initial conformation β1 = ψ = 180°, while keeping ϕ = 180° (t conformation). This (β1,ψ) energy map is shown in Figure 3.

minima corresponding to the three positions (90°, 180°, 270°) are much less pronounced. Altogether, the observed minima correspond to three distinct (ϕ, ψ) conformations. Newman representations of ϕ and ψ stable conformations are given in the Supporting Information. In the PA6I monomer, the presence of a second amide group in the meta position of the aromatic ring is responsible for a symmetry breaking, which induces that the two basins at ψ = 0° and ψ = 180° become different. This difference will be developed using MD calculation (section 3.5). 3.2.2. ψ1−ψ2 Conformations. The ψ dihedral angles (with respect to the plane of the phenyl ring) on each side of the phenyl ring in the monomer subunit may influence each other. In order to determine this influence, the symmetric phenyldibutanamide molecule has been considered for QM calculations. The energy map for the (ψ1, ψ2) angles is shown in Figure 4. The four energy minima correspond to the four possible conformations A to D. The map has axial symmetry (ψ1, ψ2) = (ψ2, ψ1), and A and D are completely symmetric and represent the same conformation. From a structural point of view, the conformation C can be seen as a linear segment, whereas conformations A (D) and B describe coil structures, with a coil angle of roughly 120° for A (D) and 90° for B. The more stable conformation is A (D); then conformation B has an excess energy of 8 kJ/mol ≃ 3 kBT as compared to A 1/ψ2 (D), with energy barriers ΔHψB→A = 14 ± 2 kJ/mol ≃ 5.5 kBT from B to A and 22 kJ/mol ≃ 9 kBT from A to B; then conformation C is the less stable conformation with a excess energy of 12 kJ/mol ≃ 5 kBT as compared to A and energy 1/ψ2 = 27 ± 2 kJ/mol ≃ 11 kBT from A to C and 15 barriers ΔHψA→C kJ/mol ≃ 6 kBT from C to A. It is interesting to note that conformation C is the most often depicted one in the literature,15,37 even though it is the least stable energetically. Note also that such alternating conformations of carbonyl groups on each sides of the aromatic rings (with ψ1 = ψ2 + 180°) have been observed in aromatic PA used in membranes.51 3.3. Charge Calculations. In order to optimize the charge distributions in PA6,6 and PA6I, representative conformers of the monomer have to be selected. In PA6,6, due to the large energy basin of (ϕ/ψ) conformations (see Figure 2), the linear conformer of the hexyldipentanamide (ω = ϕ = ψ = 180°) has been chosen for punctual charge calculations. In PA6I, due to the relatively high energy barriers between the different states (Figure 4), only the most stable conformations A (D) have been kept for charge calculations. Since these conformations are nonsymmetrical (ψ1 not equal to ψ2), charges computed in these conformations are not symmetrical, i.e., are different on amide groups 1 and 2 on either side of the phenyl ring. Thus, in order to obtain a charge distribution symmetrical with respect to the phenyl ring, we have taken the average of the charges in both A and D conformations. The complete tables of punctual Mulliken charge values for PA6,6 and PA6I, which are then used for MD simulations, are given in the Supporting Information. 3.4. Results of Molecular Dynamics Simulations in PA6,6. 3.4.1. ϕ−ψ Conformations. Only the t conformation of the ω (O−CO−NH−HN) angle (ω = 180°) is observed, in full agreement with QM calculations where the t isomer was found to be the most stable conformation. Thus, it is possible to consider that modified-OPLS force fields for PAs reproduce very well the ω angle conformation in the PA6,6 amorphous

Figure 3. Enthalpy surface in kJ/mol as a function of the two dihedral angles β1 and ψ at the HF/6-31G* level of theory of Nhexylpentanamide. The map is centrosymmetric around (180°, 180°) modulo 360°. Note that the energy level rises quite steeply for β1 below 60° and beyond 300°.

Outside the range [60°:300°] for β1 and ψ, a high-energy barrier is observed, which means that cis (c) conformations (β1 = 0 and ψ = 0) are excluded. Three potential wells appear for β1, corresponding to the g+ (β1 = 60°), t (β1 = 180°), and g− (β1 = 300°) conformations. For each β1 value, energy minima correspond to different values of ψ, which altogether gives three different isomers separated by energy barriers of about 12 ± 2 kJ/mol ≃ 5 kBT. Note that the map is centrosymmetric around (180°, 180°). 3.2. Results of ab Initio Calculations: Case of PA6I. The conformational degrees of freedom of PA6I are depicted in Figure 1B. The aromatic ring is completely rigid. Considering the symmetry of the molecule in 2D plane, six different angles appear: αi (i = 1, 2, 3), ω, ϕ, and ψ. As already mentionned, aliphatic torsion angles αi mainly adopt the t conformation (α = 180°). In order to determine the relative conformations of the ω, ϕ, and ψ angles, which may be altered by the presence of the aromatic ring, QM calculations have been performed on Nhexylphenylamide (see Figure 1B). 3.2.1. ϕ−ψ Conformations. As regards ω, c conformations are completely negligible as compared to t conformations in PA6I chains, as in PA,6. The complete energy map for the (ϕ, ψ) couple of angles in N-hexylphenylamide is shown in Figure 2B, in comparison to the map obtained in PA6,6 (Figure 2A). For ψ, relatively sharp minima correspond to c (ψ = 0°) and t (ψ = 180°) conformations, separated by an energy barrier of 25 kJ/mol ≃ 10 kBT. Because of phenyl ring symmetry, these two minima describe the same planar conformation. For ϕ, the E


Article

Macromolecules

Figure 4. Enthalpy surface in kJ/mol as a function of the two dihedral angles ψ1 and ψ2 at the HF/6-31G* level of theory of phenyldibutanamide. The different PA6I monomer conformations around the phenyl ring (illustrated on phenyldibutanamide) corresponding to the four minima are shown.

bulk. The ω torsion angle is completely rigid and cannot be involved in polymer relaxation. Conversely, both ϕ and ψ are effective degrees of freedom. The (ϕ/ψ) distribution of populations for PA6,6 is shown in Figure 5A. All observed conformations are located within a single basin in the range [60°:300°] for both ϕ and ψ. Within this region, two local maxima of population are observed, which actually correspond to the same ensemble of conformations, the map being centrosymmetric around (180°,180°). The overall fraction of conformations in these two symmetrical basins of population is ≃75%. Conversely, the linear conformation (ϕ = 180°/ψ = 180°), which is the one present in crystalline PA6,6, is observed in quite low fraction (≃25%). The ϕ angle lies within a relatively sharp interval of values ([60°:120°]) centered on 90° (due to symmetry, only the maximum on the left-hand side of the map is discussed), while the ψ angle is distributed over a broader interval [90°:240°],

with two preferred conformations around ψ = 120° and around ψ = 180°. Altogether, the map shown in Figure 5A is in excellent agreement with the results of QM calculation (see in Figure 2A). The two symmetric basins of population indeed correspond to stable conformations in QM calculations, even though QM energy maps were calculated for a fixed value β1 = 180°. Thus, this shows that the presence of other isomers of β1 does not affect significantly the distribution of (ϕ/ψ) conformations. It was shown that the flip between ϕ = 90° and ϕ = 270° involves an energy barrier of 2 kJ/mol ≃ 1 kBT (see Figure 2), which allows almost free rotation at room temperature, while the only possible large amplitude motion of ψ is the 60° flip. 3.4.2. Aliphatic Part. First, αi (i = 1, 2, 3) dihedral angles corresponding to aliphatic segments between NH groups have been investigated. MD results give predominant populations (≃75%) around 180°, indicating that all αi are mainly in t F


Article

Macromolecules

torsion angles (β1, ψ) obtained in MD is shown in Figure 6. Three local maxima of populations are observed, corresponding

Figure 6. Distribution of β1/ψ torsion angles along the MD production period of PA6,6.

to the three β1 conformations: 60° (g+), 180° (t), and 300° (g−) with relative populations of 37%, 26%, and 37%, respectively. Note that g+ and g− effectively correspond to the same ensemble of conformations, the map being centrosymmetric around (180°, 180°). The basins of populations are relatively narrow in β1 and correspond to broader distributions of ψ around 180°. Table 1 summarizes the different populations for the three angles (ϕ/ψ/β1). 3.5. Results of Molecular Dynamics Simulation in PA6I. 3.5.1. ω Conformations. As in PA6,6, the distribution of the ω angle (O−CO−NH−HN) dihedral has been determined by MD simulations. About 15% of c (cis) conformations are observed. However, as argued in the QM part, c conformations in PA6I are completely excluded. They should be eliminated Table 1. Conformations of the Angles ϕ, ψ, and β1; Populations of (ϕ/ψ/β1) Conformations

Figure 5. Distributions of populations for ϕ/ψ torsion angles along the MD production step of PA6,6 (A) and PA6I (B). The maps are centrosymmetric around (180°, 180°) modulo 360°. Note the similarity of (B) with Figure 2B.

(trans) conformations. The same is observed for β2 (central aliphatic segments between CO groups). Regarding β1, the three conformations (g+, t, and g−) appear roughly with the same ratio of population (30%), in full agreement with QM results.46 In section 3.1, an energy barrier of 10 kJ/mol (4 kBT) has been determined between g+ or g− and t. The relaxation time of the angle β1 corresponding to this high energy barrier would not be accessible in the time scale of MD. However, β1 should be considered as an important degree of freedom in PA6,6 relaxation processes, not only involving aliphatic parts but also providing a coupling between aliphatic and amide moieties. 3.4.3. β1−ψ Conformations. A possible way to investigate the role of β1 is to look at (β1, ψ) energy maps, since both angles seem to be related together. The distribution of the G

ϕ ψ β1 ϕψβ1

90° (g+) 120° (g+) 60° (g+)

180° (t) 180° (t) 180° (t) populations (%)

g+g+g+ g+g+t g+g+g− g+tg+ g+tt g+tg− g+g−g− g+g−t g+g−g+ tg+g+ tg+t tg+g− ttg+ ttt

(g−g−g−) (g−g−t) (g−g−g+) (g−tg−) (g−tt) (g−tg+) (g−g+g+) (g−g+t) (g−g+g−) (tg−g−) (tg−t) (tg−g+) (ttg−)

4.0 ± 0.2 5.8 ± 0.1 16.5 ± 0.3 10.6 ± 0.3 10.4 ± 0.2 12.6 ± 0.8 8.3 ± 0.2 3.7 ± 0.1 2.9 ± 0.1 2.6 ± 0.1 3.2 ± 0.1 8.0 ± 0.1 8.0 ± 0.2 3.5

90° (g−) 240° (g−) 300° (g−) total 75%

total 25%


Article

Macromolecules before starting MD simulations, since c conformations cannot be converted into t conformations during MD due to a very high energy barrier (110 kJ/mol). This is done in the following way. Hydrogen linked to nitrogen atoms (NH) of c conformations are exchanged with aliphatic carbon atoms (CH2) linked to the same nitrogen atom. The modified structure is optimized without LINCS constraints on hydrogen atoms. During the minimization step, most modified c conformations move off to t conformations, while some of them go back to their initial c conformations, due to a too restricted environment. Thus, it is necessary to carry out several steps. After 10 correction/optimization iterative steps have been performed, less than 1% of ω angles remain in c conformations, which provides a good starting point for MD simulations. These c and t populations are kept constant all along the MD procedure. 3.5.2. ϕ−ψ Conformations. In PA6I, as in PA6,6, MD simulations confirm that aliphatic αi dihedral angles are mainly in t conformations (population ≃ 70%). The other accessible degrees of freedom are the (ϕ/ψ) torsion angles, whose distribution is shown in Figure 5B. As regards ψ, two relatively narrow basins of population are observed at ψ = 0° and ψ = 180°, which actually correspond to the same conformation with amide groups in the same plane as phenyl rings, due to the symmetry of the monomer (see Figure 1). Conversely, in each basin, ϕ is relatively broadly distributed between ϕ = 60° and ϕ = 300°, with three maxima of population at ϕ = 90°, ϕ = 180°, and ϕ = 270° (the two most populated maxima at 90° and 270° corresponding to mirror-symmetric conformations). Thus, altogether, six (ϕ/ψ) conformations are possible around amide bonds in PA6I: (ϕ = ±90°/ψ = 0°), (ϕ = ±90°/ ψ = 180°), (ϕ = 0°/ψ = 0°), and (ϕ = 0°/ψ = 180°) (see the Newman representations in Supporting Information). The map in Figure 5B is in excellent agreement with QM calculations (see Figure 2B), where two potential wells were observed as well at ψ = 0° and ψ = 180°. There is however a significant difference that is emphasized in MD simulation. Considering only one side of the phenyl ring, conformations with ψ = 0° and ψ = 180° should be strictly equivalent, as it is indeed the case in Figure 2B. However, in PA6I, considering the relative positions on both sides of the phenyl ring, since (ψ1 = 180°,ψ2 = 180°) (conformation B in Figure 4) is favored compared to (ψ1 = 0°, ψ2 = 0°) (conformation C in Figure 4), the (ψ = 0°, ψ = 180°) symmetry is broken and the ψ = 180° energy minima are more populated than ψ = 0°, with relative populations of ≃60% and ≃40%, respectively. On the other hand, conformations with ϕ = ±90° are mirror-symmetric and remain energetically equivalent. 3.5.3. ψ1−ψ2 Conformations. To elucidate the influence of both ψ angles around the aromatic ring, we now look at the (ψ1, ψ2) energy map, in a similar way as in section 3.2. The distribution of (ψ1, ψ2) angles in PA6I is shown in Figure 7A. Four conformations are observed, corresponding to the isomers A, B, C, and D (A and D corresponding to the same conformation by symmetry) represented in Figure 4, with relative populations ∼50% for both A and D together, ∼33% for B, and ∼17% for C. The distribution in Figure 7A is in good overall agreement with the energy map obtained by QM calculations (Figure 4). As regards the C conformation, it was indeed found in QM calculations to be the less stable one (higher energy state) among A (D), B, and C, which agrees with Figure 7A, in which C has indeed the lowest relative population. As regards the B conformation, it appears to be the most populated in MD,

Figure 7. Distribution of (A) (ψ1/ψ2) torsion angles and (B) (ϕ1/ϕ2) along the MD production period of PA6I.

which does not fit to QM results in which A and D conformations were the lowest energy states, suggesting that A and D should be the most populated isomers in PA6I. However, QM calculations were performed at the fixed value ϕ = 180°, whereas in MD simulations it is observed that the main populations are located at ϕ = ±90° (Figure 5B). Such a result indicates that the (ψ1/ψ2) relative populations can be altered by the conformations of both ϕ1 and ϕ2 angles. 3.5.4. ϕ1−ϕ2 Conformations. Thus, in order to get a global overview of the possible monomer conformations in PA6I, the (ϕ1,ϕ2) population distributions are now investigated. The (ϕ1/ϕ2) distribution is shown in Figure 7B. Nine distinct basins of population appear in the range [60°,300°] for both ϕ1 and ϕ2. Given the symmetries of the molecule, these basins correspond to four distinct sets of conformations: set A (90°/90° and 270°/270°, or g+g+ and g−g−) in which the Cα−Cβ bonds on each side of the phenylamide entity are perpendicular to the plane of the H


Article

Macromolecules phenylamide entity in the same direction; B (90°/270° and 270°/90°, or g−g+ and g+g−) in which the Cα−Cβ bonds are perpendicular to the plane of the phenylamide entity in opposite directions; D (180°/180°, or tt) in which both Cα−Cβ bonds are in the plane of the phenylamide entity in t conformations; and C (90°/180°, 180°/90°, 180°/270°, and 270°/180°, or g+ t,tg+, tg−, and g−t) in which the Cα−Cβ bond on one side of the phenylamide entity is in the plane in t conformation while the one on the other side is perpendicular to the plane of the phenylamide entity. Note that the last four conformations are not strictly equivalent, and two subsets should be defined corresponding to tg+ and g−t on one side and tg− and g+t on the other side. Table 2 summarizes the different populations for the four angles (ϕ1/ψ1/ϕ2/ψ2). Table 2. Populations of (ϕ1/ψ1/ϕ2/ψ2) Conformations ϕ1/ϕ2 ψ1/ψ2 ϕ1/ψ1/ϕ2/ψ2 g+tg+t g+ttt tttt g−ttt g−tg+t g+cg+t g+ctt tctt g−ctt g−cg+t g+cg+c g+ctc tctc g−ctc g−cg+c

90° (g+) 0° (c) (g−tg−t) (ttg−t) (ttg+t) (g+tg−t) (g−cg−t; g+tg+c; g−tg−c) (tcg−t; g+ttc; ttg−c) (tttc) (tcg+t; g−ttc; ttg+c) (g+cg−t; g−tg+c; g+tg−c) (g−cg−c) (tcg−c) (tcg+c) (g+cg−c)

90° (g−) 180° (t) populations (%)

180° (t)

5.8 7.4 4.5 7.7 7.0 12.5 11.2 5.2 9.9 11.0 4.5 3.4 1.3 3.1 3.4

± 0.7 ± 0.3 ± ± ± ± ± ± ± ± ±

0.1 0.9 0.5 0.4 0.4 0.2 0.3 0.4 0.3

total 33%

total 50%

total 16%

± 0.3 ± 0.1

3.6. Metadynamics: Results in PA6,6. The metadynamic (MTD) approach is based on the indirect estimation of relatively long relaxation times through the efficient scanning of energy maps showing large energy barriers, assuming thermally activated processes. Motions of large amplitude leading to conformational changes at the scale of the PA6,6 monomer involve the three degrees of freedom ϕ, ψ, and β1 described previously. However, MTD simulations using three collective variables (CVs) are difficult to perform and analyze. Thus, to restrict the dynamics to two CVs, we shall fix one degree of freedom and let the other ones vary. We consider here the monomer in the middle of the polymer chain. 3.6.1. ψ−β1 Free Energy Barrier. The previous simulation steps have shown that the position of the ϕ dihedral angle does not affect much the ψ equilibrium positions (see Figures 2A and 5A) and that the preferred value for the ϕ dihedral angle is 90° (mod 180°). Thus, the ϕ dihedral angle may first be fixed at 90°, while ψ and β1 are chosen as CVs. The obtained free energy map as a fonction of the (ψ, β1) angles is shown in Figure 8A. This map shows energy barriers of up to 40 kJ/mol ≃ 16 kBT, which would have been impossible to explore within the time scale of conventional MD (10 ns corresponding to ∼20 kJ/mol). From the initial (ψ = 120°/β1 = 300°) conformation (conformation A in Figure 8A; see also Figure 10), three local free energy minima can be reached: (ψ = 120°/ β1 = 60°), (ψ = 270°/β1 = 300°), and (ψ = 240°/β1 = 60°),

Figure 8. Energy surfaces in kJ/mol of PA6,6 obtained by metadynamics (MTD) with two collectives variables (CVs) ψ/β1 with ϕ = 90° (A) and ϕ/β1 with ψ = 240° (B). Arrows indicate large amplitude motions discussed in section 5.

with relative d free energies as compared to initial state of respectively 17 kJ/mol ≃ 7 kBT, 8 kJ/mol ≃ 3 kBT, and 10 kJ/ mol ≃ 4 kBT. Note also the similarity of the map in Figure 8A with Figure 3. The locally stable conformation most remote from the initial state in terms of angular variation is located at (ψ = 240°/β1 = 60°) (conformation B in Figures 8A and 10), with a free energy barrier of 30 ± 5 kJ/mol ≃ 12 kBT. 3.6.2. ϕ−β1 Free Energy Barrier. However, due to the constraints applied on the ϕ angle along the first MTD simulation step, the system may not have reached its most stable conformation, in state B, even though it has been shown before that the impact of ϕ on β1 and ψ equilibrium values is rather weak. To further verify this point, in a second step, the constraint on ϕ has been released and the influence of ϕ determined. Thus, the new conformation generated from MTD (ϕ = 90°, ψ = 240°/β1 = 60°, conformation B in Figure 8A) has I


Article

Macromolecules then been taken as the starting point for a subsequent simulation step. ψ is fixed to 240°, and ϕ and β1 dihedral angles are chosen as CVs. The obtained free energy map as a fonction of the (ϕ, β1) angles is shown in Figure 8B. From the initial (ϕ = 90°, β1 = 60°) conformation (conformation B in Figure 8B), five local energy minima are observed. Interestingly, the most stable conformation is not the initial state (B) but the conformation (ϕ = 270°, β1 = 60°) denoted as C in Figures 8B and 10, which suggests that the system could move naturally toward this conformation if no constraint was imposed on the ϕ dihedral angle. It could be also interesting to perform a second MTD with conformation C as initial state in order to compare the energy barriers. The relative free energies of each conformation as compared to C can be directly read on the map in Figure 8B: 10 kJ/mol ≃ 4 kBT for the initial conformation B; 32 kJ/mol ≃ 13 kBT for (ϕ = 90°, β1 = 180°); 33 kJ/mol ≃ 13 kBT for (ϕ = 90°, β1 = 300°); 17 kJ/mol ≃ 7 kBT for (ϕ = 270°, β1 = 180°); and 23 kJ/mol ≃ 9 kBT for (ϕ = 270°, β1 = 300°). The conformational change from B to C presents a free energy barrier of 15 ± 5 kJ/mol ≃ 6 kBT. 3.6.3. Free Energy Balance. The total activation free energy ΔGMTD corresponding to the conformational change between tot conformations A (ϕ = 90°, ψ = 120°, β1 = 300°, Figure 8A) and C (ϕ = 270°, ψ = 240°, β1 = 60°, Figure 8B) may be evaluated by combining both successive MTD simulations. At first, the exchange between conformations A and B involves a free = 30 ± 5 kJ/mol ≃ 12 kBT, with an energy barrier ΔGMTD 1 /kBT). Then, the associated relaxation time τ1 = τ0 exp(ΔGMTD 1 conformational change between B (ϕ = 90°/ψ = 240°/β1 = 60°) and C (ϕ = 270°/ψ = 240°/β1 = 60°) involves a free = 15 ± 5 kJ/mol ≃ 6 kBT, with a energy barrier of ΔGMTD 2 /kBT). relaxation time τ2 = τ0 exp(ΔGMTD 2 Considering that the global reorientation from A (ϕ = 90°/ψ = 120°/β1 = 300°) to C (ϕ = 270°/ψ = 240°/β1 = 60°) corresponds to the slower motion with relaxation time τtot = τ1 + τ2 ≃ τ1, the total activation energy is given by the higher free ≃ ΔHMTD = 30 ± 5 kJ/mol ≃ 12 kBT. energy barrier ΔGMTD tot 1 3.7. Metadynamics: Results in PA6I. The four degrees of freedom described previously (ϕ1, ψ1, ϕ2, ψ2) are involved in the motion of largest amplitude in PA6I. Using the same method as in PA6,6, we shall iteratively fix two degrees of freedom and look at the evolution of the other ones. Dipoles 1 and 2 on both sides of the phenyl ring cannot move simultaneously, which would correspond to a high energy barrier (see Figure 4). We shall thus first fix one dipole and investigate the motion of the other one. 3.7.1. ϕ1−ψ1 Free Energy Barrier. First, the angles ϕ2 and ψ2 are fixed to (ϕ2 = 90°, ψ2 = 180°), which corresponds to a most probable configuration; see Figure 5B), and the angles (ϕ1, ψ1) are chosen as CVs. From the initial conformation located at (ϕ1 = 270°/ψ1 = 0°) (conformation B′ in Figures 9A and 10), five local energy minima can be reached. The global minimum (most stable conformation) is observed for (ϕ1 = 90°, ψ1 = 180°) (conformation A′ in Figures 9A and 10), which also corresponds to the most remote conformation in terms of reorientation amplitudes as compared to the initial conformation B. The relative free energies of each conformations as compared to A′ are (ϕ1 = 90°, ψ1 = 0°): 12 kJ/mol ≃ 5 kBT; (ϕ1 = 180°, ψ1 = 0°): 6 kJ/mol ≃ 2.5 kBT; (ϕ1 = 270°, ψ1 = 0°): 2 kJ/mol ≃ 1 kBT; (ϕ1 = 180°, ψ1 = 180°): 17 kJ/mol ≃ 7 kBT; (ϕ1 = 270°, ψ1 = 180°): 17 kJ/mol ≃ 7 kBT. The conformational change from A′ to B′ is not direct and takes place through the

Figure 9. Energy surfaces in kJ/mol of PA6I obtained by metadynamics (MTD) with two collectives variables (CVs) ϕ1/ψ1 with ϕ2 = 90° and ψ2 = 180° (A) and ϕ2/ψ2 with ϕ1 = 270° and ψ1 = 0° (B). Arrows indicate large-amplitude motions discussed in section 5.

intermediate conformation I (ϕ1 = 90°/ψ1 = 0°). While the free energy barrier between A′ and I is 25 kJ/mol ≃ 10 kBT and that between I and B′ is 5 kJ/mol ≃ 2 kBT, the global free energy barrier from A′ to B′ is higher: 25 ± 5 kJ/mol ≃ 10 kBT. However, this simulation does not take into account ϕ2 and ψ2 dihedral angles, and their influence has to be determined. 3.7.2. ϕ2−ψ2 Free Energy Barrier. In a second step, the same initial state B′ (ϕ2 = 90°/ψ2 = 180°) as previously is considered. ϕ1 and ψ1 are fixed to (ϕ1 = 270°/ψ1 = 0°) (point B′ in Figures 9A and 10), and (ϕ2/ψ2) dihedral angles are chosen as CVs. The obtained free energy (enthalpy) map is shown in Figure 9B. From this initial conformation B′ (which is the most stable conformation), three minima are observed in the (ϕ2, ψ2) free energy map, with relative free energies as compared to B′ of respectively: (ϕ2 = 300°, ψ2 = 180°): 26 kJ/mol ≃ 10.5 kBT; J


Article

Macromolecules

Figure 10. Largest amplitude motions in PA6,6 and PA6I.

(ϕ2 = 90°,ψ2 = 0°): 22 kJ/mol ≃ 9 kBT; (ϕ2 = 300°, ψ2 = 0°): 16 kJ/mol ≃ 6.5 kBT (point C in Figures 9B and 10). The conformational change from B′ to C′ is not direct and may take place through either the conformation J1 (ϕ2 = 90°/ψ2 = 60°) or J2 (ϕ2 = 270°/ψ2 = 120°). Energy barriers of 40 kJ/mol ≃ 16 kBT are observed from B′ to either the conformation J1 or J2. Then, from either of these intermediate conformations J1 or J2, the system reaches C′ with a free energy barrier of 15 kJ/mol ≃ 6 kBT. The change from B′ to C′ is then controlled by the higher free energy barrier 40 ± 5 kJ/mol ≃ 16 kBT. 3.7.3. Free Energy Balance. Both MTD simulations can be combined to get a value of the total activation free energy, ΔHGMTD of the overall largest motion which occurs in the tot PA6I monomer, i.e., the conformational change from A′ to C′. Considering that the global reorientation from A′ to C′ is controlled by the slowest motion, (i.e., the highest free energy barrier), the overall activation free energy is ΔGMTD ≃ ΔHMTD tot 2 = 40 ± 5 kJ/mol ≃ 16 kBT.

associated with the considered motion is implicitly included in the value of the prefactor τ∞ as τ∞ = τ0 exp(−ΔS/R). Thus, experimentally, measuring a prefactor much smaller than 10−13−10−12 s indicates that significant entropy change (associated with complex conformational changes or change in the confinement state of short chain segments) is associated with the considered relaxation process. More precisely, the entropy generally decreases when the considered molecular variable is confined in a free energy well. Comparing MTD and experimental results may thus potentially allow discrimininating enthalpic and entropic contributions to local motions. Of course, as already mentioned, the proposed method is restricted to the study of the amorphous phase. In semicrystalline polymers, it is well-known that the dynamics of the amorphous phase may be strongly affected in the vicinity of crystallites. We shall see that some discrepancy between simulation and experimental results reported in semicrystalline materials may indeed come from this fact. Implementing the presence of a crystalline phase is a possible perspective of this work. Crystalline lamellae might indeed be modeled. However, modeling the confinement induced by the presence of crystalline lamellae and the related boundary conditions certainly involves a considerable further work. 4.1. γ-Relaxation in PA6,6 and PA6I. First, the dynamics of aliphatic segments have been investigated. Both the αi and β2 (PA6,6) dihedral angles are mainly in trans conformations in both aliphatic and semiaromatic PA6,6. The autocorrelation functions of all aliphatic dihedral angles (αi and β2 in PA6,6) have been calculated in MD simulations and fitted with monoexponential decays, which gives average values of the aliphatic relaxation times 10 ± 2 ns for PA6,6 and 12 ± 2 ns for PA6I. All aliphatic dihedral angles give similar values of the relaxation time. On the other hand, several experimental studies have determined, mostly by broadband dielectric spectroscopy, that the activation energy (enthalpy) ΔHγ of the γ relaxation process is between 28 kJ/mol ≃ 11 kBT and 46 kJ/mol ≃ 18.5 kBT, with an Arrhenius prefactor of the order 10−13−10−12 s.3,15,52−54 This would correspond to relaxation times between 10−7 and 10−8 s in normal conditions (T = 300 K, P = 1 bar), assuming an Arrhenius-like temperature variation τ = τ0 exp(−ΔG/kBT) with τ0 ≃ 10−12 s. Thus, the relaxation times

4. DISCUSSION In dry PAs, it is considered to be well established that the γ and β relaxation processes involve the aliphatic parts and the amide groups of the polymer, respectively.16−19 Here we discuss these assignments by looking at the different motions which have been studied using our multiscaling modeling approach. In MTD simulations, free energy (in fact enthalpy) barriers ΔG = ΔH − TΔS associated with the considered type of motion are measured, since MTD includes all local configuration changes associated with the transition between both equilibrium states. Indeed, even though they are very local, secondary relaxation processes, such as the β relaxation, may contain a significant entropic contribution. Then, from (molar) free enthalpy barriers ΔG estimated by MTD simulations, molecular relaxation times are deduced by assuming that relaxation processes are thermally activated, with a relaxation time varying as τ = τ0 exp(ΔG/RT) (where R is the gas constant), in which τ0 ≃ 10−12 s is the characteristic trial frequency at the molecular scale. On the other hand, the experimentally measured Arrhenius plot τ = τ∞ exp(ΔH/RT) gives the molar enthalpy barrier ΔH associated with the considered relaxation process, while the entropy change K


Article

Macromolecules

= 300 K, which is significantly faster than experimental data. Assuming τ0 to be unaffected by the crystalline state, this difference would correspond to a difference in free enthalpy barriers associated with the motion. One possible explanation would be that the mobility of the amorphous phase is affected in the vicinity of crystallites in PA6 and PA6,6 due to confinement or constraint effects. 4.3. β-Relaxation in PA6I. The same approach can be adopted in PA6I. In PA6I, the possible motion of largest amplitude accessed indirectly through MTD simulations corresponds to the conformational change from A′ to C′ in Figure 9, which was found to be controlled by an overall energy ≃ 40 ± 5 kJ/mol ≃ 16 kBT. The overall barrier ΔHMTD tot relaxation mechanism from A′ to B′ is schematized in Figure 10. In QM calculations, the energy barrier between A′ and (ϕ1 = 270°, ψ1 = 0°, ϕ2 = 90°, ψ2 = 180°) is 14 ± 2 kJ/mol ≃ 5.5 kBT, which is much lower than determined by MTD. As well as in PA6,6, dihedral motions can be controlled by hydrogen bond relaxation. A relaxation time of hydrogen bonds of 20 ps has been estimated by MD, corresponding to an energy ΔHMD H‑bonds = 11 kJ/mol ≃ 4.5 kBT. This relaxation time is smaller than in PA6,6, which may be due to the fact that CONH dipoles are smaller in PA6I than in PA6,6. By combining both conformaMD tional (ΔHQM 1 ) and interactional (ΔH1 ) parts, a global free QM MD energy barrier ΔGtot = ΔH + ΔH = 25 kJ/mol ≃ 10 kBT, 1 1 1 exactly the same as ΔGMTD = 25 ± 5 kJ/mol ≃ 10 kBT, is 1 obtained. Regarding the motion from conformation (ϕ1 = 270°, ψ1 = 0°, ϕ2 = 90°, ψ2 = 180°) to B′, the QM energy barrier is ΔHQM 2 = 27 ± 2 kJ/mol ≃ 11 kBT. By combining both conformational MD (ΔHQM 2 ) and interaction (hydrogen bonds) (ΔH1 ) parts, a tot QM MD global free energy barrier ΔG2 = ΔH2 + ΔH2 = 38 kJ/mol ≃ 15 kBT, comparable to ΔGMTD = 40 ± 5 kJ/mol ≃ 16 kBT, is 2 obtained. In PA6I, aromatic moieties may interact through π−π stacking.57,58 This interaction plays an important role in aromatic polymers, especially for thermomechanical properties. The autocorrelation function of π−π interactions along MD simulations has been calculated and fitted with a monoexponential decay. A relaxation time of 2 ps has been found, corresponding to an energy barrier ΔHMD Π−Π = 5 kJ/mol ≃ 2 kBT, with an estimated distance of 4 Å between the centers of mass of both aromatic ring. The barrier ΔHMD π−π may possibly contribute to the difference between ΔGtot and ΔGMTD . 2 2 However, the value 5 kJ/mol lies within error bars on those values, and thus, results are not precise enough to discriminate a possible contribution of π−π interactions. Whereas the global motion is governed by the slower ≃ ΔHMTD = movement, the total free energy barrier is ΔHMTD tot 2 40 ± 5 kJ/mol ≃ 16 kBT. The overall relaxation mechanism from A′ to C′ is schematized in Figure 10. Broadband dielectric spectroscopy performed on PA6I/6T (70/30) have shown that the β relaxation process follows an Arrhenius law with an energy Eβ = 100 kJ/mol ≃ 40 kBT, an Arrhenius prefactor τ∞ ≈ 10−18 s, and a relaxation time τ ≃ 10−3 s at T = 300 K.59 In the same way as in PA6,6, the difference between experimental and simulated barriers and the very small τ∞ value may correspond to the entropic contribution. Assuming τ0 = 10−12 s, τ ≃ 10−5 s at T = 300 K, which is faster than experimental data. Here again, the potential presence of small nuclei of crystalline phase may possibly shift the β relaxation.

of aliphatic carbons in PAs obtained by MD are of the same order of magnitude as experimental γ relaxation times. This suggests that it is indeed relevant to assign the γ process in PAs to the relaxation of aliphatic dihedral angles. Also, QM calculations indeed show that aliphatic parts have nonzero charges, thus effectively giving rise to a nonzero dielectric response. Finally, this also indicates that the entropic contribution to the γ process is negligible. 4.2. β-Relaxation in PA6,6. In PA6,6, the possible motion of largest amplitude involving the surrounding of an amide group is indirectly accessible through MTD simulations and corresponds to the conformational change from A to C in Figure 8. As discussed previously, this global motion is controlled by an overall free enthalpy barrier ΔGMTD = 30 ± tot 5 kJ/mol ≃ 12 kBT. The overall mechanism of A to C relaxation in PA6,6 is schematized in Figure 10. On the other hand, QM calculations show that the (intramolecular) energy barrier between conformations A and (ϕ = 90°,ψ = 180°, β1 = 180°) is 12 ± 2 kJ/mol ≃ 5 kBT. Similarly the energy barrier between (ϕ = 90°, ψ = 180°, β1 = 180°) and B is 12 ± 2 kJ/mol ≃ 5 kBT (Figure 3). Therefore, the QM energy barrier from A to B is ΔHQM = 12 ± 5 kJ/mol 1 ≃ 5 kBT. However, this motion may be hindered by the formation of (intermolecular) hydrogen bonds. The autocorrelation function of hydrogen bonds has been calculated along MD calculation and fitted with a monoexponential decay, giving a relaxation time of 100 ps, which would correspond to = 15 kJ/mol ≃ 6 kBT for hydrogen an energy barrier ΔGMD 1 bonds (estimated with τ0 ≃ 10−12 s and no entropic contribution). By combining both the conformational (or intramolecular), ΔHQM 1 , and interaction (or intermolecular), tot QM + ΔGMD 1 , parts, a global free energy barrier of ΔG1 = ΔH1 MD ΔG1 = 27 kJ/mol ≃ 11 kBT is obtained, which is indeed = 30 ± 5 kJ/mol ≃ 12 kBT. Then, the comparable to ΔGMTD 1 = 2 ± 1 kJ/ QM energy barrier for the B to C motion is ΔHQM 2 mol ≃ 1 kBT. This value, combined with the energy barrier of hydrogen bonds ΔHMD 1 = 15 kJ/mol ≃ 6 kBT, gives a global free QM + ΔHMD = 17 kJ/mol ≃ 7 kBT, energy barrier ΔGtot 2 = ΔH2 2 = 15 ± 5 kJ/mol ≃ 6 which is indeed comparable to ΔGMTD 2 kBT. Altogether, this shows that the global free energy barrier determined in MTD is mostly determined by interactional free energies, related to hydrogen bonds, rather than conformational energy barriers. From a structural point of view, the A to C conformational change in PA6,6 corresponds to a reorientation of the dipole moment of the amide group (CONH) (Figure 10). Other simulation results also described the CONH reorientation motion with an activation energy of 61 ± 5 kJ/mol ≃ 24.5 kBT.55 On the other hand, broadband dielectric spectroscopy performed on PA6 reveals that the β relaxation process follows an apparent Arrhenius law with an activation energy Eβ = 60 kJ/ mol ≃ 24 kBT, an Arrhenius prefactor τ∞ ≈ 10−15−10−16 s, and a relaxation time τ ≃ 10−5 s at T = 300 K.52 The β relaxation processes have been characterized recently in PA6,6 in the dry state as well as in the presence of water under various chemical activities.56 The same value τ ≃ 10−5 s is found at T = 300 K in the dry state. The difference between the experimentally measured energy barrier and the simulated free energy barrier, as well as the small experimental value of τ∞, would then correspond to the presence of a significant entropic contribution in the motion. Assuming an Arrhenius parameter τ0 = 10−12 s, corresponding to C−H bond vibration, we would find here τ ≃ 10−7 s at T L


Article

Macromolecules 4.4. Comparison between PA6,6 and PA6I. A difference of 10 kJ/mol ≃ 4 kBT has been found between the energy barriers of the slowest motion at the monomer scale in PA6,6 and PA6I. While interaction energies are roughly similar MD (ΔHMD H‑bonds = 15 kJ/mol ≃ 6 kBT in PA6,6 and ΔHH‑bonds + MD ΔHΠ−Π = 16 kJ/mol ≃ 6.5 kBT), the main difference comes from the conformational energy. According to this result, we therefore infer that the fact that the β relaxation process is slower in semiaromatic PAs (PA6I) as compared to aliphatic PAs (PA6,6) can be attributed to the increased backbone rigidity in semiaromatic PAs, which shifts relaxation times by about 2 decades between PA6,6 and PA6I. This result is in good agreement with broadband dielectric spectroscopy results.52,59 In aliphatic PA6, a relaxation time of 10−5 s is observed for the β relaxation at room temperature.52 In the case of semiaromatic PA6I/6T copolymer, a relaxation time of 10−3 s is obtained for the β relaxation at room temperature.59 Thus, the shift of relaxation time between aliphatic and semiaromatic PAs that we obtained is in good concordance with experimental one, showing that β relaxation is indeed slower in semiaromatic PAs. It therefore suggests that the β relaxation may indeed correspond to CONH dipole reorientation in PAs. The case of aliphatic/semiaromatic copolyamides may also be clarified by this method. It is important to notice that MTD simulations generally underestimate free energy (enthalpy) barriers due to the fact that punctual charges are kept constant along the simulations. While being aware of the limitations of such simulations, we claim that this method may be extremely powerful to understand semiquantitatively the dynamics of relatively largeamplitude motions in amorphous polymer materials. The new approach proposed here which cover multiple time scales may potentially allow going beyond existing simulation works. We therefore suggest that the β relaxation process in PA6,6 may be a combination of both (ϕ/ψ/β1) free rotations and hydrogen bond relaxation, while in PA6I, it may be a combination of both free rotations on each side of the phenyl ring (ϕ1/ψ1 and ϕ2/ ψ2) and relaxation of interactions consisting of both hydrogen bonding and Π−Π stacking. While this relaxation was previously attributed to γ relaxation processes due to the limitation of the MD time scale,55 we emphasize that CONH reorientation could rather correspond to β relaxation processes. This distinction is important, since γ and β relaxations correspond to two very distinct processes measured experimentally, which have different impact on the material properties.

around 10 ns in both aliphatic and semiaromatic PAs and have been proposed to describe the γ relaxation process. Amide relaxations show a relatively high activation energy for both aliphatic and semiaromatic PAs and have been proposed to correspond to the β relaxation process. These results are in good agreement with several experimental studies. 16−19 However, the amide relaxation energies are significantly different between PA6,6 and PA6I. Such a difference is assigned to the larger rigidity of polymer backbone in semiaromatic PAs. Also, in PA6I, the amide relaxation involves both amide groups cooperatively around aromatic rings. This study also emphasizes the important role of the aliphatic carbon in α position of the carbonyl groups in PA6,6. Whereas in PA6,6 the major part of intermolecular interactions are attributed to hydrogen bonds, in PA6I those are shared into hydrogen bonds and Π−Π stacking between aromatic groups. The different intermolecular interactions do not contribute to the difference observed between PA6,6 and PA6I regarding the activation energy of amide part (CONH), related to β relaxation. Nevertheless, the balance between the different types of intermolecular interactions involved in PAs (hydrogen bonds and Π−Π stacking) can play a significant role in terms of long time-scale dynamical process such as α relaxation. This particular aspect will be developed in a future work. Note finally that the method strongly depends on the energy path initially selected. In systems described by more than two degrees of freedom, fixing one parameter can lead to overestimating activation energies. It is therefore important to well select initial states and well determine all possible conformations at the MD step. A good knowledge of the studied system is required to identify conformations potentially involved in the different transitions. For a given transition, a large number of conformations in a system shall lead to a distribution of the total energy barrier. Indeed, in the literature, the β and γ relaxations are often described by broad relaxation time distributions.15 In our case, we have tried to identify the motion of largest amplitude in order to obtain the highest associated relaxation time. While being aware of the limitations of the method to reproduce the complex properties of polymers, we argue that it provides a predictive tool in order to understand secondary relaxations in amorphous polymers. Specifically, it may be extremely helpfull to compare secondary (subglassy) relaxations in different polymers of the same family. It is also important to mention that phenomena involving other, even longer, time scales, such as diffusion of small molecules in polymer matrices or the α relaxation, may be potentially modeled with this method. The method is of course limited to relatively local motions. Further theoretical development would be needed before being able to address large scale, cooperative motions associated with the glass transition (or equivalently α relaxation). Nevertheless, the approach proposed here may be a powerful predictive tool in order to tailor the relaxations in polymer materials, thus playing a key role in the design of innovative polymer materials with enhanced thermomechanical properties.

5. CONCLUSION The highlight of this innovative approach in polymers is the ability to simulating long time scale motions while retaining an all-atomistic description of the polymer material. Also, conformational (intramolecular, entropic) and interactional (intermolecular, enthalpic) contributions to the free energy landscape, and thus to energy barriers and time scales, can be discriminated clearly. Amorphous polyamides with different backbone flexibility, namely PA66 and PA6,I, have been modeled to estimate how chain conformations and intermolecular interactions vary between aliphatic and semiaromatic PAs. Aliphatic and amide relaxations have been highlighted using a multiscaling modeling approach combining ab initio calculations (QM), molecular dynamics simulations (MD), and metadynamics simulations (MTD). Aliphatic relaxations occur

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.macromol.5b01963. M


Article

Macromolecules

■

(23) Ahlrichs, R.; Bar, M.; Haser, M.; Horn, H.; Kolmel, C. Electronic structure calculations on workstation computers: The program system turbomole. Chem. Phys. Lett. 1989, 162, 165−169. (24) Fischer, C. The Hartree-Fock Method for Atoms: A Numerical Approach; John Wiley and Sons, Inc.: New York, 1977. (25) Hariharan, P.; Pople, J. The influence of polarization functions on molecular orbital hydrogenation energies. Theoretical Chemistry Accounts: Theory, Computation, and Modeling (Theoretica Chimica Acta) 1973, 28, 213−222. (26) Mulliken, R. S. Electronic Population Analysis on LCAOMO Molecular Wave Functions. III. Effects of Hybridization on Overlap and Gross AO Populations. J. Chem. Phys. 1955, 23, 2338−2343. (27) Goudeau, S.; Charlot, M.; Vergelati, C.; Müller-Plathe, F. Atomistic simulation of the water influence on the local structure of polyamide 6, 6. Macromolecules 2004, 37, 8072−8081. (28) Theodorou, D.; Suter, U. Detailed molecular structure of a vinyl polymer glass. Macromolecules 1985, 18, 1467−1478. (29) Berendsen, H.; van der Spoel, D.; Van Drunen, R. GROMACS: A message-passing parallel molecular dynamics implementation. Comput. Phys. Commun. 1995, 91, 43−56. (30) Lindahl, E.; Hess, B.; van der Spoel, D. GROMACS 3.0: a package for molecular simulation and trajectory analysis. J. Mol. Model. 2001, 7, 306−317. (31) Van Der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A.; Berendsen, H. GROMACS: fast, flexible, and free. J. Comput. Chem. 2005, 26, 1701. (32) Hess, B.; Kutzner, C.; Van Der Spoel, D.; Lindahl, E. GROMACS 4: Algorithms for highly efficient, load-balanced, and scalable molecular simulation. J. Chem. Theory Comput. 2008, 4, 435− 447. (33) Kaminski, G.; Friesner, R.; Tirado-Rives, J.; Jorgensen, W. Evaluation and reparametrization of the OPLS-AA force field for proteins via comparison with accurate quantum chemical calculations on peptides. J. Phys. Chem. B 2001, 105, 6474−6487. (34) Essmann, U.; Perera, L.; Berkowitz, M.; Darden, T.; Lee, H.; Pedersen, L. A smooth particle mesh Ewald method. J. Chem. Phys. 1995, 103, 8577−8593. (35) Hess, B.; Bekker, H.; Berendsen, H.; Fraaije, J. LINCS: a linear constraint solver for molecular simulations. J. Comput. Chem. 1997, 18, 1463−1472. (36) Kirkpatrick, S.; Gelatt, C., Jr.; Vecchi, M. Optimization by simulated annealing. Science 1983, 220, 671−680. (37) Cousin, T.; Galy, J.; Dupuy, J. Molecular modelling of polyphthalamides thermal properties: Comparison between modelling and experimental results. Polymer 2012, 53, 3203−3210. (38) Flory, P. J. Principles of Polymer Chemistry; Cornell University Press: Ithaca, NY, 1953. (39) Flory, P. J. Statistical Mechanics of Chain Molecules; WileyInterscience: New York, 1969. (40) Fox, T. G.; Flory, P. J. Second-order transition temperatures and related properties of polystyrene. I. Influence of molecular weight. J. Appl. Phys. 1950, 21, 581−591. (41) Saunders, P. The unperturbed dimensions of nylon 66. J. Polym. Sci., Part A: Gen. Pap. 1964, 2, 3765−3770. (42) Flory, P. J.; Williams, A. D. Configurational statistics of polyamide chains. Journal of Polymer Science Part A-2: Polymer Physics 1967, 5, 399−415. (43) Fetters, L. J.; Lohse, D. J.; Richter, D.; Witten, T. A.; Zirkel, A. Connection between Polymer Molecular Weight, Density, Chain Dimensions, and Melt Viscoelastic Properties. Macromolecules 1994, 27, 4639−4647. (44) Bonomi, M.; Branduardi, D.; Bussi, G.; Camilloni, C.; Provasi, D.; Raiteri, P.; Donadio, D.; Marinelli, F.; Pietrucci, F.; Broglia, R.; et al. PLUMED: A portable plugin for free-energy calculations with molecular dynamics. Comput. Phys. Commun. 2009, 180, 1961−1972. (45) Jorgensen, W. L.; Madura, J. D.; Swenson, C. J. Optimized intermolecular potential functions for liquid hydrocarbons. J. Am. Chem. Soc. 1984, 106, 6638−6646.

Table of atomic Mulliken charges from ab initio calculations; Newman projections of stable conformations in PA6,6 and PA6I; characteristic Flory’s ratio CM and radius of gyration as a function of monomer number M (PDF)

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected] (P.S.). Notes

The authors declare no competing financial interest.

■

REFERENCES

(1) Monnerie, L.; Halary, J. L.; Kausch, H.-H. Intrinsic Molecular Mobility and Toughness of Polymers I; Springer: 2005; pp 215−372. (2) Starkweather, H. W., Jr. Distribution of activation enthalpies in viscoelastic relaxations. Macromolecules 1990, 23, 328−332. (3) McCrum, N. G.; Read, B. E.; Williams, G. Anelastic and Dielectric Effects in Polymeric Solids; John Wiley & Sons Ltd.: New York, 1967. (4) Kremer, F.; Schönhals, A. Broadband Dielectric Spectroscopy; Springer: Berlin, 2003. (5) Williams, M. L.; Landel, R. F.; Ferry, J. D. The temperature dependence of relaxation mechanisms in amorphous polymers and other glass-forming liquids. J. Am. Chem. Soc. 1955, 77, 3701−3707. (6) Colmenero, J.; Arbe, A.; Alegria, A. The dynamics of the α-and βrelaxations in glass-forming polymers studied by quasielastic neutron scattering and dielectric spectroscopy. J. Non-Cryst. Solids 1994, 172, 126−137. (7) Schmidt-Rohr, K.; Spiess, H. W. Multidimensional Solid-State NMR and Polymers; Academic Press: London, 1994. (8) Monnerie, L.; Lauprêtre, F.; Halary, J. L. Intrinsic Molecular Mobility and Toughness of Polymers I; Springer: 2005; pp 35−213. (9) Bocahut, A.; Bernad, S.; Sebban, P.; Sacquin-Mora, S. Relating the diffusion of small ligands in human neuroglobin to its structural and mechanical properties. J. Phys. Chem. B 2009, 113, 16257−16267. (10) Leone, V.; Marinelli, F.; Carloni, P.; Parrinello, M. Targeting biomolecular flexibility with metadynamics. Curr. Opin. Struct. Biol. 2010, 20, 148−154. (11) Barducci, A.; Bonomi, M.; Parrinello, M. Metadynamics. Wiley Interdisciplinary Reviews: Computational Molecular Science 2011, 1, 826−843. (12) Elliott, J. Novel approaches to multiscale modelling in materials science. Int. Mater. Rev. 2011, 56, 207−225. (13) Laio, A.; Parrinello, M. Escaping free-energy minima. Proc. Natl. Acad. Sci. U. S. A. 2002, 99, 12562−12566. (14) Laio, A.; Gervasio, F. L. Metadynamics: a method to simulate rare events and reconstruct the free energy in biophysics, chemistry and material science. Rep. Prog. Phys. 2008, 71, 126601. (15) Laurati, M.; Sotta, P.; Long, D.; Fillot, L.-A.; Arbe, A.; Alegrı ́a, A.; Embs, J.; Unruh, T.; Schneider, G.; Colmenero, J. Dynamics of water absorbed in polyamides. Macromolecules 2012, 45, 1676−1687. (16) Owen, A.; Bonart, R. Cooperative relaxation processes in polymers. Polymer 1985, 26, 1034−1038. (17) Le Huy, H.; Rault, J. Remarks on the α and β transitions in swollen polyamides. Polymer 1994, 35, 136−139. (18) Diaz-Calleja, R.; Riande, E. Comparative study of mechanical and dielectric relaxations in polymers. Mater. Sci. Eng., A 2004, 370, 21−33. (19) Brandrup, J.; Immergut, E.; Grulke, E. A. Polymer Handbook, 4th ed.; CRC Press: 2003. (20) Kohan, M. I. Nylon Plastics Handbook; Hanser: Wilmington, 1995. (21) Scienomics, MAPS 3.3 platform; Paris, France, 2012. (22) Mayo, S.; Olafson, B.; Goddard, W. DREIDING: a generic force field for molecular simulations. J. Phys. Chem. 1990, 94, 8897−8909. N


Article

Macromolecules (46) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (47) Jorgensen, W. L.; Gao, J. Cis-trans energy difference for the peptide bond in the gas phase and in aqueous solution. J. Am. Chem. Soc. 1988, 110, 4212−4216. (48) Jorgensen, W.; Tirado-Rives, J. The OPLS [optimized potentials for liquid simulations] potential functions for proteins, energy minimizations for crystals of cyclic peptides and crambin. J. Am. Chem. Soc. 1988, 110, 1657−1666. (49) Roterman, I.; Lambert, M.; Gibson, K.; Scheraga, H. A comparison of the CHARMM, AMBER and ECEPP potentials for peptides. II. φ-ψ maps for N-acetyl alanine N-methyl amide: comparisons, contrasts and simple experimental tests. J. Biomol. Struct. Dyn. 1989, 7, 421−453. (50) Linge, J. P.; Williams, M. A.; Spronk, C. A.; Bonvin, A. M.; Nilges, M. Refinement of protein structures in explicit solvent. Proteins: Struct., Funct., Genet. 2003, 50, 496−506. (51) Xiang, Y.; Liu, Y.; Mi, B.; Leng, Y. Hydrated Polyamide Membrane and Its Interaction with Alginate: A Molecular Dynamics Study. Langmuir 2013, 29, 11600−11608. (52) Laredo, E.; Grimau, M.; Sanchez, F.; Bello, A. Water absorption effect on the dynamic properties of nylon-6 by dielectric spectroscopy. Macromolecules 2003, 36, 9840−9850. (53) Hirschinger, J.; Miura, H.; Gardner, K.; English, A. D. Segmental dynamics in the crystalline phase of nylon 66: solid state deuterium NMR. Macromolecules 1990, 23, 2153−2169. (54) Pathmanathan, K.; Cavaillé, J.-Y.; Johari, G. The dielectric properties of dry and water-saturated nylon-12. J. Polym. Sci., Part B: Polym. Phys. 1992, 30, 341−348. (55) Goudeau, S.; Charlot, M.; Müller-Plathe, F. Mobility enhancement in amorphous polyamide 6, 6 induced by water sorption: A molecular dynamics simulation study. J. Phys. Chem. B 2004, 108, 18779−18788. (56) Preda, F.-M.; Alegría, A.; Bocahut, A.; Fillot, L.-A.; Long, D. R.; Sotta, P. Investigation of Water Diffusion Mechanisms in Relation to Polymer Relaxations in Polyamides. Macromolecules 2015, 48, 5730− 5741. (57) Hunter, C. A.; Sanders, J. K. The nature of pi-pi interactions. J. Am. Chem. Soc. 1990, 112, 5525−5534. (58) Hunter, C. A.; Singh, J.; Thornton, J. M. π-π interactions: The geometry and energetics of phenylalanine-phenylalanine interactions in proteins. J. Mol. Biol. 1991, 218, 837−846. (59) Avakian, P.; Matheson, R., Jr.; Starkweather, H., Jr. Implications of dielectric response for the mechanism of changes in oxygen transport due to traces of moisture in amorphous nylons. Macromolecules 1991, 24, 4698−4700.

O


Modeling Molecular Relaxation Mechanisms in ... - ACS Publications

Recommend Documents