On the Size Effect of Additives in Amorphous Shape Memory Polymers

Small additive molecules often enhance structural relaxation in polymers. We explore this effect in a thermoplastic shape memory polymer via molecular dynamics simulations. The additive-to-monomer size ratio is shown to play a key role here. While the effect of additive-concentration on the rate of shape recovery is found to be monotonic in the investigated range, a non-monotonic dependence on the size-ratio emerges at temperatures close to the glass transition. This work thus identifies the additives’ size to be a qualitatively novel parameter for controlling the recovery process in polymer-based shape memory materials.


Introduction
Shape memory polymers (SMPs) are a promising class of functional materials with interesting properties and a wide range of applications [1,2]. Some of these features such as biocompatibility [3][4][5] and large strain recovery make them better candidates particularly in biomedical field compared to metallic shape memory alloys [6,7]. The thermo-responsive shape memory effect (SME) is an intrinsic feature among almost all the polymers due to their entropic elasticity, that is, a tendency of the stretched chains to deform back towards the most probable conformation [8]. However, exploiting the entropic elasticity in order to generate a usable SME depends on the molecular architecture of polymer as well as an appropriate programming procedure [7,9].
Depending on the underlying mechanism, the shape recovery in SMPs can be triggered in different ways such as by thermal, electrical and magnetic stimuli, as well as by light, ultrasound energy, moisture and change in PH [10][11][12][13][14][15][16][17]. Crucial to shape memory effect in polymers is the change of chains' conformation upon deformation and the tendency to recover a high entropic state during structural relaxation. These high entropic conformations have been recently studied for quiescent copolymers of different architectures at low temperatures and under poor-solvent conditions [18][19][20].
In most of the biomedical applications for SMPs, the sample shall operate and remain effective in a chemical environment [21][22][23][24][25][26]. Therefore, a question of high importance concerning the functionality of SMPs in physiological environments is how efficiently they can perform when small molecules are absorbed into their structure from the environment.
Small molecules are reported to influence morphology [27] and relaxation dynamics [28][29][30][31][32][33][34][35][36] in polymers. The glass transition temperature, T g , of the polymers decreases upon absorbing small molecules [31]. It is reported that this effect can also lead to an antiplasticization by enhancing the local stiffness [32]. Furthermore, antiplasticizers are found to alter the nature of glass formation via enhancing the packing efficiency in polymers [28] which leads to the formation of stronger glass-forming materials [30]. In thermoplastic SMPs, the additive-effect often manifests itself in a shift of the T g and the onset of shape recovery at a lower temperature as compared to the pure sample [37]. This solvent-induced effect on SME has been investigated experimentally by many researchers [37][38][39][40][41][42][43][44][45][46][47][48]. A recent experimental study reported on a new strategy to tune the switching temperature through sequence-rearranged cocrystallization of copolymer blocks [49].
The shape memory behavior in polymers has also been investigated by means of molecular dynamics (MD) computer simulations [50][51][52][53][54][55][56]. In a recent MD-work [54], different models of polymers were analyzed in a shape memory cycle, highlighting the importance of temperature-dependent chain mobility. Another study [56] employed dissipative particle dynamics (DPD) to investigate phase-separation in the polymer network and individual functions of different segments in a polyurethane-based shape memory polymer. As to additive-induced shape recovery, it has been shown that plasticizing effect of additive molecules on structural relaxation in polymers leads to an enhancement of shape recovery rate [53,57,58].
In the present work, MD simulations are employed to explore additive effects on shape memory behavior with a qualitatively new perspective, showing that the size of absorbed molecules provides a new independent parameter to control the shape recovery process. Importantly, at a constant number concentration, the temperature at which significant relaxation towards the original shape takes place changes in a non-monotonic way with the additive diameter.

Model
In this study, we use a toy model made of a sequence of strongly and weakly interacting spherical particles which we loosely name 'hard' and 'soft', arranged in a way as to mimic the sequence of hard and soft segments in ESTANE ® . The model consists of ar block copolymers with two monomer types A and B. The connectivity of the chains' backbone is ensured via a FENE potential [59,60], where k = 30 is the strength factor and R 0 = 1.5 is the breaking limit of covalent bonds. In addition to the FENE-potential, all particle pairs (including spherical additive molecules) interact via a Lennard-Jones (LJ) potential, In Equation (2), α, β ∈ {A, B and S}, where the letter S stands for spherical additive particles. r αβ is a short hand notation for the distance between a particle, i, of type α and another one, j, of type β: r αβ = | r i,α − r j,β |. The LJ potential is truncated at a cutoff radius of r c,αβ = 2 × 2 1/6 σ αβ . The monomer diameter is set to σ AA = σ BB ≡ 1 and serves as unit of length. For α = β, the parameter σ αβ is chosen as the arithmetic mean, σ αβ = 0.5(σ αα + σ ββ ).
As mentioned above, A and B monomers follow the same sequence as the hard and soft segments in ESTANE ® (Figure 1a,b). Using the LJ-interaction parameter of AA-particle pairs as unit of energy, ε AA = 1, the BB-interaction strength is set to ε BB = ε AA /2 = 0.5. The corresponding cross term of interaction energy is chosen as ε AB = 0.3ε AA = 0.3. The justification for this choice of interaction parameters is that they lead to a phase separation between hard and soft segments on the nano-scale (Figure 1c-e), and thus reproduce an important feature-which is well-known from experiments [61][62][63]-at least on a qualitative level. We also remark that bond-angle and torsional potentials are not included in the present model. Results discussed in this manuscript thus correspond to a generic flexible co-polymer model without reference to a specific chemistry.
The phase separation of hard and soft regions is a slow process, whose simulation requires advanced techniques. Here, we equilibrate the system using an approach proposed by Parker and Rottler [64], where the excluded-volume LJ-interaction is replaced by a soft potential. A main feature of this approach is that the chains become crossable and exhibit a Rouse-like dynamics which leads to a considerably faster relaxation.
Since the thus obtained dynamics is a biased one, we use this method only for the sample preparation step (from a random initial configuration towards a nanostructured sample). All the simulations of the shape memory cycle are performed using the full interaction model described through Equations (1) and (2). If additive particles are included, the energy scales pertaining to the corresponding spherical particles are all set to that of the soft monomers (ε SS = ε AS = ε BS = ε BB ). The employed parameters of LJ potential for different components are summarized in Table 1. The mass of an additive particle is set to be equal to that of the monomers, m S = m A = m B = 1. Temperature is measured in units of ε AA /k B with the Boltzmann constant k B (≡ 1). All other quantities are given as a combination of the above described units. The unit of time, for example, is given by τ LJ = (mσ 2 AA /ε AA ) 1/2 and that of pressure is ε AA /σ 3 AA .
All quantities addressed here are expressed in this set of reduced LJ units. The equations of motion are integrated using the Velocity-Verlet algorithm with a time step of δt = 0.005 for the preparation of initial configurations and δt = 0.003 for the actual simulations of shape memory behavior. Unless otherwise stated, all the simulations of the SMP-behavior are conducted at a constant pressure of p = 0. For the simulation software, we use the open source molecular dynamics simulator LAMMPS [65]. All the 3D snapshots of particle configurations in this work are prepared by the visualization software OVITO [66].

Sample Preparation
For the preparation of initial configurations, first, a regular arrangement of the chains is simulated within a constant volume and at a fixed temperature (NVT-ensemble) for a duration of 150,000 τ LJ at a monomer number density of ρ 0 = 0.93 σ −3 AA and temperature T = 1 ε AA /k B employing the above-mentioned soft repulsive potential [64]. This allows a considerably faster equilibration of our long-chain polymeric system since entanglements are avoided and the monomers follow a Rouse-like dynamics. The mean-squared displacements and end-to-end autocorrelation function plotted in Figure 2a-d demonstrate that the soft potential version of the system is completely equilibrated. Then the soft potential is reversed to the full LJ model with excluded-volume effect following the procedure explained in Reference [64]. Finally, the sample is cooled down to T = 0.35 with a cooling rate ofṪ = 10 −4 using N pT-ensemble at zero pressure.
The polymer-additive mixture for additives of σ SS = 0.5 is prepared by inserting, one after another in regular time intervals, spherical particles into the polymer matrix, making use of the soft potential U soft (r αS ) = ε αS 1 + cos πr αS σ αS , with ε αS = 0.5 and σ αS = 0.5(σ αα + σ SS ), where α ∈ {A,B}. The number concentration is defined via φ s = N s /(N s + N p ), with N s and N p being the number of additives and monomers. Particle insertion is performed during a time interval of 1500 τ LJ at T = 0.35 in N pTensemble at zero pressure. Then, the sample is simulated for 500 τ LJ which is sufficient for the additives to travel a relatively large distance of 20 particle diameters thanks to their soft interaction (cosine potential) with the surroundings. By surveying concentration profiles along different spatial directions and also within various planes, we have verified that additives reach a uniform distribution before switching on the full LJ-potential (for brevity, data are not shown). After these steps, the energy parameter of the soft potential is gradually increased to a high value of ε αS = 30 during a time interval of 100 τ LJ . This reduction in softness of interactions serves to prevent the additives from overlapping with the surrounding particles and is necessary for avoiding high energy collisions when the LJ is switched on. Subsequently, the cosine potential is replaced by the LJ potential and the system is simulated for a time span of 9000 τ LJ . To generate mixtures with additives of other diameters (e.g., σ SS = 0.3 and 0.8), we do not repeat the above procedure but use the σ SS = 0.5-configurations as the starting point and gradually change, at a constant pressure of p = 0, the additive diameter σ SS from 0.5 to the desired value. (e) Mean-squared displacements of different components plotted for the full model in the absence of additive at T = 0.23 (<T g ≈ 0.25, see Figure 3). The emergence of a quasi-Rouse-like dynamics, is visible (∼t 0.63 ). The symbols show dynamic measurements using t 0 = 10 5 τ LJ (note that MSD are plotted versus t − t 0 ). The observed difference in the dynamics for t 0 = 0 (lines) and t 0 = 10 5 τ LJ (Symbols) is a signature of an ongoing aging. This effect is, however, weak and can be neglected as compared to the effects of temperature and additive molecules. (f) Pair distribution function for hard and soft subchains, averaged over a first time interval of t ∈ [0, 10 5 τ LJ ] (solid lines) and a second one t ∈ [10 5 τ LJ , 2 × 10 5 τ LJ ] (dashed lines), highlighting that the structure is hardly changed during aging. An anti-correlation is visible between the self and cross terms. This can be regarded as a consequence of the phase separation between the regions of hard and soft segments (Figure 1c-e). The thus obtained samples are simulated for another 10,000 τ LJ in N pT-ensemble at zero pressure. This procedure, of course, does not lead to a perfectly equilibrated configuration of the full model . Rather, it slows down the structural relaxation and leads to a plateau in mean square displacements, signaling temporal arrest in the nearest neighbor cage (see, e.g., Figure 2e for a pure sample quenched further to T = 0.23 atṪ = 10 −4 ). The structural data further underlines the disordered character of the system (Figure 2f). As a crucial aspect, we make sure that the duration of temporal arrest for all the samples prepared this way is beyond the time relevant for the main part of strain recovery, which is ≤1000 τ LJ . This way, the shape memory effect is not deteriorated by structural relaxation.

Simulating the Shape Memory Effect
The samples obtained in the previous section are used as initial configurations for studying the shape memory effect via the following procedure (see Figure 4 for an illustration). First, the simulation box is stretched along the z-direction with a constant engineering strain rate (˙ = 10 −3 ), defined via = dL z /L z (t = 0), at a temperature above the glass transition point (T deform > T g ) of the sample under investigation.
It is noteworthy that, due to the presence of two distinct domains with different self-interactions, that is, hard and soft segments, the model exhibits two glass transition temperatures (Figure 3). In a simplified view, the higher T g can be assigned to a transition from a rubbery state to a partly glassy state, where hard segments are essentially immobile but soft segments keep significant mobility. As a consequence, shape memory effect does not occur at temperatures around this high T g . The lower T g , on the other hand, represents a temperature below which both hard and soft domains have negligible dynamics (are temporally arrested) within the deformation time scale. Therefore, we performed all our simulations at temperatures around and below the lower T g . Since the higher T g plays no role throughout this work, we refer to the lower value of T g ≈ 0.25 as the glass transition point of the pure polymer system. Similarly, glass transition temperatures of polymeradditive mixtures refer to the lower value of the two T g s.
During deformation, the two dimensions transverse to the direction of deformation, L x and L y , are independently coupled to the Nosé-Hover barostat at zero pressure (p ⊥ = p xx = p yy = 0) and can vary in response to internal forces. Once a maximum strain of max = 100% is reached, the sample is cooled at fixed strain to a temperature of T = 0.1 with a cooling rate ofṪ cooling = 10 −2 . Applying such a fast cooling helps to freeze chain conformations in the elongated state and thus minimize the undesired effects which would arise from structural relaxation. At the end of this deep and fast cooling process, the size of the simulation box along the direction of deformation, L z , is coupled to a barostat at p = p zz = 0, while the lateral dimensions (L x and L y ) remain independently under stress-free conditions. This leads to the release of the stress and is, therefore, referred to as the unloading stage. The above discussed temporal arrest ensures that, in the absence of external load, the system still remains in the elongated state (with the exception of a small spring-back effect due to ordinary elasticity). From the beginning of deformation process up to this point is often referred to as the 'programming' step. Shape recovery process is triggered via heating to a temperature above the sample's glass transition point, T recovery > T g , where structural relaxation becomes sufficiently fast to activate entropic forces, which drive the shape recovery process in polymers. For all the simulations, this heating is applied at a relatively high rate (Ṫ heating = 10 −2 ) to suppress effects of relaxation processes, which start to become active as temperature increases. Figure 4. A shape memory cycle for the present bead-spring copolymer model. The sample is first deformed at a temperature above its glass transition, T g . Then it is cooled down to a temperature below T g followed by unloading. The sample first deforms back slightly due to its energetic elasticity but the major part of deformation is stored as it is connected to entropic elasticity, the latter becoming active at temperatures above T g . Therefore, shape recovery takes place only if the sample is heated to a temperature T > T g .

Results and Discussion
In this section, the following findings shall be demonstrated: (1) The simple copolymer model presented above is capable of entropy-driven shape memory behavior, (2) adding small molecules to the sample plays a role similar to increasing temperature with regard to shape recovery and (3) the temperature, at which significant shape recovery takes place, first decreases for smaller additive particles but then increases as the particle size falls below roughly half a monomer diameter.

Effect of Additive Concentration on Shape Recovery
A way to highlight the entropic effect in shape memory polymers is to monitor the chain conformation. For this purpose, we survey the dynamics of the so-called gyration tensor of a chain, defined as S αβ = 1 Here, Greek letters denote spacial directions, α, β ∈ {x, y, z}, N p is the number of monomers (or segments) in a single chain (including both hard and soft segments) and R cm,α = 1 N p ∑ N p i=1 r i,α is the αth-component of chain's center of mass position. The symbol · · · stands for statistical average.
A first set of results obtained within the present model is shown in Figure 5, for the cases of a pure sample and a polymer-additive mixture. In the both cases investigated, deformation starts at a temperature slightly above the glass transition point ( where the system is soft enough to prevent it from crazing [67,68]. The plot depicts projections of the chains' gyration tensor onto the directions parallel and perpendicular to deformation, R g,⊥ ≡ 1 2 S xx + S yy and R g, ≡ S zz . One can clearly see from Figure 5 the increase in anisotropy of chain conformation during the externally imposed deformation (see data for t < 1000 τ LJ ). As the load is removed, the chain conformation still remains in this out-of-equilibrium state as long as temperature is kept sufficiently low (T = 0.1 T g ). Upon heating to a temperature above T g (T recovery = 0.27), however, R g, starts to decrease while at the same time R g,⊥ increases, both converging towards one another with time, thus signaling the evolution of chain conformation towards an isotropic distribution. It is, in fact, this spontaneous process of regaining isotropic conformation, which drives shape recovery. Importantly, Figure 5 provides evidence that, within a shift in recovery temperature, the same entropy driven scenario is prevalent both in a pure polymer sample and in a probe containing additive molecules. and perpendicular (⊥, dashed lines) to the direction of deformation. Thick blue lines correspond to a pure sample and green lines to a probe containing 20% additives with σ SS = 0.8. End of deformation is marked by the leftmost vertical line. After the cooling stage (between the 1st and 2nd vertical lines), L z is coupled to a barostat with p zz = 0. This leads to a small partial recovery, followed by a plateau (best visible in the inset) indicative of temporal arrest (shape fixation). As temperature raises above T g , chains become sufficiently mobile to acquire isotropic conformation and the probe relaxes to its original shape. The third and fourth vertical lines in the inset mark the beginning and end of the heating process. For better visibility, T is multiplied by a factor of 40. The upper right image is a snapshot of an undeformed configuration at T = 0.35 with 20% additive particles (yellow). Monomers (blue) are made semi-transparent.
In view of the close connection between shape recovery and chains' mobility, the upper temperature to which the sample is reheated plays a major role for the recovery rate. To address the influence of this 'recovery temperature', T recovery , statistically equivalent copies are prepared in exactly the same way until the end of programming step. Then, temperature is raised to different T recovery -values, all above T g . As shown in Figure 6a, rate of shape recovery is enhanced at higher T recovery , consistent with the idea that kinetics of chain relaxation is faster at a higher temperature. Similarly, since adding small molecules also affects the chains' relaxation dynamics, results for a number of additive concentrations are shown in the same plot. Similar to a raise in T recovery , one observes an enhancement of recovery rate upon increasing the number concentration of additive molecules. A closer scrutiny of the data shown in Figure 6a reveals that recovery profile of a pure sample at T recovery = 0.32 is quite similar to that a polymer-additive blend at a concentration of φ s = 20% but a lower temperature of T recovery = 0.30. Clearly, the initial strain relaxation is very similar in both cases and extends down to a strain value of roughly 0.2, below which slight deviations become visible. If one is interested in the first 70-80% of strain recovery, additive essentially play the role of a higher temperature. This interpretation is confirmed further by observing that a larger T-reduction (from T recovery = 0.34 to T recovery = 0.30) is possible if a higher additive concentration (φ s = 0.30) is used. These findings are in qualitative agreement with experimental report on the effects of moisture and organic solvents on shape recovery in SMPs [38][39][40]42,58,69].
A more quantitative analysis of this observation is provided in Figure 6b, where the initial recovery rate-defined via˙ 0, = (t) − (t 0 ) (t − t 0 ) with t − t 0 = 150τ LJ -is depicted for a number of concentrations and three different values of T recovery . The plot also highlights the above mentioned idea that the same recovery rate can be realized by various combinations of temperature and additive-concentration.
An important consequence of this observation is that additive molecules provide an independent way to trigger the shape recovery process, since one can keep T constant but allow for penetration of additive molecules into the sample. While this finding is in line with previous studies [36,47,53], the effect of additive-to-monomer size ratio on shape recovery process has not yet been explored in the literature. The next subsection is devoted to this novel aspect.

Effect of the Additive Size
It has been shown above that the plasticizing effect of small molecules, which penetrate into the polymer sample, remains intact in a shape memory cycle to the extent that, additives' concentration can be used as a new control parameter to influence the shape recovery process. Here, we demonstrate a further and qualitatively new way of using additive molecules to control the shape memory effect. This second independent route is related to the molecules' size relative to a monomer unit.
To focus on this new effect, we follow the same procedure as described on page 6 and prepare polymer-additive mixtures for different sizes of additive particles. In all the samples investigated below, the number concentration is kept constant, φ s = 0.20. First, a series of simulations are performed for different probes using exactly the same shape memory protocol with T deform = T recovery = 0.23. This temperature is slightly below the glass transition point of the pure sample (T g (φ s = 0) ≈ 0.25).
Strain profiles obtained from these studies are shown in Figure 7. As expected, under load, the parallel component of strain, , evolves in identical way for all the cases shown. As the load is released, strain follows different paths depending on the size of molecules added to the polymer. A close survey of the data reveals that the presence of additive molecules reduces the early elastic response, which occurs right after load removal. This effect becomes stronger between additive sizes σ SS = 0.3 and σ SS = 0.5 but then seems to saturate with further increase of σ SS . Figure 7. Additive size effect on shape recovery. The plot shows strain versus time in a shape memory cycle for a pure polymer and polymer-additive mixtures at φ s = 0.20 for three choices of additive's diameter, σ SS , as indicated. The recovery temperature is T recovery = 0.23. A close scrutiny of the data reveals a non-monotonic dependence of the shape recovery process on particle size, the largest recovery occurring for the sample with σ SS = 0.5. However, the most prominent effect of molecular size shows itself in the early stage of shape recovery, where a non-monotonic dependence of recovery rate on σ SS becomes visible (Figure 8a). Noteworthy, this non-monotonic effect is a specific feature of low temperature response, as the rate of strain decrease is a monotonic function of σ SS at higher temperatures. It is, therefore, intriguing to check whether this effect extrapolates to lower temperatures to the point of being observed in the glass transition temperature. To examine this idea, we have determined T g via cooling simulations at a constant additive concentration but different additive-diameters. These simulations provide evidence that the glass transition temperature depends on the additive diameter in a non-monotonic way (Figure 8b).
(a) (b) Figure 8. (a) Initial shape recovery rate versus σ SS for three values of recovery temperature as indicated and for a constant number concentration of φ s = 0.20. At the highest temperature shown here (T recovery = 0.35), the rate of strain recovery is a monotonic function of σ SS but a non-monotonic trend develops for lower T recovery and becomes distinguishable at a temperature of T recovery = 0.23 < T g (φ s = 0) ≈ 0.25. Filled symbols on the right vertical axis give the recovery rate of the pure sample at the corresponding T recovery . (b) The glass transition temperature, T g , obtained from cooling simulations (Figure 3) versus additive-diameter. In agreement with low-T data in panel (a), a non-monotonic effect is clearly visible here. All quantities are given in Lennard-Jones units.

Summary and Outlook
In this work, a simple bead-spring type copolymer model is used to study shape memory effect in AB-copolymers. The model is inspired by ESTANE ® , but has no spe-cific chemistry. Rather, A and B monomers are arranged in the same sequence as the so-called hard and soft segments occur in the real material. Regardless of this simplification, the model captures well the shape memory behavior. This fact is traced back to the central role of entropic elasticity and the possibility to freeze non-equilibrium (elongated) chain conformations via cooling below the glass transition.
Using this model, the effect of added spherical particles on shape recovery process is investigated and is shown to be qualitatively similar to that of temperature. Plasticizing effects of additives on shape memory polymers are also reported from experiments [38][39][40]42,58,69]. The fact that a simple copolymer model reproduces such interesting effects hints towards the generic character of the phenomenon, encoded in frozen entropic elasticity.
In addition to these findings, it is shown that at temperatures close to the glass transition point, the rate of shape recovery first increases but then decreases again for progressively larger additive particles. This non-monotonic effect is shown to be also present in the dependence of the glass transition temperature on additive diameter.
Non-monotonic effects result from competing trends which develop in opposite ways as a control parameter, for example, temperature, is varied [70][71][72]. Recently, a mechanism was proposed to rationalize such non-monotonic behavior based on the competition between higher mobility of smaller additive molecules versus stronger additive-host coupling strength for larger ones [73]. A key aspect here is the strong sensitivity of structural relaxation dynamics to tiny changes of local particle packing in the proximity of T g [74]. This issue is also highlighted in recent theoretical developments [75,76] and experiments [77,78] on binary colloidal mixtures-where packing effects are most prominent. The connection between local packing and structural relaxation dynamics is, however, a subtle one. Importantly, we do not observe a non-monotonic variation of free volume with additivesize, which could serve as a direct cause of the observed phenomenon. For more details, the interested reader is referred to the related literature [73][74][75][76][77][78].
Regarding further theoretical work, the non-monotonic size effect shown in this work deserves special notice, as it occurs in polymers' response to large deformation strains. This is a highly non-trivial finding and motivates new theoretical work. In view of the evidences obtained from the studies of various quiescent systems, it seems promising if one could develop a theory for polymers-additive systems which accounts both for local packing features and at the same time considers large deformation strains.
It must also be emphasized that the reported size dependence is a subtle effect and needs being examined in the case of a specific chemistry. It would be very interesting to check this issue using more realistic atomistic models. This would also open the way to investigate molecular coupling mechanisms between additive molecules and polymer and thus better understand how small molecules influence the structure and dynamics of the polymer matrix and hence the shape memory effect.
Our results also motivate new experiments with a focus on particle size effects. A possibility here would be to use small organic molecules of similar structure but variable sizes. Another option could be the use of different ideal gases-as good model systems for spherical particles with variable sizes-whose absorption into the polymer sample is driven by entropy [79] and which have relatively weak energetic interactions with the polymer matrix.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.