Nucleation theory of polymer crystallization with conformation entropy

Polymer xxx (xxxx) xxx

Contents lists available at ScienceDirect

Polymer journal homepage: http://www.elsevier.com/locate/polymer

Nucleation theory of polymer crystallization with conformation entropy Hiroshi Yokota a, *, Toshihiro Kawakatsu b a b

RIKEN Interdisciplinary Theoretical and Mathematical Sciences (iTHEMS), Wako, Saitama, 351-0198, Japan Department of Physics, Graduate School of Science, Tohoku University, Aramaki, Aoba, Sendai, 980-8578, Japan

A R T I C L E I N F O

A B S T R A C T

Keywords: Polymer crystallization Nucleation theory

Based on classical nucleation theory, we propose a couple of theoretical models for the nucleation of polymer crystallization, i.e. one for a single chain system (Model S) and the other for a multi-chain system (Model M). In these models, we assume that the nucleus is composed of tails, loops and a cylindrical ordered region, and we evaluate the conformation entropy explicitly by introducing a transfer matrix. Using these two models, we evaluate the occurrence probability of critical nucleus as a function of the polymer chain stiffness. We found that the critical nucleus in Model M is easier to occur than in Model S because, for semi-flexible chains, the nucleus in Model M can grow by adding a new polymer chain into the nucleus rather than to diminish the loop and tail parts as in the case of Model S.

1. Introduction For many years, a lot of researchers were attracted by the isothermal crystallization process of polymers, which is composed of multi-step ordering processes [1–3]. Especially, early stage of the crystallization has been extensively investigated both experimentally [2–7] and theo­ retically [8–11]. In the early stage, some scenarios of ordering mechanisms have been proposed based on experimental results, for example, spinodal-assisted crystallization [4,12,13], nucleation and growth [3,5] and appearance of a mesomorphic phase [1,14]. The mesomorphic phase, which is an intermediate state between the liquid and the solid phases, has the same symmetry as the liquid crystal. For example, in the case of poly (butylene-2, 6-naphthalate) (PBN), the structure of the mesomorphic phase has a smectic periodicity which is composed of a one dimensional periodic order along a certain direction and a liquid-like structure in the plane perpendicular to the direction [14]. In the present research, we focus on the time regime after induction period, i.e., the time regime where a critical nucleus appears. Here, the critical nucleus is defined as the minimum-size nucleus which can grow stably. Recently, the existence of the nuclei in polymer crystallization was directly confirmed by small angle and wide angle X-ray scattering (SAXS and WAXS) experiments, Fourier transform infrared spectroscopy (FTIR) and so on [2,3]. Moreover, these experimental results were supported by recent particle simulations using molecular dynamics and Monte Carlo

techniques [10,15]. Many researchers have tried to uncover the mechanism of such nucleation process by proposing simple but essential models for the polymer crystallization [15]. One of the most primitive theoretical models is the classical nucleation theory (CNT) [16], where a competi­ tion between the bulk free energy difference and the surface excess free energy determines whether the nucleus grows or shrinks. In the CNT, by calculating the free energy difference before and after a nucleus appears, we obtain the size and the occurrence probability of the critical nucleus. A nucleus shrinks if it is smaller than the critical size, while a nucleus larger than the critical size grows. It is noted that the CNT accounts for both of homogeneous and heterogeneous nucleation processes where the homogeneous nucleation is driven by thermal fluctuation in a uni­ form system while the heterogeneous nucleation is initiated by the im­ purity such as a solid wall. In 1960, Lauritzen and Hoffman applied the CNT for both of homogeneous and heterogeneous nucleations in poly­ mer crystallization, where they assumed an anisotropic shape such as a cylindrical shape for the nucleus [8]. In the homogeneous nucleation process, the free energy difference before and after the nucleation, Δf, of an isolated nucleus is given by Δf ðr; lÞ ¼

π r2 lΔμ þ 2πr2 σ t þ 2πrlσ s ;

(1)

where r and l are the radius of the top surface and the height of the cylindrical shaped nucleus. Moreover, Δμ; σ s and σ t are the bulk energy difference per volume between liquid and solid, the surface tension of

* Correponding author. E-mail address: [email protected] (H. Yokota). https://doi.org/10.1016/j.polymer.2019.121975 Received 4 July 2019; Received in revised form 8 October 2019; Accepted 4 November 2019 Available online 9 November 2019 0032-3861/© 2019 Elsevier Ltd. All rights reserved.

Please cite this article as: Hiroshi Yokota, Toshihiro Kawakatsu, Polymer, https://doi.org/10.1016/j.polymer.2019.121975

H. Yokota and T. Kawakatsu

Polymer xxx (xxxx) xxx

but also of tails and loops (see Fig. 1 (a)), where the tails and the loops are described by bead spring model. On the other hand, the effect of the conformation of single chain system is different from the one of the multi-chain system whose crystallization has been one of the main tar­ gets in the experiments and the simulations. In addition to the number of chains participating in the nucleus, the chemical detail of the polymer chains is another important factor for the nucleation. Although results of the experiments and of the simulations depend on the chemical details [5–7,10], Muthukumar’s model cannot reflect such a dependence because Muthukumar assumed that the polymers are simple flexible bead-spring chains. In this paper, we construct nucleation theory which is applicable to both of single chain and multi-chain systems. Our model includes a chain stiffness, which reflects an important chemical detail of the polymer chain, through the calculation of the conformation entropy, where the stiffness is described by the persistence length or by the en­ ergy difference between trans and gauche conformations. Here, the chain stiffness affects the important chemical details of the polymer chains in the samples. This paper is organized as follows. In the next section, we discuss our theoretical models of nucleation in single-chain and multi-chain sys­ tems. First, we construct a theoretical model for single-chain systems. Next, a model for multi-chain systems is constructed based on the singlechain model introduced above. In Sec. 3, we compare these two models by discussing the effect of the stiffness of the polymer chain on the size of the critical nucleus. We further discuss the relationship between the induction time and the stiffness of the polymer chain. Finally, we conclude our results in Sec. 4.

Fig. 1. (a)A schematic picture of Model S of the nucleus. The tails and loops are constructed by a semi-flexible chain. (b)If the ordered region is negligibly thin, the conformation entropy of the tails and of the loops are approximated as the one of the flower micelle because the number of micro states of the ordered region is just unity.

the side surface and, that of the top/bottom surfaces. We note that σ t implicitly includes the effect of the conformation entropy of the chains. Let us consider a free energy landscape defined by a curved surface Δf ðr; lÞ on the ðr; lÞ plane. We can easily confirm that this free energy landscape has a saddle point at ðr� ; l� Þ, which specifies the size of the critical nucleus [8]. By using the value of Δfðr� ;l� Þ, we can evaluate the activation energy of the critical nucleus which determines its occurrence probability. Here, it should be noted that, among the two types of studies performed by Lauritzen and Hoffman, only the CNT for the heteroge­ neous nucleation of polymer crystallization is called Lauritzen-Hoffman (LH) theory [8,17]. In LH theory, the nucleation is assumed to start from a solid wall, where the chains in the nucleus are aligned parallel to the solid wall. In this case, the area of the side surface of the nucleus is less than the one in the case of the homogeneous nucleation, which leads to a change in the contribution from the surface tension of the side surface in eqn. (1). LH theory was successfully used in the analysis of the growth rate of the spherulite, where the values of the surface tensions are regarded as fitting parameters. Despite of the importance of the critical size of the nucleus in the theoretical studies, it is difficult to measure it experimentally due to the limitation of the spatial resolution of experiments. Then, computer simulation is a powerful tool to investigate the characteristics of the critical nucleus. In fact, many computer simulations unveiled the dy­ namics of the nucleation and the structural characteristics of the nucleus [10,11]. In these studies, however, the sizes of the critical nuclei have ambiguity because the crystalline region is not well-defined [15,18]. The crystalline region is usually defined by introducing threshold values of orientational order parameters which describe the interchain and intrachain bond correlations [11]. Here, a problem comes from the fact that the threshold values in many simulation studies are different from each other. Although Welch reported a technique to evaluate these threshold values autonomously by using machine learning, there still remains a problem on the precision of the calculation [11]. Some simulation studies reported that chain conformation in the nucleus affects the condition of the crystallization [15,19]. The nucle­ ation from a large solid wall of a bulk crystal leads to the appearance of a hairpin-like structure of a chain [19]. On the other hand, such a hairpin-like structure does not appear in the simulation study on the homogeneous nucleation [15]. Although Anwar et al. concluded that the difference between these behavior is due to the condition of the nucle­ ation [15], the detail of the physical mechanism has not been clarified. To explain these results of the simulations, we need to incorporate the effect of the conformation entropy into the analytical model such as CNT. Such an attempt has first been done by Muthukumar in 2003, where an extension of CNT was proposed for crystallization of a single polymer chain by adding the effect of conformation entropy [20]. In this model, a nucleus is assumed to be composed not only of straight parts

2. Model 2.1. Nucleation theory for single polymer chain (model S) In this subsection, we construct a theoretical model of nucleation of a single polymer chain (Model S) which serves as a basic component of the multi-chain model discussed later. We calculate the free energy differ­ ence before and after the nucleation takes place, where we assume that the nucleus is formed by a single semi-flexible polymer chain. We also assume that the shape of the nucleus is as shown in Fig. 1 (a). This nu­ cleus is composed of straight parts in the ordered region, and tails and loops outside of the ordered region. We assume that the two ends of each loop are on the same side of the ordered region. The free energy difference Δf is evaluated as a function of the number of loops α and the height of the ordered region m � b where b is the size of a segment of the polymer chain (l ¼ m � b, where l appeared in eqn. (1)). Then, we obtain pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Δf ðα; mÞ ¼ ðα þ 1Þmb3 Δμ þ 2 πðα þ 1Þmb2 σs þ 2ðα þ 1Þb2 σ t � � �� (2) Δε kB T ln Z α; m; ; N þ 1; M ¼ 1 : kB T In eqn. (2), Δμ is the bulk free energy difference per unit volume, σs the surface tension of the side surface, σt the surface tension of the top/ bottom surfaces, Δε the energy difference between trans and gauche conformations (the chain stiffness), kB Boltzman constant, T the actual temperature, M the number of chains participating in the nucleus, N þ 1 the total number of segments in a polymer chain and Zðα; m; Δε =ðkB TÞ; N þ1; M ¼ 1Þ the partition function of conformations of the tails and the loops of the single chain. In deriving eqn. (2), we assumed that the segments are closely packed in the ordered region, which leads to the condition on the total volume of the ordered region as

πmbr2 ¼ ðα þ 1Þmb3 :

(3)

The reference value of the free energy is that of a single chain in the free space, which corresponds to the free energy in the liquid state composed of theta solvents and of a single chain with N þ 1 mðα þ1Þ 2

H. Yokota and T. Kawakatsu

Polymer xxx (xxxx) xxx

segments. The bulk energy difference Δμ � b3 in eqn. (2) is given by Δμb3 ’

T ð0Þ T m Δh; T ð0Þ m

(4)

where T m , T and Δh are the equilibrium melting temperature, the actual ð0Þ

temperature of the system (thus, T m T is the degree of supercooling) and the heat of fusion per segment, respectively [21]. We ð0Þ

non-dimensionalize eqn. (2) by using kB T m as ð0Þ

Δf ðα; mÞ ¼ kB T ð0Þ m

ðα þ 1Þm

T ð0Þ T Δh m T ð0Þ kB T ð0Þ m m

pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi b2 σs b2 σ t þ2 π ðα þ 1Þm þ 2ðα þ 1Þ ð0Þ kB T m kB T ð0Þ m � � �� ð0Þ T Δε T m ln Z α; m; ; N þ 1; M ¼ 1 : ð0Þ ð0Þ T Tm kB T m

(5)

Fig. 2. Our model of nucleus (Model M) composed of multi-chains. In this example, the number of chains participating in the nucleus M is 3. Moreover, the number of loops in the nucleus α is 3.

Hereafter, we use the non-dimensional parameters where the energy, ð0Þ temperature and length are normalized by kB Tð0Þ m , T m and b, respec­ ð0Þ tively. For example, we simply refer to Δh=ðkB T m Þ as Δh. Thus, eqn. (5) is rewritten in the following manner pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Δf ðα; mÞ ¼ ðα þ 1Þmð1 TÞΔh þ 2 πðα þ 1Þmσ s þ 2ðα þ 1Þσ t �i h � (6) Δε T ln Z α; m; ; N þ 1; M ¼ 1 : T

� T~ ηξ ðqÞ ¼ exp

� i q ⋅ eðηÞ T 2

ηξ

� exp

� i q ⋅ eðξÞ : 2

(11)

here r and q are the position vector and the wave number vector. As shown in Appendix A, the expression of the statistical weight of conformation of the flower micelle Zðα; m; Δε =T; N þ1; M ¼ 1Þ is written as � � Δε Z α; m; ; N þ 1; M ¼ 1 T

For simplicity, we assume that the ordered region is negligibly thin, therefore the conformations of the tails and the loops of the nucleus are assumed to be the same as those without the ordered region. Under this assumption, our model nucleus is just like a flower micelle as shown in Fig. 1 (b). The statistical weight of the conformation is evaluated by using the transfer matrix T where the polymer chain is regarded as a sequence of rod-like segments [22]. Hereafter we will refer the rod-like segment as simply ‘segment’ unless otherwise noted. The orientation vector of each segment is chosen from the 12 basis vectors of a face centered cubic lattice eðηÞ (η ¼ 1; 2; ⋯12), where the basis vectors are pffiffiffi pffiffiffi pffiffiffi defined as eð1Þ ¼ ðb; b; 0Þ= 2; eð2Þ ¼ ð b; b; 0Þ= 2; eð3Þ ¼ ðb; 0; bÞ = 2; pffiffiffi ð5Þ pffiffiffi ð6Þ pffiffiffi ð4Þ ðξÞ e ¼ ð b; 0; bÞ= 2; e ¼ ð0; b; bÞ= 2; e ¼ ð0; b; bÞ= 2 and e ¼ eðξ 6Þ (ξ ¼ 7; 8; ⋯12). The transfer matrix is defined as 8 < 1 ðif neighboring 2 bonds are in trans conformationÞ T ηξ ¼ δ ðif neighboring two bonds are in gauche conformationÞ ; : 0 ðotherwiseÞ

¼

1 1 12 α!

​ �U

N X mðαþ1ÞN X mðαþ1Þ

N X mðαþ1Þ

X

⋯ n0 ¼0

n1 ¼0

ηη0 ðn0 ÞJ η0 η1 ðn1 ÞJ η1 η2 ðn2 Þ⋯J ηα

​ � δ0;N

(12)

nαþ1 ¼0 η;η’;η0 ;η1 ;⋯;ηαþ1 1 ηα

ðnα ÞU

ηα ηαþ1 ðnαþ1 Þ

mðαþ1Þ ðn0 þn1 þ⋯nαþ1 Þ ;

where the prefactor 1=12 is a normalization factor according to the 12 possible orientations of the initial segment and 1=α! is the correction of the overcounting of the micro states for loops. U ηξ ðnÞ and J ηξ ðnÞ are the statistical weights of a tail and of a loop, respectively, where both of these are composed of n segments. The expressions of these statistical weights are as follows:

(7)

tail; ​ ​ U ηξ ðnÞ ¼

where η and ξ specify the orientations of 2 consecutive segments, respectively. Moreover, 1 and δ are the statistical weights of the trans and the gauche conformations. The value of δ is described by using the energy difference between trans and gauche conformations Δε as h Δεi δ ¼ exp : (8) T

loop; ​ ​ J

ηξ ðnÞ

1 ðT n Þηξ ; ð1 þ 4δÞn

¼

1 ð1 þ 4δÞn

(13)

Z dq½ðT~ ðqÞÞn �ηξ :

(14)

Using these definitions, eqn. (12) is rewritten as � Δε Z α; m; ; N; M ¼ 1 T � ��� � ��α � � � N X mðαþ1Þ� 1 1 1 ¼ U� ηα ηαþ1 p U� η0 η1 p J� p 12 α! N mðα þ 1Þ þ 1 p¼0 η1 ηα �

We will refer Δε as simply ‘gauche energy’. The persistence length of the chain is related to Δε through hΔεi 1 lp ¼ ∝exp : (9) δ T

� �� 2π pðN mðα þ 1ÞÞ �exp i ; N mðα þ 1Þ þ 1

By using the transfer matrix, the path integral for a polymer chain composed of N þ 1 segments (N 6¼ 0) Qð0; η; r; N; ξ; r’Þ is defined as follows: Z N dqT~ ηξ ðqÞexp½iq ⋅ ðr’ rÞ�; (10) Qð0; η; r; N; ξ; r’Þ ¼

where � means discrete Fourier transformation with respect to the segment index, for example,

where

U� ηη0 ðpÞ ¼

(15)

N X mðαþ1Þ

U n¼0

3

� ηη0 ðnÞexp

i

N

� 2π pn : mðα þ 1Þ þ 1

(16)

H. Yokota and T. Kawakatsu

Polymer xxx (xxxx) xxx

Fig. 3. The excess free energy of a nucleus obtained using Model S for various values of the gauche energy. The vertical axis is the height of the nucleus m and the horizontal axis is the number of loops α. The parameters are set as N ¼ 127, Δh ¼ 23:75, σ s ¼ σ t ¼ 1:0 and T ¼ 0:90. It should be noted that the ranges of color bars for (a) and (b) are different from that for (c). (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

This Model S gives us the information about the critical nucleus composed of a single polymer chain. In the next subsection, we construct a theoretical model of the nucleus of multi-chain system based on this Model S.

� � Δε Z α; m; ; N þ 1; M T � � α � X α! Δε ¼ Z α1 ; m; ; N þ 1; M ¼ 1 M! α1 ;α2 ;⋯;αM ¼0 T

2.2. Nucleation theory for multi-chain system (model M)

� � � � Δε Δε �Z α2 ; m; ; N þ 1; M ¼ 1 � ⋯ � Z αM ; m; ; N þ 1; M ¼ 1 T T � �δ0;N�M mðαþMÞ FðN mðα1 þ1ÞþN mðα2 þ1Þþ⋯N mðαM þ1ÞÞ

In this subsection, we model the nucleus composed of monodispersed polymer chains which are statistically independent with each other (Model M). A schematic picture of the model of nucleus in our model M is shown in Fig. 2. In Model M, the number of chains participating in the nucleus can fluctuate. Thus, we should introduce the chemical potential conjugate to the number of chains as an independent variable. This choice of the variable means that we use the grand canonical ensemble to evaluate the conformation entropy. In the case of the nucleus with α loops and with length of ordered region m � b, the total grand potential difference before and after the nucleation is as follows: ΔΩðα; m; μc Þ ¼ ðα þ hMiÞmb3 Δμ pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 þ2 π ðα þ hMiÞmb σ s þ 2ðα þ hMiÞb2 σ t � � Δε T lnΞ α; m; ; N þ 1; μc ; T

(19)

� � α � X α! Δε ¼ Z α1 ; m; ; N þ 1; M ¼ 1 M! α1 ;α2 ;⋯;αM ¼0 T � � � � Δε Δε �Z α2 ; m; ; N þ 1; M ¼ 1 � ⋯ � Z αM ; m; ; N þ 1; M ¼ 1 T T � PM � �δ 0; α

¼ (17)

i¼1

αi

� � α h � �iM 1 X Δε 2 π sα Z s; m; ; N þ 1; M ¼ 1 exp i ; αþ1 M! α þ 1 s¼0 T

α!

(20)

where means the discrete Fourier transformation with respect to α and, s is the conjugate variable to α. The expression of Ξðα; m; Δε =T; N þ1; μc Þ is as follows: � � ∞ � X � � Δε μM � Δε Ξ α; m; ; N þ 1; μc ¼ exp c Z α; m; ; N þ 1; M ; (21) T T T M¼0

where μc is the chemical potential conjugate to the number of chains and Ξðα; m; Δε =T; N þ1; μc Þ is the grand partition function of the nucleus including the effect of the conformation entropy of the tails and loops. h Mi is the average number of the chains participating in the nucleus: � � � �i ∂ h Δε lnΞ α; m; ; N þ 1; μc : (18) M ¼ ∂μc T

where exp½μc =T� is the so-called ‘fugacity’. It is noted that the interaction between chains are described not only by bulk energy difference Δμ but also by the chemical potential μc . By calculating eqn. (17) for each α þhMi and m, we can obtain the information about the critical nucleus composed of multiple chains. In present Model M, the number of chains participating in the nucleus can fluctuate. Such a fluctuation effect is essentially important in describing critical nucleus, where the number of chains in the nucleus is not a fixed value. As the chain conformation obeys the ideal chain statistics due to the screening effect in the melt state, we can evaluate the conformation entropies of the tails and the loops using the Gaussian chain statistics. Even if the Gaussian chain statistics can be used, it is not completely

Here we ignore the contribution from CNT terms to the value of hMi to simplify the calculation. To evaluate Ξðα; m; Δε =T; N þ 1; μc Þ, we start from the partition function of the nucleus composed of M chains, where the number of loops and the hight of the nucleus are α and m, respectively (For example, in the case of Fig. 2, M ¼ 3 and α ¼ 3.). We denote the ca­ nonical partition function of a nucleus with M chains as Zðα; m; Δε =T; N þ 1; MÞ. As in the case of Model S, the width of the ordered region is also assumed to be negligibly thin also in the present Model M. Then, we obtain

4

H. Yokota and T. Kawakatsu

Polymer xxx (xxxx) xxx

Fig. 4. (left panel)The dependence of α� on Δε=T. (right panel)The dependence of m� on Δε=T. In both panels, the parameters are the same as the ones in Fig. 3. m� tends to increase with the gauche energy, while α� fluctuates.

realistic due to the topological constraint. For real polymer chains, the statistics of their conformations is affected by the non-crossing nature of the chains due to the excluded volume interaction. Such a topological constraint is especially noticeable for ring polymers [23], which can be a model for loop chains. As it is difficult to introduce the topological effect in our path integral in Model M, this effect will be considered in our future work, and then, in the present study, the conformation entropies of the loops are calculated based on the ideal chain statistics of ring polymers.

The position of the critical size ðα� ; m� Þ depends on the gauche en­ ergy Δε=T. In the case of the parameters Δh ¼ 23:75, T ¼ 0:90 and σ s ¼ σt ¼ 1:0, the values of α� and m� depend on the gauche energy as shown in Fig. 4. The height of the critical nucleus m� tends to increase with the gauche energy, while α� fluctuates with Δε=T. The fluctuation of α� is due to the competition between the energy gain of bulk energy difference in the CNT terms and the conformation entropy. When α increases, the energy gain from the bulk energy dif­ ference increases, while the conformation entropy loss also increases. According to the increase in Δε=T, the conformation entropy in eqn. (2) changes. As a result, α� and m� change with the gauche energy Δε=T. It should be noted that, in the case of the CNT parameters that are the same as those in Fig. 4, the contributions from the CNT terms and the conformation entropy term are almost the same when α changes. These contributions lead to the fluctuation of α� . We show the behavior of the activation energy Δf � as a function of the gauche energy Δε=T in Fig. 5. The behavior of Δf � shown in Fig. 5 is explained as follows. Since the parameters of the CNT are fixed, the activation energy is mainly affected by the conformation entropy term in eqn. (2). When the value of Δε=T is large, the loss of the conformation entropy is relatively large compared with that in the case of small Δε=T. Thus, the activation energy is a monotonically increasing function of Δε=T as is shown in Fig. 5. Δf � is related to the induction time τ of the nucleation, which is defined as the period before the critical nucleus is generated, as

3. Result and discussion 3.1. Result of model S In this subsection we discuss the size and the occurrence probability of the critical nucleus in Model S. In the present work, as we are inter­ ested in the dependence of the effect of the conformation entropy on Δε= T, we evaluate the free energy difference eqn. (6) when the parameters Δh, σ s and σ t are fixed. In the present work, we mainly discuss the nucleation in a shallow quench case, i.e. 1 T ¼ 0:10. First we show the excess free energy of a nucleus in Fig. 3. The size of the critical nucleus is specified by the saddle point of the excess free energy Δf ðα; mÞ shown in Fig. 3. In general, the saddle point of a continuous function gðx; yÞ is obtained by ∂g=∂x ¼ ∂g=∂y ¼ 0 and by λ1 � λ2 < 0 where λ1 and λ2 are eigen values of Hessian matrix of gðx;yÞ, i.e., 0 2 1 ∂ g ∂2 g B 2 C B ∂x ∂x∂y C C: H ¼B (22) B 2 C @ ∂ g ∂2 g A ∂x∂y ∂y2

τ ¼ τ0 exp½Δf � �;

(24)

where τ0 is the atomistic time scale. The dependence of the induction

In our model, the variables α and m in Δf are discrete values, thus, the saddle point (α� ;m� ) cannot be calculated by using the derivatives of Δf with respect to α and m. Instead, the activation energy in our model is defined as the lowest energy barrier among the all paths from ðα; mÞ ¼ ð0; 0Þ to ðα’; m’Þ, where ðα’; m’Þ are the values on the other side of the energy barrier. The values of ðα� ; m� Þ which specify the activation en­ ergy give the size of critical nucleus. To obtain the activation energy and the critical size, we should choose the lowest energy barrier in the paths from ðα; mÞ ¼ ð0; 0Þ to ðα’; m’Þ. Since ðα’; m’Þ specify the values on the other side of the energy barriers, ðα’; m’Þ satisfy 8 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi < ðα’Þ2 þ ðm’Þ2 < L : (23) : α’ > 0 m’ > 0 Here, L > 0 is chosen as an appropriate value that the size of the critical nucleus does not change. Here, to search all paths from ðα; mÞ ¼ ð0; 0Þ to ðα’; m’Þ, the backtracking algorithm [24] is applied. It is noted that, by choosing the small value of L, we can reduce the computational cost for searching the all paths from ðα; mÞ to ðα’; m’Þ.

Fig. 5. The dependence of Δf � on Δε=T. The parameters of the CNT are the same as those in Fig. 4. 5

H. Yokota and T. Kawakatsu

Polymer xxx (xxxx) xxx

Fig. 6. The grand potential difference ΔΩ for the various values of the gauche energy Δε=T in Model M. The vertical axis is the height of the nucleus m and the horizontal axis is the number of straight parts in the cylindrical ordered region α þ hMi. The other parameters are set as N ¼ 127, Δh ¼ 23:75, σs ¼ σ t ¼ 1:0, T ¼ 0:90 and μc =T ¼ 0:0.

Fig. 7. (left panel)The dependence of the number of the straight parts in the critical nucleus α� þhMi on the dimensionless gauche energy Δε= T in Model M. (right panel) The dependence of the height of the critical nucleus m� on the dimensionless gauche energy in Model M. The parameters are the same as the ones in Fig. 6.

time on the gauche energy is qualitatively the same as the one of the activation energy. In experimental studies, the induction time is usually obtained by using DSC (Differential Scanning Calorimetry) technique by observation of the latent heat, where the observed induction time cor­ responds to the generation of several nuclei. On the other hand, the induction time in our model is the time by a single critical nucleus ap­ pears. Although this induction time is not perfectly the same as the one

obtained with DSC technique, it qualitatively corresponds to the experimental values. Moreover, the gauche energy depends on the chemical species of polymers. The result shown in Fig. 5 leads to the relationship between the observable “induction time” and the micro­ scopic parameter Δε=T.

Fig. 8. The dependence of the activation energy ΔΩ� on the gauche energy Δε= T in Model M. The vertical axis is ΔΩ� and the horizontal axis is the Δε= T. The parameters are the same as the ones in Fig. 6.

Fig. 9. The grand potential difference for μc =T ¼ 0:10. The parameters of classical nucleation terms are the same as the ones in Fig. 6 and Δε=T ’ 1:5. 6

H. Yokota and T. Kawakatsu

Polymer xxx (xxxx) xxx

Fig. 10. (left panel)The dependence of the number of straight parts of the critical nucleus α� þhMi on Δε=T for μc =T ¼ 0:00 (purple), 0.10 (green) and 0.20 (blue) in Model M and (right panel)the dependence of the critical nucleus m� . The parameters of CNT are the same as the ones in Fig. 6. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.)

3.2. Result of model M

small μc =T region. The activation energy in Model M is expected to depend on μc =T since hMi increases with μc =T as is shown in Appendix B. The increase in hMi with μc =T is explained by the contribution from the fugacity in eqn. (21). As a large value of μc =T leads to a large value of α þhMi (see Fig. 10), the decrease in bulk free energy term leads to a decrease in the activation energy (see Fig. 11). The behavior of the size of the critical nucleus for each μc =T is shown in Fig. 10. In the left panel of Fig. 10, α� þhMi varies with Δε=T due to the fluctuation of the number of loops. On the other hand, α� þhMi increases with μc =T as long as α� is constant (for example, see the value of α� þhMi at Δε=T ’ 1:5). The contribution to the increase in α� þhMi from hMi is dominant, which means that μc =T can be treated as the driving force for a polymer chain to participate in the nucleus. The gauche energy which changes m� increases with the chemical potential μc =T as shown in the right panel of Fig. 10. The behaviors of the activation energy for each μc = T are shown in Fig. 11. According to Fig. 11, we recognize that ΔΩ� monotonically decreases with μc =T. In Model M, the size of the critical nucleus is determined by Δε=T and by μc =T, while in Model S, the critical size is determined only by Δε=T.

In this subsection, the results obtained from Model M are discussed. First, we discuss the grand potential difference before and after the nucleation, ΔΩ, for the case Δh ¼ 23:75, T ¼ 0:90, σs ¼ 1:0, σ t ¼ 1:0 and μc =T ¼ 0:0 in Fig. 6. The size of the critical nucleus ðα� þhMi; m� Þ is expected to depend on the gauche energy Δε=T and the chemical potential μc = T. We show the relationship between the critical sizes and the gauche energy in Fig. 7 in the case of μc =T ¼ 0. Fig. 7 implies that α� þhMi varies with gauche energy due to the change in the number of loops. Moreover, the height of the critical nu­ cleus m� tends to increase with gauche energy. These dependences of α� þhMi and m� on Δε=T are due to the same mechanism as in Model S. We show the dependence of the activation energy of Model M on the gauche energy in Fig. 8. The plateau of ΔΩ� corresponds to α� ¼ 0 ðα� þ hMi ’ 1Þ. It is noted that, as shown in Figs. 5 and 8, ΔΩ� in Model M is smaller than Δf � in Model S, which is due to the competition between the CNT terms and the conformation entropy. For example, the conformation entropy in Model M is larger than the one in Model S because the chain does not neces­ sarily generate a high-energy hairpin loop in Model M. Next, we discuss the effect of the chemical potential μc = T. We show in Fig. 9 the grand potential difference for the cases of μc = T ¼ 0:10 and Δε=T ’ 1:5. It is noted that the free energy landscape does not change in

4. Conclusion We extended the classical nucleation theory (CNT) which does not contain the conformation entropy and the information of the number of chains participating in the nucleus. We proposed 2 theoretical models of nucleation behavior. One is the model for a nucleus composed of a single chain (Model S), and the other is the model for a nucleus composed of multiple chains (Model M). In these models, the nucleus is composed of tails, loops and an ordered region whose shape is assumed to be a cyl­ inder, where the width of cylindrical ordered region is assumed to be negligibly thin. By using these models, we discuss the behavior of the critical nucleus in terms of the chain stiffness and of the number of chains participating in the nucleus which are not taken into account in Muthukumar’s model [20]. In Model S, we obtain the information about the critical nucleus by calculating the free energy difference before and after the nucleation. The size of the critical nucleus is determined by the competition between the CNT terms and the conformation entropy term. In our model, the number of loops α and the height of the ordered region m are discrete values. Thus, the critical size and the activation energy cannot be determined by using derivative of the free energy difference with respect to α and m. Then, the activation energy is defined by the lowest free energy among all paths from ðα; mÞ ¼ ð0; 0Þ to ðα’; m’Þ where the latter is the value on the other side of the energy barrier. Moreover, the size of the critical nucleus is defined as the values of the number of loops and the height of the nucleus which correspond to the saddle point. The activation energy is a monotonically increasing function of the chain

Fig. 11. The dependence of the activation energy ΔΩ� on the dimensionless gauche energy Δε=T for μc =T ¼ 0:00 (purple), 0.10 (green) and 0.20 (blue) in Model M. The parameters of CNT are the same as the ones in Fig. 6. (For interpretation of the references to color in this figure legend, the reader is referred to the Web version of this article.) 7

Polymer xxx (xxxx) xxx

H. Yokota and T. Kawakatsu

stiffness. This result is interpreted by using the fact that the conforma­ tion entropy loss increases with the chain stiffness. In Model M, we introduced the chemical potential μc which is con­ jugate to the number of chains participating in the nucleus. The grand potential difference discussed in Model M for μc =T ¼ 0:0 is smaller than the free energy difference evaluated in Model S. This means that the critical nucleus can be formed more easily in the multi-chain system than in the single chain system due to the difference in the conformation entropy loss associated with the creation of loops. In Model M, the size of the ordered region is proportional to the average of the number of chains participating in the nucleus hMi. As the chemical potential plays a role of a driving force for a polymer chain to participate in the nucleus, hMi is an increasing function of the chemical potential. Thus, the high μc leads to a large radius of the cylindrical ordered region despite of our assumption that the size of the ordered region is negligibly small. Therefore, when we study the crystallization process for high μc , we cannot approximate the conformation entropy of the nucleus as the one in free space due to the effect of the excluded

volume of the ordered region. The number of chains participating in the nucleus is determined by the following 2 contributions; The first one is the contribution from the competition between surface tensions and bulk energy difference which is controlled by the temperature. The other contribution is the chemical potential conjugate to the number of chains. If the state of the nucleus is specified, e.g. the size of the ordered region, hMi is obtained from the derivative of the grand potential difference with respect to μc . Declaration of competing interest We have no conflict of interest, financial or otherwise. Acknowledgement This work is partially supported by the Grant-in-Aid for Scientific Research (Grant Number 26287096 and 16K13844) from The Ministry of Education, Culture, Sports, Science and Technology(MEXT), Japan.

Appendix A. Derivation of Zðα; m; Δε =T; N þ 1; M ¼ 1Þ The conformation of a flower micelle is calculated by using the path integral Q defined in eqn. (10) as Z1 ðN þ 1 mðα þ 1ÞÞ Z N X mðαþ1Þ N X mðαþ1Þ N X mðαþ1Þ ¼ drdxdr’ ⋯ n0 ¼0

n1 ¼0

nαþ1 ¼0

X η;η’;η0 ;η1 ;⋯ηαþ1

​ � Qð0; η; r; n0 ; η0 ; xÞ

(A.1)

α 1 α �X � X ​ � Qðn0 ; η0 ; x; n0 þ n1 ; η1 ; xÞ⋯Q ni ; ηα 1 ; x; n i ; ηα ; x i¼0

i¼0

α αþ1 �X � X ​ �Q ni ; ηα ; x; ni ; ηαþ1 ; r’ i¼0

​ � δ0;N

i¼0

mðαþ1Þ ðn0 þn1 þ⋯nαþ1 Þ ;

where η; η’ and ηi are indices of the segment orientation. Variables r and r’ are the positions of two ends of the tails, and x is the position of the branching point of the flower micelle. By using Fourier transformation with respect to the positions of segments, we derive mðα þ 1ÞÞ

Z1 ðN þ 1 ¼

1 1 mðα þ 1Þ þ 1 12ðt þ 4gÞN

N

mðαþ1Þ

Z drdxdr’

​ �

​ � exp½iq0 ⋅ðx

N X mðαþ1ÞN X mðαþ1Þ

N X mðαþ1Þ

X

Z dq0

⋯

​ � n0 ¼0

Z

N

1 1 mðα þ 1Þ þ 1 12ðt þ 4gÞN N X mðαþ1ÞN X mðαþ1Þ n0 ¼0

​ � ðT ​ � δ0;N

n0

n1 ¼0

mðαþ1Þ ðn0 þn1 þ⋯nαþ1 Þ

dq1 ⋯

dqαþ1

N X mðαþ1Þ

1 ηα

½ðT~ ðqαþ1 ÞÞnαþ1 �ηα ηαþ1

(A.2)

mðαþ1Þ

X

Z

Z dq1

⋯

​ �

xÞ� ​ � δ0;N

nαþ1 ¼0 η;η0 ;η1 ;⋯ηαþ1

n1 ¼0

​ � ½ðT~ ðq0 ÞÞn0 �ηη0 ½ðT~ ðq1 ÞÞn1 �η0 η1 ½ðT~ ðq2 ÞÞn2 �η1 η2 ⋯½ðT~ ðqα ÞÞnα 1 �ηα

¼

rÞ�exp½iqαþ1 ⋅ðr’

Z

Z dq2 ⋯ dqα

(A.3)

nαþ1 ¼0 η;η0 ;η1 ;⋯ηαþ1

Þηη0 ½ðT~ ðq1 ÞÞn1 �η0 η1 ½ðT~ ðq2 ÞÞn2 �η1 η2 ⋯½ðT~ ðqα ÞÞnα 1 �ηα

1 ηα

ðT

nα

Þηα ηαþ1

mðαþ1Þ ðn0 þn1 þ⋯nαþ1 Þ :

By using eqns. (13) and (14), we can derive eqn. (12). Appendix B. Behavior of hMi and driving force μc In this appendix, we discuss the behavior of hMi which depends on μc =T, α and m. We show the relationship between hMi and μc = T in the cases with ðα; mÞ ¼ ð0; 3Þ; ð1; 3Þ and ð2; 3Þ, and with Δε=T ¼ 0:0. As the value of α is relatively small, the behavior of hMi does not depend on Δε= T in these sizes 8

H. Yokota and T. Kawakatsu

Polymer xxx (xxxx) xxx

of the nuclei. Although hMi depends on Δε=T in the case of the large value of α, this situation corresponds to the thick ordered region which is not consistent with the assumption of Model M.

Fig. B.12. The dependence of the average number of chains hMi participating in the nucleus on the chemical potential μc =T. The vertical axis is logarithmic scale and the horizontal axis is linear scale. The sizes of the nuclei are ðα; mÞ ¼ ð0; 3Þ (purple solid line), ðα; mÞ ¼ ð1; 3Þ (green dashed line) and ðα; mÞ ¼ ð2; 3Þ (blue dotted line). The gauche energies are Δε=T ¼ 0:0 for all sizes.

We notice the relationship between hMi and μc =T. While hMi increases with μc =T, it does not depend on the size of the nucleus. hMi is roughly an exponential function of μc =T. Thus, in the high μc =T regime, the assumption of the small ordered region in the nucleus can not be justified. We focus on the behavior of hMi at μc =T ¼ 0. As we mentioned in subsection 3.2, μc =T is regarded as a driving force for the polymer chain to migrate from the bulk phase to the nucleus. Therefore, μc =T ¼ 0 means that the driving force does not exist except for the contribution from the CNT terms. In Fig. B.12, at μc =T ¼ 0 with ðα; mÞ ¼ ð0; 3Þ (the purple solid line), we obtain hMi ’ 1 which is the minimum value for creating the nucleus. The behavior of hMi around μc =T ¼ 0 does not depend on m but on α. In the case of μc =T ¼ 0 with α ¼ 0 (see the purple solid line in Fig. B.12), hMi ’ 1 which is the minimum value for M for creating the nucleus. In the case of α ¼ 1 and 2 specified by the green dashed line and the blue dotted line, respectively, in Fig. B.12, we obtain hMi > α at μc =T ¼ 0. This result is explained as follows; When the total number of the loops is specified by α in Model M, the single chain tends not to include multiple loops due to the large loss of the conformation entropy. In other words, taking the confor­ mation entropy loss into account, the nucleus with α loops tend to be composed of at least α chains where a chain may or may not include a single loop. Therefore, the average number of chains is larger than α at μc =T ¼ 0.

References [10]

[1] G. Strobl, From the melt via mesomorphic and granular crystalline layers to lamellar crystallites: a major route followed in polymer crystallization? Eur. Phys. J. E 3 (2) (2000) 165–183, https://doi.org/10.1007/s101890070030. [2] W.-T. Chuang, W.-B. Su, U.-S. Jeng, P.-D. Hong, C.-J. Su, C.-H. Su, Y.-C. Huang, K.F. Laio, A.-C. Su, Formation of mesomorphic domains and subsequent structural evolution during cold crystallization of poly(trimethylene terephthalate), Macromolecules 44 (5) (2011) 1140–1148, https://doi.org/10.1021/ma1025175, arXiv. [3] J.-B. Jheng, W.-T. Chuang, P.-D. Hong, Y.-C. Huang, U.-S. Jeng, C.-J. Su, G.-R. Pan, Formation of mesomorphic domains associated with dimer aggregates of phenyl rings in cold crystallization of poly(trimethylene terephthalate), Polymer 54 (22) (2013) 6242–6252, https://doi.org/10.1016/j.polymer.2013.09.006. http://www. sciencedirect.com/science/article/pii/S0032386113008586. [4] M. Imai, K. Mori, T. Mizukami, K. Kaji, T. Kanaya, Structural formation of poly (ethylene terephthalate) during the induction period of crystallization: 1. ordered structure appearing before crystal nucleation, Polymer 33 (21) (1992) 4451–4456, https://doi.org/10.1016/0032-3861(92)90399-H. http://www.sciencedirect.com/ science/article/pii/003238619290399H. [5] Z.-G. Wang, B.S. Hsiao, E.B. Sirota, P. Agarwal, S. Srinivas, Probing the early stages of melt crystallization in polypropylene by simultaneous small- and wide-angle xray scattering and laser light scattering, Macromolecules 33 (3) (2000) 978–989, https://doi.org/10.1021/ma991468t, arXiv. [6] Z.-G. Wang, B. Hsiao, E. Sirota, S. Srinivas, A simultaneous small- and wide-angle x-ray scattering study of the early stages of melt crystallization in polyethylene, Polymer 41 (25) (2000) 8825–8832, https://doi.org/10.1016/S0032-3861(00) 00225-1. http://www.sciencedirect.com/science/article/pii/S0032386 100002251. [7] P. Panine, E.D. Cola, M. Sztucki, T. Narayanan, Early stages of polymer melt crystallization, Polymer 49 (3) (2008) 676–680, https://doi.org/10.1016/j. polymer.2007.12.026. http://www.sciencedirect.com/science/article/pii/ S0032386107011561. [8] J. Lauritzen, J.D. Hoffman, Theory of formation of polymer crystals with folded chains in dilute solution, J. Res. Natl. Bur. Stand. 64 (1960) 73–102. [9] N. Waheed, M. Ko, G. Rutledge, Molecular simulation of crystal growth in long alkanes, Polymer 46 (20) (2005) 8689–8702, https://doi.org/10.1016/j.

[11] [12]

[13]

[14] [15]

[16] [17]

[18]

[19]

9

polymer.2005.02.130. http://www.sciencedirect.com/science/article/pii/ S003238610500710X. M. Anwar, F. Turci, T. Schilling, Crystallization mechanism in melts of short nalkane chains, J. Chem. Phys. 139 (21) (2013), 214904, https://doi.org/10.1063/ 1.4835015 arXiv. P.M. Welch, Examining the role of fluctuations in the early stages of homogenous polymer crystallization with simulation and statistical learning, J. Chem. Phys. 146 (4) (2017), 044901, https://doi.org/10.1063/1.4973346 arXiv. M. Imai, K. Mori, T. Mizukami, K. Kaji, T. Kanaya, Structural formation of poly (ethylene terephthalate) during the induction period of crystallization: 2. kinetic analysis based on the theories of phase separation, Polymer 33 (21) (1992) 4457–4462, https://doi.org/10.1016/0032-3861(92)90400-Q. http://www.sc iencedirect.com/science/article/pii/003238619290400Q. M. Imai, K. Kaji, T. Kanaya, Structural formation of poly(ethylene terephthalate) during the induction period of crystallization. 3. evolution of density fluctuations to lamellar crystal, Macromolecules 27 (24) (1994) 7103–7108, https://doi.org/ 10.1021/ma00102a016, arXiv. T. Konishi, K. Nishida, G. Matsuba, T. Kanaya, Mesomorphic phase of poly (butylene-2,6-naphthalate), Macromolecules 41 (9) (2008) 3157–3161, https:// doi.org/10.1021/ma702383b, arXiv. M. Anwar, T. Schilling, Crystallization of polyethylene: a molecular dynamics simulation study of the nucleation and growth mechanisms, Polymer 76 (2015) 307–312, https://doi.org/10.1016/j.polymer.2015.08.041. http://www.sciencedi rect.com/science/article/pii/S0032386115301804. S. Auer, D. Frenkel, Prediction of absolute crystal-nucleation rate in hard-sphere colloids, Nature 409 (6823) (2001), 1020. J. Huang, F. Chang, Crystallization kinetics of poly(trimethylene terephthalate), arXiv, J. Polym. Sci. B Polym. Phys. 38 (7) (2004) 934–941, 10.1002/(SICI)10990488(20000401)38:7<934::AID-POLB4>3.0.CO;2-R, https://onlinelibrary.wiley. com/doi/pdf/10.1002/, https://onlinelibrary.wiley.com/doi/abs/10.1002/. J.-U. Sommer, C. Luo, Molecular dynamics simulations of semicrystalline polymers: crystallization, melting, and reorganization, J. Polym. Sci. B Polym. Phys. 48 (21) (2010) 2222–2232, https://doi.org/10.1002/polb.22104, arXiv, http s://onlinelibrary.wiley.com/doi/pdf/10.1002/polb.22104, https://onlinelibrary. wiley.com/doi/abs/10.1002/polb.22104. T. Yamamoto, Molecular dynamics of polymer crystallization revisited: crystallization from the melt and the glass in longer polyethylene, arXiv, J. Chem.

H. Yokota and T. Kawakatsu

Polymer xxx (xxxx) xxx

Phys. 139 (5) (2013), 054903, https://doi.org/10.1063/1.4816707, 10.1063/ 1.4816707. [20] M. Muthukumar, Molecular modelling of nucleation in polymers, Mathematical, Physical and Engineering Sciences, Phil. Trans. R. Soc. Lond. A 361 (1804) (2003) 539–556, https://doi.org/10.1098/rsta.2002.1149. arXiv, http://rsta.royalsoci etypublishing.org/content/361/1804/539.full.pdf, http://rsta.royalsocietypublish ing.org/content/361/1804/539. [21] G. Strobl, The Physcs of Polymers: Concepts for Understanding Their Structure and Behavior, Springer-Verlag, 2007.

[22] H. Yokota, T. Kawakatsu, Modeling induction period of polymer crystallization, Polymer 129 (2017) 189–200, https://doi.org/10.1016/j.polymer.2017.09.022. http://www.sciencedirect.com/science/article/pii/S0032386117308844. [23] J. Suzuki, A. Takano, T. Deguchi, Y. Matsushita, Dimension of ring polymers in bulk studied by monte-carlo simulation and self-consistent theory, J. Chem. Phys. 131 (14) (2009), 144902, https://doi.org/10.1063/1.3247190, 10.1063/ 1.3247190. [24] S.W. Golomb, L.D. Baumert, Backtrack programming, J. ACM 12 (4) (1965) 516–524, https://doi.org/10.1145/321296.321300.

10