Coarse-Grained Molecular Dynamics Modeling of Segmented Block Copolymers: Impact of the Chain Architecture on Crystallization and Morphology

: We extend our recent coarse-grained model describing semicrystalline homopolymers to simulate the morphol-ogy and phase transitions of thermoplastic elastomers made of segmented (hard/soft) block copolymers. The generic model is adapted to match the physical characteristics of the two chemical units involved in the copolymer chains by using classic scaling rules. We investigate the crystallization kinetics of the hard segments as well as their phase separation from the soft units in either triblock or pentablock copolymers. We identify the soft segment molecular weight as a key parameter resulting in the following observations when decreasing the temperature from a homogeneous state. On the one hand, the phase separation preceding the crystallization process in triblock copolymers results in a constant temperature of crystallization when varying the soft segment length. On the other hand, the limited phase separation achieved in pentablock copolymers constrains them to crystallize at progressively lower temperatures while increasing the soft segment length. Finally, increasing the soft segment molecular weight was found to lead to a higher relative crystallinity which can be interestingly related to a rise of the loop segment’s content. bridges two di ﬀ erent crystallites, “ loop ” if it is connected twice to the same crystallite, “ dangling end ” if it has only one extremity connected to a crystal, and “ purely amorphous ” if the whole chain has no crystallized segment. This analysis is performed every 0.1 ε u / k B during the slow cooling process. The MSID of the di ﬀ erent types of segments is also computed during this ramp. Finally, the morphology of the multiphasic material is studied by using the radial distribution function (RDF) computed with the Ovito software. 78 To analyze the phase separation process, the RDF of HS beads is computed during the cooling ramp (having removed the SS beads from the box to keep the hard phase only). The resulting morphology after crystallization of the hard segments is also analyzed by computing the RDF on the crystallized HS beads only, at a ﬁ xed temperature T = 1 ε u / k B .


INTRODUCTION
Thermoplastic elastomers (TPE) are widely used in the industry since the 1960s for their rubbery properties and ease of processing. They can be found in various applications such as dashboard elements, wires sheaths, and bitumen modifiers. 1−4 Their outstanding properties come from the specific chain morphology of these segmented block copolymers consisting of an alternation of soft and hard segments (respectively denoted SS and HS below) along the chain. 5 While SS are made of an amorphous rubbery polymer at room temperature, HS are characterized by a higher glass transition temperature 6 and a possible ability to associate with each other by crystallizing. 7 This particular chain morphology leads to a microphase-separated structure of soft and hard domains, the latter being glassy regions or crystallites. These hard domains behave like physical nodes linking the soft segments, in contrast to chemical cross-links in vulcanized rubbers. However, the structure−properties relationship of thermoplastic elastomers is still unclear. Experiments have pointed out that these materials can exhibit quite different morphologies according to their chemical composition and processing route, 5,8 suggesting that interaction strengths between the blocks, chain rigidity, and thermal history are determining factors.
Actually, while the microphase separation and ordering processes in amorphous diblocks and triblocks are well understood with clear conclusions in both experimental 9,10 and theoretical 11−13 approaches, increasing the number of blocks in the chains makes their behavior significantly more complex. Theoretical works show that multiblock copolymers can present the same types of phases as in diblocks. 14,15 Nevertheless, structural characterizations also show that changing the number of blocks can be a useful tool to tune the phase diagram 4,16−18 and the order−disordered transition, sometimes called the microphase separation transition (MST). 19−21 This transition is known to depend on four parameters: the interaction parameter, the fractions of each component, the chain length, and the number of blocks. 18,21 However, interpreting the consequences of additional blocks in the chain can sometimes be difficult because in most of these experimental studies both the number of blocks and the chain length change simultaneously. While increasing the chain length promotes the phase separation at higher temperatures, increasing the number of blocks tends to lower the MST temperature. This can be understood from the thermodynamic point of view and has been modeled by Krause considering surface energies at the junction of two consecutive blocks. 22 The full understanding of these ordering mechanisms becomes even more intricate when dealing with semicrystalline multiblock copolymers. The thermal behavior of such systems becomes quite complex with several phase transitions: glass transition of each block, melting/crystallization of each crystallizable block, and microphase separation transition. The morphologies obtained after crystallization are highly dependent on the relative values of these phase transitions temperatures. 23 When the melting point of the crystallizable block exceeds the MST temperature, the crystallization can proceed from the homogeneous melt. 24 In this case, crystallization can even induce the phase separation. 25 On the contrary, when the microphase separation transition occurs above the melting point, ordering happens before crystallization. If the glass transition of the amorphous segment is higher than the crystallization temperature of the other block, then the crystallization is forced to occur between the hard glassy domains. 26−29 However, this is not the case for TPE. Thus, crystallization can reorder the melt morphology if the segregation strength is too weak. 30 A crystallization confined to the phase-separated domains requires a higher segregation strength. 31,32 In addition, the link between the structure and the mechanical behavior still remains elusive. Very few models can explain the role of the SS conformations on the rubbery modulus of these materials. 33 A reason for that is the difficulty to get such reliable topological information. To get a better understanding of this complex structuring behavior, one can use a simple coarse-grained molecular dynamics modeling. As mentioned above, TPE are made of segmented copolymers with two types of blocks (SS and HS), each one having specific chain rigidity, density, interaction strength, and ability to crystallize. These four criteria must be taken into account to correctly describe the TPE morphology and thermomechanical behavior. Generally, coarse-grained molecular dynamics simulations of block copolymers concern amorphous systems, with glassy and rubbery domains. 34 These can be obtained by defining several regions containing each type of block 35 or by starting from a homogeneous melt which will phase separate when decreasing the temperature. 36−39 It has been demonstrated that such models can provide a good description of the nonlinear tensile response of amorphous thermoplastic elastomers like SIS (styrene−isoprene) or SBS (styrene− butadiene) triblock copolymers. 40 Similar numerical models have been used to study the microphase separation transition and the resulting patterns for directed self-assembly systems. 41−43 Molecular dynamics simulations also demonstrated the key role of the copolymer architecture on the resulting morphology in solvent. The influence of the SS length 44 and grafting 45 have notably been investigated. Single c h a i ni nm e a nfield simulations have led to similar conclusions. 46 Besides, a recent study has shown the modification of the chain conformations when adding associating groups to the copolymer, resulting in the formation of intramolecular loops. 47 Apart from that, no molecular dynamics model describing the morphology and phase transitions in crystallizable segmented block copolymers has been reported to our knowledge.
Nevertheless, several coarse-grained models enable the simulation of semicrystalline homopolymers. Some of them are directly derived from atomistic simulations using a systematic coarse-graining procedure, as proposed for poly-(vinyl alcohol) (CG-PVA). 48−51 The relationship between the semicrystalline structure obtained with this CG-PVA model and its mechanical behavior has been extensively studied by Jabbari-Farouji and co-workers. 52−54 Another approach consists of using a generic coarse-grained model, such as the wellknown Kremer−Grest (KG) bead−spring model, 55,56 and modifying its parameters to obtain an energetically favorable crystalline structure. 57,58 In such a model, beads represent groups of atoms interacting through two types of forces: covalent bonds and weak interactions. Each interaction is modeled by a potential leading to two equilibrium distances: the bond length l b and the distance between nonbonded beads r eq (see details in section 2.1). Hoy et al. showed that crystallization could be obtained by equalizing these two equilibrium lengths: l b = r eq . 57 However, this results in a very coarse crystallization, in which the beads are aligned in space to form crystals, but with unrealistic chain conformations. A very recent publication based on a similar approach to build semicrystalline linear and star homopolymers has shown that increasing the complexity of the chain architecture (increasing the number of arms) can inhibit the polymer crystallization. 59 Another approach consisting of changing the ratio r eq /l b to obtain a perfect crystalline structure has been recently implemented by some of us to model semicrystalline homopolymers. 58 The optimal set of parameters leading to a stable body-centered orthorhombic structure has notably been determined. This method has demonstrated its ability to model semicrystalline polymers with realistic morphologies such as folded chains, ties, and loop molecules. 60,61 In the present work, we propose to combine this technique to the simulation of segmented copolymers to model the structuring behavior of semicrystalline TPE previously studied from the experimental point of view. 8,62,63 The article is organized as follows: After explaining the different methods and modeling aspects in section 2,w e present and discuss the results in three sections. First, results concerning the properties of the neat soft and hard homopolymers are presented in section 3.1. Then, section 3.2 details how the chain architecture influences the crystallization kinetics. Finally, the resulting morphologies are presented in section 3.3, and the relationship with the copolymer architecture is discussed. Section 4 is devoted to concluding remarks.

METHODS
2.1. Potentials. In the KG model, polymer chains are made of beads representing groups of atoms linked by nonlinear springs. The weak interactions between nonbonded beads separated by a distance r are given by a Lennard-Jones (LJ) potential: i where ε is the depth of the potential well and σ corresponds to the distance at which the potential is zero. The bonds between neighbors separated by a distance l are modeled by the combination of a finitely extensible nonlinear elastic (FENE) potential and a Weeks−Chandler−Anderson (WCA) potential (purely repulsive LJ potential): where ε′ is the depth of the potential well and σ′ the distance at which the potential is ε′.
The parameters defining these potentials are chosen in a careful way which is described in sections 2.2 and 2.3. In the following, physical quantities related to the simulations are expressed in the LJ unit system. The energy, length, and time units are given by ε u , σ u , and τ u , respectively (with m / uu u , where m u is the mass unit).
2.2. Copolymer Model. Modeling block copolymer chains requires using two types of beads, corresponding to the two types of blocks (SS and HS). This involves using two types of bonds (U bond SS and U bond HS ) and three types of weak interactions (U LJ SS , U LJ HS , and U LJ SS−HS ), defined by a total of 14 parameters. Two additional parameters have to be considered, the bead mass m b and the number of beads per segment N b , resulting in a total of 18 parameters. Because we wish the HS to be able to crystallize, the parameters for U bond HS and U LJ HS were directly set to the values found by Morthomas et al. 58 Then, the corresponding parameters for SS must be determined relative to the HS ones by taking into account the differences of chain rigidity, length, density, and interaction strength between the two types of blocks. This has been achieved in the same vein as the work of Svenaborg et al., who adapted the KG model to specific polymers by tuning the chain stiffness. 64,65 The density can be expressed as chain length per unit volume. This quantity, expressed as the bond length l b divided by the volume of a bead V b in the model, is wished to match with the ratio of the Kuhn length l k over the volume occupied by a Kuhn monomer V k : The volume occupied by a Kuhn monomer can be expressed as a function of the Kuhn molecular mass and the specific mass: Also, the volume of a bead is related to the bond length in one direction and depends on the equilibrium distance between the nonbonded bodies in the two directions perpendicular to the axis of the chain: V b ∝ l b σ 2 . Combining these three relations leads to eq 4a, which can be used to fix the density of each type of block to a target value.
The chain stiffness is usually expressed as the Flory's characteristic ratio C R Nl 2 bb 2 = ∞ (ratio between the mean-square end-to-end distance of the actual chain and the one resulting from a random walk). 66 The freely rotating chain (FRC) model, which seems quite relevant in this case since all torsion angles are equally probable (no torsion potential in our model), provides the relation between C ∞ and the bond angle In consequence, the stiffness of the chain can be tuned by modifying the bond angle θ, which is usually achieved by adding a bending potential. 48,49,64,65 However, we chose another way to tune the stiffness that does not require any additional potential and consists of changing the ratio r eq /l b . 58 Indeed, the bond angle is related to this ratio according to . Because we have already fixed σ in the LJ potential (i.e., r eq ) to control the density, we can only change l b to modify the angle. To do that, one can modify the minimum position of the bond potential by changing σ′ in the WCA part. To predict the effect of such a change on the chain stiffness, one can use the previous relation between r eq , l b , and θ to express the total potential U tot = U LJ + U bond as a function of the bond length l b (keeping the bond angle θ constant) or, in return, as a function of the bond angle θ (keeping the bond length l b constant). That way, the mean bond angle and bond length can be predicted by using the partition function: By using the expressions 3a and 3b, one can then find the optimal values of the bond potential parameters which lead to the target stiffness (expressed as Kuhn length l k ) which must be compared with the product l b C ∞ in the simulation. This is the meaning of eq 4b.
The value of C ∞ is also related to the number of beads included in a Kuhn monomer: Because the number of Kuhn monomers N k can be simply expressed from the ratio between the mass of a segment M s and the Kuhn molar mass M k , one can obtain the total number of beads in a segment N b as shown in eq 4c. Following the same idea, the mass of a bead m b can be defined from eq 4d.
Finally, the interaction strength between the beads is controlled by the depth of the LJ potential ε. It is estimated from the solubility parameter δ of the polymer. Indeed, the latter is related to the square root of the interaction energy per volume unit which is described by eq 4e. These five considerations constitute the following system of eqs 4a−4e, which need to be solved to determine the simulation parameters corresponding to specific polymer properties.
Once the parameters for HS and SS are calculated, the values for cross-interactions (HS and SS interacting together with U LJ SS−HS ) can be determined from the usual Lorentz−Berthelot combining rules expressed in eqs 5a and 5b.

Systems.
In this paper, the model is applied to the specific case of PBT−PTHF (poly(butylene terephthalate)−polytetrahydrofuran) copolymers whose target properties are provided in Tables 1 and 2. One can note that some properties of PTHF and PBT are estimated from values concerning polymers of close chemistry (respectively PEO and PET) because of the lack of data in the literature. As mentioned before, the HS parameters were fixed to the values determined by Morthomas et al. 58 Only the size of the segment N b HS was adjusted according to eq 4c (see Table 3). Then the corresponding SS parameters were determined so that the ratios between the HS and SS physical properties comply with the five scaling eqs 4a−4e. Bonds parameters k SS , R 0 SS , and ε′ SS were fixed to their classic values, leaving only one degree of freedom for varying the product l b C ∞ in eq 4b which was controlled by adjusting σ′ SS . This method results in a five equations system with five unkowns: m b SS , N b SS , ε SS , σ SS , and σ′ SS . All parameters calculations have been performed at a temperature of 4 ε u /k B , resulting in the numerical values presented in Table 4. As explained in section 2.2, the parameters values for crossinteractions ε HS−SS and σ HS−SS are deduced from eqs 5a and 5b.
Although real TPE (such as PBT/PTHF) are often made of multiblock copolymers with 7 or 8 HS per chain on average, we restricted our simulations to the case of shorter chains to limit the computing time: pentablock copolymers (containing 3 SS and 2 HS) and dumbbells (triblock copolymers with 1 SS and 2 HS). The HS length was fixed to 21 beads, corresponding to a molar mass of 970 g mol −1 . For each copolymer type, seven samples were made with different SS lengths (between 5 and 114 beads) corresponding to molecular weights ranging from 130 to 3000 g mol −1 . In the following, they are denoted SSN b SS where N b SS is replaced by the number of beads per SS. The limit case of 0 bead per SS is also treated as a neat HS of 42 beads (denoted neat HS42), namely twice the length of a HS. Neat hard and soft polymer chains of 100 beads are also created and denoted respectively neat HS100 and neat SS100. The number of chains in each box is fixed to 500. The different chain architectures are detailed in Figure 1 and Table 5.

Simulation
Procedure. All the molecular dynamics simulations are performed using the classic LAMMPS code. 71 The copolymer chains are built using a radical-like polymerization algorithm (RLP) previously developed by Perez and co-workers. 72,73 The simulation is performed in the NPT ensemble, using periodic boundary conditions, at P = 0.5 ε u /σ u 3 and T =4ε u /k B . It starts from a bath containing two types of loner beads (SS and HS types). 500 beads are randomly chosen to become radicals which will try to find a loner in their vicinity to bind to. The radical is then transferred to the new extremity of the molecule while the old one becomes an inner monomer. This procedure is repeated until the chain reached its target size, building the different blocks in a stepwise manner. A short delay of 5 τ u is introduced between each growth step so that the chain can relax all through its growth. At the end of this process, the remaining loner beads are removed and the system is equilibrated during 2 × 10 5 τ u .
Once the copolymer melt is generated and equilibrated, cooling ramps are performed until T = 0.01 ε u /k B at two different rates: 5 × 10 −5 τ u −1 and 5 × 10 −6 τ u −1 for the study of glass transition and nonisothermal crystallization, respectively.
To improve the statistics of the simulations, the whole equilibration and cooling ramp has been repeated on one system (pentablock SS76) replicated eight times in space (two times in each direction). This new system has been equilibrated during 5 × 10 5 τ u at T =4ε u / k B to obtain different (new) initial configurations in the whole box before starting the cooling ramp. A comparison with the results obtained on the small system is provided in Supporting Information section 1, showing the good reliability of these results.
2.5. Post-Treatment Analysis. Post-treatment analysis is performed right after equilibration or during the cooling ramp. To look at the chain conformations in the molten state, the mean-square internal distance (MSID) ⟨R 2 (n)⟩ is computed at T =4ε u /k B .I n addition, an entanglement analysis is performed using the Z1 code of Kroger et al. which uses geometrical considerations to extract the primitive path of a dense polymer melt. 74−77 This algorithm is used to compute the molecular weight between entanglements M e .
The study of crystallization during the cooling ramp is performed using an algorithm recently developed by Perez and co-workers and based on hierarchical clustering. 58 It enables us to detect all the crystallites in the simulation box and estimate their size. Previous works using a similar coarse-grained model have shown that crystal clusters made of fewer than 50 beads are not stable over time. 60 Therefore, a threshold has been fixed so that every organized cluster below this size is considered as amorphous. The analysis of the amorphous segments is also part of the algorithm that provides a global insight into the chains conformation and network topology. It classifies the amorphous segments into four categories: "tie" if it       bridges two different crystallites, "loop" if it is connected twice to the same crystallite, "dangling end" if it has only one extremity connected to a crystal, and "purely amorphous" if the whole chain has no crystallized segment. This analysis is performed every 0.1 ε u /k B during the slow cooling process. The MSID of the different types of segments is also computed during this ramp. Finally, the morphology of the multiphasic material is studied by using the radial distribution function (RDF) computed with the Ovito software. 78 To analyze the phase separation process, the RDF of HS beads is computed during the cooling ramp (having removed the SS beads from the box to keep the hard phase only). The resulting morphology after crystallization of the hard segments is also analyzed by computing the RDF on the crystallized HS beads only, at a fixed temperature T =1ε u /k B .

RESULTS AND DISCUSSION
3.1. Impact of Chain Rigidity. The main difference between the two types of blocks constituting the copolymers is their chain stiffness. Considering the potentials parameters values, a FRC model predicts a C ∞ value of 6.19 for the HS and 4.47 for the SS at the reference temperature T =4ε u /k B . Figure 2a shows the normalized MSID C(n)=⟨R 2 (n)⟩/(nl b 2 ) which is expected to tend to these values for high n. Indeed, the graph shows C(100) ≈ 5.9 for the HS and C(100) ≈ 4.2 for the SS. Moreover, one can see that the full FRC prediction given by eq 6 is in good agreement with the effective C(n) curve. This indicates that the Kuhn lengths of our modeled systems are well controlled.
Another conformational property resulting from the chain rigidity and density is the tendency to form entanglements, quantified by the molecular weight between entanglements M e = m b N e = m b C ∞ N Ke , where N e and N Ke are respectively the number of beads and number of Kuhn segments separating two consecutive entanglements. The primitive path analysis of the two neat polymer chains (N b = 100) at T =4ε u /k B gives N e SS = 44.7 (N Ke SS = 10) and N e HS = 23.2 (N Ke SS = 3.7). To compare with the literature, these values can be expressed in terms of molecular weights (eq 4c), resulting in M e SS = 1180 g mol −1 and M e HS = 1070 g mol −1 , in good agreement with the available data: 1080−1420 g mol −1 for PTHF 79 and 1160 g mol −1 for PBT. 80 One can note that the numbers of Kuhn segments per entanglement strand are also qualitatively consistent with the data concerning similar chemistries, PEO and PET, for which N Ke equals respectively 14.7 and 5.5. 68 Such a result is of primary importance for a good modeling of the mechanical properties, in particular the plateau modulus.
Two critical properties for TPE are the glass transition temperature T g and the crystallization (or melting) temperature T c (or T m ) of each block. Typical values of these temperatures follow T g SS < T m SS < T room < T g HS < T m HS to have a rubbery mechanical behavior at room temperature, T room .I n the case of PBT−PTHF systems, T g SS = −60°C, T m SS =0°C, T g HS =6 0°C, and T m HS = 220°C. 8,63,81 The model is expected to reproduce (at least in a qualitative way) this thermal behavior with the four phase transitions occurring at different temperatures. Indeed, the glass transition temperature is expected to be related to the chain stiffness C ∞ and the interaction energy ε, which have been both carefully adjusted for each block. Empirical models predict a linear dependence of the glass transition temperature with the cohesive energy per flexible unit, E coh . 82 Considering a Kuhn segment containing C ∞ beads whose interaction energy is ε consequently leads to E coh ∝ C ∞ ε. Given the present parameters values, this quantity should be 3 times larger for HS compared to SS. To check the resulting difference of T g between the two modeled polymer chains, one needs to follow the evolution of the volume per bead during a quick cooling ramp (see section 2.4 for details). One can note that the glass transition temperature is expected to be rate-dependent (which is not studied here) because of its second-order nature, contrary to the melting/crystallization process. The results for both neat polymer chains are shown in Figure 2b. The glass transition temperature is extracted from the break in slope position, giving T g SS = 0.46 ε u /k B and T g HS = 1.5 ε u /k B . This seems to be in good agreement with the previous estimation of cohesive energy in HS being 3 times higher than the SS one. The crystallization temperature is also expected to be affected by both the chain mobility (which depends on the chain stiffness C ∞ ) and the driving force related to the interaction energy ε. The crystallization can be observed in a slower cooling process (see section 2.4 for details), and T c can be taken as the temperature at which the time derivative of enthalpy is maximum. Figure 2c gives T c SS = 0.9 ε u /k B and T c HS = 2.0 ε u /k B . These values are in qualitative agreement with our previous studies in which we reported that the isothermal crystallization of neat HS can start at ca. 2.3 ε u / k B . 58,60,61 As expected, the model leads to four distinct temperatures which qualitatively match the reality. A simple mapping of these results with the real temperatures cited above indicates that the room temperatureusual operating temperature for TPEshould be around 1 ε u /k B .

Copolymer Architecture Drives Crystallization.
In the previous part, we showed how the chain rigidity controls the conformation of the neat polymer chains. Thus, changing the copolymer morphology should impact the chain conformation. Figure 3 provides the MSID of both copolymer types evidencing that modifying the chain architectureby changing either the SS length or the number of segments enables to tune the chain stiffness between the two extreme cases corresponding to neat HS and neat SS. Indeed, the global rigidity of the chain increases when decreasing the SS length, i.e., when increasing the HS fraction. However, it is clear that the MSID of the copolymer does not follow a classic FRC model whose bond length and bond angle were taken as the averaged values of each neat polymer. In fact, the MSID curves are nonmonotonic, which is not related to any lack of equilibration but rather to the alternation of two rigidities in the copolymer chain (see the convergence of MSID over equilibration in Supporting Information section 2). In particular, one can notice the pronounced higher rigidity at the extremities of the dumbbellsending with HScompared to pentablocks whose ends are made of SS.
As explained in section 2.4, the nonisothermal crystallization of the copolymers is studied during a cooling ramp starting from the melt and ending in the SS glassy state. The rate (5 × 10 −6 τ u −1 ) has been adjusted so that the crystallization can easily occur within a reasonable time. The phase transitions can be observed by looking at the evolution with temperature of the time derivative of enthalpy. The so-called thermograms of all copolymers are plotted in Figure 4. All the curves exhibit a peak between T = 1.5 and 2 ε u /k B corresponding to the HS crystallization. Its position can be taken as the crystallization temperature, T c HS . One can note the decreasing amplitude of the peak when increasing the SS length. This is due to the fact that the area under the peak is proportional to the degree of crystallinity of the whole system. Therefore, lowering the HS fraction obviously leads to smaller degree of crystallinity. A second peak corresponding to the SS crystallization can be observed at lower temperatures for the systems containing the highest fraction of SS (see the second peak at 1.2 ε u /k B for pentablock SS114 on inset of Figure 4a).
These results come together with the crystal analysis described in section 2.5. The algorithm quantifies the number of crystallized beads, distinguishing HS and SS, from which we can extract the relative fractions of crystallized HS and SS: φ HS and φ SS . Their evolution with temperature is represented in Figure 5 for both pentablock copolymers and dumbbells. The crystallinity remains zero down to a temperature close to T c HS determined on the thermograms. Then φ HS starts to increase up to a pseudoplateau ranging between 70% and 90% for 2 )a tT =4ε u /k B for pentablock copolymers (a) and dumbbells (b). Lines stand for the FRC predictions defined in eq 6 for neat soft and hard polymer chains. temperatures below 1 ε u /k B . One can note the good agreement between this value and experimental data obtained for multiblock copolymers made of short and monodisperse HS. 62 The increase of φ SS is much slower and leads to lower relative crystalline fractions (between 10% and 30% for temperatures below 0.5 ε u /k B ). This will be commented further.
An interesting result is the evolution of the HS crystallization temperature when changing the chain architecture. Figure 4a shows that increasing the SS length in pentablock copolymers tends to lower the value of T c HS .I n fact, the crystallization temperature ranges from ca. 1.6 ε u /k B for SS114 up to 1.9 ε u /k B for SS5. It is even higher for neat HS42 (T c HS = 2.0 ε u /k B ). This is also clear in Figure 5a when looking at the onset of the φ HS upturn which is postponed to lower temperatures when increasing the SS length. On the contrary, the crystallization kinetics of the pentablock copolymer tends to the one of neat HS when decreasing the SS length. These results are consistent with previous experimental studies on PBT−PTHF segmented block copolymers revealing that sufficiently short SS resulted in similar kinetics as for neat PBT. 63 However, one can notice that this dependence on SS length is not present for dumbbells. T c HS ranges between 1.9 and 2 ε u /k B whatever the SS length (see Figure 4b), and the onset of the crystallinity increase is always at 2 ε u /k B (see Figure 5b). In other words, the HS crystallization is only impacted in the multiblock systems, not in the dumbbells (given the fixed HS length). Figure 5 also shows the impact of the copolymer architecture on the SS crystallization kinetics. While neat SS crystallize at T c SS = 0.9 ε u /k B , SS inside the copolymers start crystallizing at temperatures much closer to T c HS , around 1.7 ε u /k B . In fact, it seems that the SS crystallization is initiated by the HS one. This is not surprising since the lattice parameters of both crystalline phases are very close. This is a consequence of the very similar LJ σ parameters for both SS and HS interactions (see Tables 3 and 4) and the way we choose the value for cross-interactions (see Lorentz rule in eq 5a), promoting the SS crystallization from the HS crystalline walls. This is quite clear when looking at the snapshots of Figure 10 showing the crystalline structures only. There is no isolated SS crystallite: all of them are attached to their HS counterparts. However, we know from experimental studies that the real crystalline structures of PBT and PTHF have different lattice parameters preventing this promoted SS crystallization. 83,84 This artifact, generated by the choice of the potentials, should not affect the following results because of the relatively low crystalline SS fractions.
Focusing on the understanding of the HS crystallization kinetics, we propose below to follow the evolution of the systems during the cooling ramp above the crystallization point, namely between 3 and 2 ε u /k B . We know from scattering 8,62 and rheology 85 experiments performed on such segmented copolymers that a remaining phase separation is expected in the molten state, likely to play a major role on the crystallization kinetics and resulting network topology. The ordering process while cooling can be monitored by following the RDF of the HS beads. Its evolution with temperature is represented in Figure 6 for two couples of comparable systems having the same HS fraction and same chain length, but different numbers of blocks: (a) pentablock SS19 vs dumbbell SS57 (x HS = 42%) and (b) pentablock SS38 vs dumbbell SS114 (x HS = 27%). At high temperature (T =3ε u /k B ), the RDF curve is almost completely flat in the high-r range (for distances greater than the cutoff radius 2.5σ), suggesting there is no structuring at long distances. When decreasing the temperature, the shape of the RDF changes. A correlation peak progressively grows, which is the sign of structuring at this length scale. The HS start to form domains, separated from the SS by a distance r ≈ 30 σ u . It can be noted that this happens at temperatures higher than T c HS , meaning that the phase separation occurs before crystallization. Comparing then the pentablocks and the dumbbells clearly evidences that the phase separation occurs at higher temperatures in the latter. Therefore, an increase of the number of blocks seems to significantly delay the phase separation process, in good agreement with experimental studies comparing in a systematic way the melt structure in multiblock and triblock copolymers. 86 We propose to rationalize these differences between triblocks and pentablocks using the Krause theory, 22 which provides an expression of the free energy change on microphase separation of multiblocks copolymers made of A and B units using the lattice theory of polymers. Unlike classic approaches for diblocks, this model includes an additional entropy decrease caused by the immobilization of the A and B units which are linked at the interface between two blocks. The Krause prediction for free energy change on microphase separation is then   (7) where N r is the total number of sites, ν A and ν B are the volumes fractions of A and B units, respectively, χ AB is the Flory interaction parameter (which is a function of temperature, see Supporting Information section 3), z is the coordination number, N 0 is the number of chains, m is the number of blocks per chain, and ΔS dis is the entropy loss when one segment of a chain is immobilized. This model can be used to determine the MST temperature defined as the temperature at which the free energy change is zero. This temperature depends on the chain length (N r /N 0 ), fractions of A and B units (ν A and ν B ), interaction strengths (χ AB ), and the number of blocks m. In our case, the comparison between pentablock SS19 and dumbbell SS57 or pentablock SS38 and dumbbell SS57involves only the change of the parameter m. The application of this model to our systems is presented in Figure 7 where we can see the evolution of the free energy change on microphase separation with temperature for the four systems listed above (see Supporting Information section 3 for details concerning the application of the Krause model to our systems). Although we do not claim to determine the absolute values of the MST temperatures from this graph, we note the consistent order of magnitude for the four values: 1.5, 2.1, 2.7, and 3.8 ε u /k B for pentablocks SS19 and SS38, dumbbells SS57 and SS114, respectively. Moreover, the differences in terms of MST temperatures qualitatively match with the phase separation kinetics observed on the RDF ( Figure 6). Basically, increasing the number of blocks makes the phase separation more difficult. In other words, the phase separation occurs at lower temperatures when increasing the number of blocks.
Once shown the differences in terms of phase separation between dumbbells and pentablocks, we can now try to explain the different crystallization kinetics observed for both types of system. We have previously observed that all the dumbbells crystallize at T c HS =2ε u /k B , which is also the crystallization temperature of the neat HS system. On the contrary, the pentablocks have a lower crystallization temperature, which further decreases when increasing the SS length. Looking at the RDF and using the simple thermodynamic model of Krause, we have found that dumbbells phase separate at higher temperatures compared to pentablocks. In fact, we have demonstrated the dumbbells are already well phase-separated at T =2ε u /k B , whereas it is not the case for pentablocks. With this picture in mind, one may assume that the crystallization temperature is driven by the phase separation. In the case of pentablocks, the system is not enough phase-separated to crystallize when reaching the crystallization temperature of the neat HS (T c HS =2ε u /k B ). It requires a lower temperature, favoring the HS gathering, to subsequently crystallize. On the contrary, since dumbbells are already well phase-separated when reaching the crystallization temperature of the neat HS, they crystallize immediately.
However, Krause theory predicts that increasing the SS length from 19 beads to 38 beads will increase the MST temperature (see Figure 7), while the simulation exhibits a decrease on the crystallization temperature (see Figure 5a). In fact, increasing the SS length decreases the HS fraction in the system and, thus, the probability for a HS to find a pair. Therefore, the system requires a high level of phase separation to crystallize. On the contrary, a system with very short SS does not require such a high degree of phase separation to crystallize because of the high HS fraction. Even though the phase separation is difficult in such a systemlow MST temperature because of the low SS fractionit can easily crystallize due to the high probability for a HS to find another  . Normalized free energy of phase separation predicted by the Krause model as a function of temperature. Comparisons for two couples of systems at iso-HS fractions (x HS = 42% for pentablock SS19 and dumbbell SS57, x HS = 27% for pentablock SS38 and dumbbell SS114). The dotted line shows the limit between the onephase (ΔG > 0) and phase-separated (ΔG < 0) states. one even at high temperatures. In this way, it appears that the phase separation is not the critical parameter for the HS crystallization which is mainly driven by the probability for a HS to find a pair in its vicinity. Therefore, one can expect the crystallization to be controlled by the local HS number density: 8) where ρ n HS is the global HS density defined as the ratio between the number of HS beads and the volume of the box, g(r) is the RDF of HS beads, and L is a cutoff radius defining the vicinity of a HS bead. Equation 8 has been used to obtain the evolution of the local HS number density before and at the crystallization onset for all the systems. A cutoff radius L = 2.5 σ u was used to take into account the nearest neighbors only (see RDF of HS beads in Supporting Information section 4 for the choice of L). In the case of neat HS, the number density being about 0.24 σ u −3 , this cutoff radius encompasses πρ −= neighboring beads. Figure 8 represents the evolution of the relative HS crystalline fraction with this local HS number density on a temperature range going from 3 to 1.5 ε u /k B . While the global HS densities of all these systems are widely distributed (ranging from 0.03 σ u −3 for pentablock SS114 to 0.23 σ u −3 for dumbbell SS5), the onset of crystallization always occurs at the same local HS density: 0.14 ± 0.02 σ u −3 . The fact that the crystallization onset occurs at an almost constant local HS density independently of the chain architecture shows that crystallization is mostly driven by the probability for a HS to find a pair.
3.3. Resulting Semicrystalline Morphologies. Apart from kinetics, the copolymer architecture is also found to impact the final crystallinity, as evidenced in Figure 5. Modifying the SS length appears to affect the HS crystallinity reached at low temperature: systems with longer SS reach a higher relative crystallinity of both hard and soft segments. φ HS goes from 70% for neat HS polymer chain to more than 90% for copolymers with the longest SS. The increase is even slightly more pronounced for dumbbell systems, going from 65% to more than 95% (probably due to the higher mobility of the HS located at the chain ends). This result is not intuitive. In fact, one could expect that decreasing the HS content would lead to smaller relative HS crystalline fractions due to the reduced probability for a HS to find a pair. However, looking at the topology of the amorphous phase can help to understand this result. As explained in section 2.5, amorphous segments surrounding the crystals can be classified in three categories: dangling ends, loops, and ties. Figure 9 represents the numbers of each type of amorphous segments per crystallized segment for all the systems at T =1ε u /k B . It should be read as the tendency for amorphous segments to preferentially form loops (N L ), ties (N T ), or dangling ends (N DE ). Interestingly, while the number of ties per crystallized segment remains almost independent of the chain architecture (N T ranging between 0.3 and 0.4), the number of loops per crystallized segment significantly increases when increasing the SS length (by a factor 2 when comparing the extreme SS sizes). Then, the higher relative HS crystallinty observed when increasing the SS length can be explained by a higher ability to form loops. Furthermore, we can see their number increasing more quickly for N b SS > 38 beads. This size corresponds to 1 kg mol −1 , which remarkably matches the experimental literature regarding the minimal size for a PTHF segment to fold, being between 1 and 2 kg mol −1 8,33,87 (see Supporting Information section 5 for additional analysis of the MSID in SS loops). This ability to fold has also an impact on the relative SS crystallinity φ SS . One can see in Figure 5 that the final value of φ SS reached at low temperature is almost constant for the systems having a SS length smaller than 38 beads (φ SS = 0.22 for pentablocks and 0.11 for dumbbells). However, when N b SS exceeds this value, φ SS increases quite significantly, by a factor 1.5 for pentablock  SS114 and a factor 2 for dumbbell SS114. This effect is once again consistent with the experimental literature showing that PTHF blocks crystallize when their molecular weight exceeds 1 kg mol −1 . 8, 63 The crystal's arrangement in the box and its dependence on the chain architecture is also interesting to analyze because it can be easily compared with experiments such as scattering or microscopy. 8,62 For that purpose, the RDF of the crystalline HS beads have been computed at T =1ε u /k B for all systems and are represented in Figure 10. The curves present a maximum corresponding to the correlation length d* between aggregated domains corresponding here to the HS crystal-litesone can note that d* is much smaller than the box length (last point of each curve in Figure 10), meaning that finite size effects are negligible here (see also Supporting Information section 6 for more information about finite size effects). This distance (d* ≈ 30 σ u ) is very close to the correlation length observed at higher temperature for systems being phase-separated but still amorphous (see Figure 6). This corroborates the idea of a crystallization occurring in a preordered melt whose morphology is preserved after crystallization (see section 3.2). Moreover, one can note the slight but significant change on the correlation length when varying the SS length. As one could have expected, decreasing the SS length leads to closer crystallites, i.e., a smaller value of d*. A comparison with experiments can be done on the pentablock SS76 system, which has the same block sizes (but fewer blocks) as a PBT−PTHF segmented copolymer studied with X-ray scattering and atomic force microscopy (AFM) at room temperature in a previous study. 8 Converting the correlation length d* =3 4σ u from LJ units to standard unitsusing eq 4bgives a distance between the crystals about 8 nm. Experimentally, it was found that such a system elaborated through hot-pressing crystallizes forming long ribbons spaced by a distance being approximately 38 nm. These ribbons, not perfectly regular, are actually formed of smaller crystals, forming a necklace structure with crystalline nodules spaced by a distance of 10 nm. Considering the relatively small size of the simulation box, it is reasonable to guess that the larger length scale (distance between ribbons) is not accessible with the simulation. However, it seems that the smaller length scale (distance between the nodules) matches quite well with the experiments, reinforcing the consistency of the coarse-grained model presented here (see Supporting Information section 7 for a comparison with AFM imaging). In addition, the shape of the crystals, almost isotropic for pentablocks (see snapshot in Figure 10), is in good agreement with the nodules observed in experiments.
Interestingly, dumbbells are found to form a different structure, with much longer crystallites having one preferential direction of growth. This can be verified looking at the snapshots of Figure 10 representing the crystals of a pentablock and a dumbbell system having the same HS fraction. Quantitatively, we also show in Supporting Information section 8 that pentablocks tend to form further, albeit smaller, crystallites than triblocks. Experimentally, this change from a necklace structure to regular ribbons has been observed in the same study of a unique PBT−PTHF copolymer elaborated through two different processing routes. 8 Hot-pressing was found to lead to a necklace structure while solvent-casting provided continuous (and closer) ribbons. This difference was attributed to a much higher mobility of the chains in the solvent compared to the melt. It is possible to explain the difference observed here between the pentablocks (forming nodules) and the dumbbells (forming long ribbons) by using the same argument. The mobility of the HS in dumbbells is increased for two reasons: first, the shorter diffusion time of the whole chain because of the smaller chain length compared to pentablocks and, second, the higher mobility of blocks located on chain ends.
Such information concerning the semicrystalline structures in TPE and chain topologies is of prime importance to understand their mechanical properties. A recent work compiling several experimental data about the structure and dynamical mechanical properties in TPE has shown that the rubbery modulus is driven by two key parameters: the volume fraction and width of the crystallites, 88 as predicted by an analytical model developed by Baeza 33 based on the chains topology. The results presented in this work could be a useful starting point for further investigations on the mechanical properties of TPE.

CONCLUSIONS
By adapting a generic coarse-grained model originally designed for semicrystalline homopolymers to the specific case of PBT− PTHF TPE, we have implemented a molecular dynamics simulation of semicrystalline segmented block copolymers for the first time to our knowledge. We have demonstrated the key role of the chain architecture in the crystallization process, corroborating the insights from experimental studies. On the Figure 10. Radial distribution function of the crystallized HS beads at T =1ε u /k B in pentablock copolymers (a) and dumbbells (b). The arrows indicate the position of the maximum (long period d*) corresponding to the mean distance between crystals. Snapshots show the crystalline phases of pentablock SS38 (a) and dumbbell SS114 (b). Crystalline HS beads are represented in blue and crystalline SS beads are in orange.
basis of the fact that the crystallization process is controlled by the probability for a HS to find a pair, we have studied the impact of both the number of blocks and the SS length on the crystallization kinetics. Consistently with thermodynamics predictions, reducing the number of blocks has been found to promote the microphase separation, increasing this probability and thus enhancing the crystallization. When the phase separation is weaker because of a higher number of blocks in the chain, we have evidenced the decrease of the crystallization temperature when increasing the SS length. Consideration of the HS mobility has enabled explaining the different crystalline arrangements observed in pentablocks and dumbbells. The resulting morphologies obtained after crystallization have revealed a good agreement with structural experimental characterizations, reinforcing the relevance of the model. Finally, we believe the development of this coarsegrained model can serve as a basis for further investigations of the static and dynamical properties of TPE.