A computational study of short-range surface-directed phase separation in polymer blends under a linear temperature gradient

Chemical Engineering Science 137 (2015) 884–895

Contents lists available at ScienceDirect

Chemical Engineering Science journal homepage: www.elsevier.com/locate/ces

A computational study of short-range surface-directed phase separation in polymer blends under a linear temperature gradient Mohammad Tabatabaieyazdi, Philip K. Chan n, Jiangning Wu Department of Chemical Engineering, Ryerson University, 350 Victoria Street, Toronto, Ontario, Canada M5B 2K3

H I G H L I G H T S

G R A P H I C A L

A B S T R A C T


Model composed of Cahn–Hilliard and Flory-Huggins-deGennes theories is solved.
Growth is fastest in the early stage of phase separation.
Thickness of the wetting layer increased with increasing temperature.
Transition time occurs earlier with higher surface potential.
Spinodal wave became more visible in the bulk as the surface potential increased.

art ic l e i nf o

a b s t r a c t

Article history: Received 29 January 2015 Received in revised form 13 June 2015 Accepted 14 July 2015 Available online 26 July 2015

The nonlinear Cahn–Hilliard theory and the Flory-Huggins-de-Gennes theory were used to study numerically the surface-directed phase separation phenomena of a model binary polymer blend quenched into the unstable region of its binary symmetric upper critical solution temperature phase diagram. Short-range surface potential within a square geometry, where one side of the binary polymer blend is exposed to a surface with preferential attraction to one component of the blend that is under a linear temperature gradient along the direction perpendicular to the surface, was integrated into the model. The structure factor analysis showed a faster exponential growth at the early stage of phase separation and a slower growth rate at the intermediate stage with a slope of 0.31 within the bulk, which is consistent with the Lifshitz–Slyozov growth law. The investigation of surface enrichment rate at the surface wall demonstrated faster growth rate at the early stage with the slope of 0.5. This growth rate became slower at the intermediate stage with a slope of 0.13 near the surface. The effect of various temperature gradient values on the surface enrichment rate with constant temperature T1* at the surface preferentially attracting one of the polymer components and different temperature T2* at the opposite surface, where T1*4T2*, was studied for the first time. The results showed that the thickness of the wetting layer increased with increasing temperature difference ΔT*, where ΔT*¼ T1*  T2*. The structure factor analysis of the surface potential h1 effect on the phase separation within the bulk close to the surface showed earlier transition time for higher values of h1. However, there was no difference observed for transition time within the bulk at distances farther away from the surface. As the surface potential increased, spinodal wave became more visible in the bulk and the transition time from complete wetting to partial wetting occurred at a later time on the surface. & 2015 Elsevier Ltd. All rights reserved.

Keywords: Surface-directed phase separation Spinodal decomposition Polymer blend Surface potential Wetting Morphology development

n

Corresponding author. E-mail address: [email protected] (P.K. Chan).

http://dx.doi.org/10.1016/j.ces.2015.07.023 0009-2509/& 2015 Elsevier Ltd. All rights reserved.

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

1. Introduction The surface-directed phase separation (SDPS) mechanism has attracted much attention and has been intensively studied (Jones et al., 1991; Krausch et al., 1994; Geoghegan et al., 1995, 2000; Brown and Chakrabarti, 1992; Puri and Binder, 1992; Henderson and Clarke, 2005; Yan and Xie, 2006a,b). The presence of a surface may alter the course of phase separation by spinodal decomposition in polymer blends by breaking translational and rotational symmetry. SDPS in polymer blends at a range of conditions could lead to a variety of structural morphology. Therefore, it is important to understand the phase separation behavior of polymeric materials and also the effect of different external fields (such as surface effect and temperature gradient) that could lead to the formation of structural anisotropy in polymer blends and could help produce new products with enhanced property and functionality. The adsorption of polymers onto surfaces, whether preferred or not, has great consequences in polymer product formation. Moreover, understanding and controlling such processes is significant and is necessary in many technological features varying from paper industry and paint formulation to pharmaceutical applications (Dalmoro et al., 2012; Evans and Wennerstrom, 1999), biophysics (Norde, 2003; Gray, 2004; Yaseen et al., 2008), and nanocomposite materials (Mishra et al., 2013; Ramanathan et al., 2008). Phase separation induced through a temperature jump (for polymer blends with a lower critical solution temperature) or quench (for polymer blends with an upper critical solution temperature) into the unstable region of the phase diagram is known as thermallyinduced phase separation (TIPS). Phase separation in the TIPS method occurs via spinodal decomposition (SD); this particular process of phase separation does not require an activation energy and proceeds spontaneously in the presence of minimal concentration fluctuations or thermal noise (Cahn, 1965). Phase separation of polymer blends can lead to different morphologies, such as the bi-continuous interconnected structure and the droplet-type morphology when altering the system characteristics such as the composition, molecular weight and structure, film thickness, solvent, or changes in the exterior environment, including the substrate, pressure, temperature, and external fields. This offers a means to create patterns in polymeric materials by controlling the phase separation conditions in thin polymer blend films. The first experimental observation of surface-directed phase separation (SDPS) via spinodal decomposition (SDPS) was reported by Jones et al. (1991). They found that in a poly(ethylenepropylene) and perdeuterated poly(ethylenepropylene) blend there is a preferential attraction of the latter component to the surface, and that phase separation by spinodal decomposition occurs with a wavevector normal to the surface which extends for some distance into the bulk. Since then, there have been numerous experimental (Krausch et al., 1994; Geoghegan et al., 1995, 2000) and numerical (Krausch et al., 1994; Henderson and Clarke, 2005; Yan and Xie, 2006a,b) work published on the SDPS of polymeric materials with various polymers and processing conditions. Krausch et al. (1994) studied the SDPS in a symmetric poly (ethylenepropylene) and perdeuterated poly(ethylenepropylene) blend with off-critical concentrations. They found experimentally that the growth rate of the wetting layer grows slower if the minority phase wets the surface, while the reverse is true if the majority phase wets the surface. Furthermore, their numerical studies indicate that this growth rate follows the t1/3 diffusive scaling law when the majority phase wets the surface. Geoghegan et al., (1995, 2000) studied the quench depth effect of the SDPS in a deuterated polystyrene and poly (α-methylstyrene) blend. They found that the growth law of the surface wetting layer follows t1/3 for the deepest quenches, while the growth law is logarithm with time for the shallower quenches. Brown and Chakrabarti (1992) studied numerically the SDPS in a two-dimensional model incorporating the Cahn-Hillard-Cook

885

theory for spinodal decomposition and a free energy functional composing of the Ginzburg-Landau free energy and a long range surface interaction term. They found that the thickness of the wetting layer varies with t1/3. Moreover, although the length scales l(t) in the direction parallel and perpendicular to the surface are different, they both scale as l(t) t1/3. Puri and Binder (1992) used the Cahn-Hilliard equation derived from a semi-infinite model with Kawasaki spin exchange dynamics to study the SDPS. Their numerical results confirm that there exists a concentration wave growing with wavevector perpendicular to the surface. Henderson and Clarke (2005) model the SDPS in a symmetric polymer blend using a model composed of the Cahn-Hilliard-Cook theory for phase separation and the Flory-Huggins-de Gennes free energy. They concluded that to model the SDPS and to be capable of replicating experimental observations, only one length scale should be used. Yan and Xie (2006a,b) simulated the SDPS in a polymer system by cell dynamic system. For the case of short range potential, they found that the wetting layer grows according to the logarithmic growth law if the noise term is absent, while it will grow according to the Lifshitz–Slyozov 1/3 growth law if the noise term is present. Furthermore, the length scale parallel to the surface also obeys the Lifshitz–Slyozov 1/3 growth law as well. They also found that the thickness of the wetting layer and degree and speed of phase separation all increase with surface potential. According to the authors’ knowledge of published work, the temperature has always remained constant during the SDPS in polymeric materials. Lee et al. (2002, 2003, 2004) and Hong and Chan (2010) studied numerically the effect of a temperature gradient on the thermal-induced and polymerization-induced phase separation processes in polymer solutions. The presence of a temperature gradient during phase separation by spinodal decomposition significantly alters the transient morphology being formed. Temperature gradients inherently exist during the fabrication of polymeric products, such as in the extrusion, and injection molding and blow molding processes. Moreover, temperature gradients may be introduced deliberately in the phase separation process to strategically achieve a graded morphology which is required in product specifications. In this present paper, we present and solve a mathematical model composed of the Cahn-Hilliard (CH) theory for spinodal decomposition (Cahn, 1965) and the Flory-Huggins-de Gennes (FHdG) free energy function (deGennes, 1980; Flory, 1953). The model also incorporates the following two external forces: (1) a surface with short-range surface potential, and (2) a linear temperature gradient normal to the surface. The wetting layer formation mechanisms on the surface and morphology development and evolution are examined under different conditions. The results are presented and discussed in the form of morphology formation and surface enrichment growth rate. The effects of diffusion coefficient, quench depth, temperature gradient and surface potential on the surface enrichment are as well investigated.

2. Model development This section explains the model development for the thermalinduced SDPS method in a binary polymer blend, involving a short-range surface potential field and an externally imposed spatial linear temperature gradient. This model is developed using the nonlinear CH theory for spinodal decomposition (Cahn, 1965) and the FHdG free energy (deGennes, 1980; Flory, 1953). The CH theory is derived using the following continuity equation: ∂c ¼ ∇Uj ∂t

ð1Þ

where c is the concentration of the solvent taken as volume fraction in this paper, and j is the interdiffusional flux. The flux

886

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

binary polymer blend is obtained:

j is related to the gradient in chemical potential through:

δF j ¼  M∇ μ2  μ1 ¼  M δc




ð2Þ

where M is the concentration dependent mobility, and m1 and m2 are the chemical potentials of components 1 and 2, respectively. The total free energy F is obtained from the FHdG theory: Z
 ð3Þ F¼ f ðcÞ þ κ ‖∇c‖2 dV where f(c) is the homogeneous Flory-Huggins (FH) free energy (Flory, 1953) and κ is the interfacial energy parameter. The FH free energy is expressed as:   kB T c 1c ln c þ lnð1  cÞ þ χ cð1  cÞ ð4Þ f ðc Þ ¼ N2 υ N1 where kB is Boltzmann’s constant, T is the temperature, υ is the volume of a cell or segment, N1 and N2 are the degrees of polymerization of components 1 and 2, respectively, and χ is Flory’s interaction parameter. deGennes (1980) proposed that κ is the sum of an enthalpic term relating to the effective range of the interactions and an entropic term whose origin is the configurational entropy of the Gaussian coils:

κ ðcÞ ¼ κ entropic þ κ enthalpic ¼ In polymer blends, becomes:

κ ðc Þ ffi

a2 36cð1 cÞ

a2 þ a2 χ 36cð1  cÞ

ð5Þ

κentropic c κenthalpic (deGennes, 1980)and Eq. (5) ð6Þ

Eq. (6) can be used in the CH equation to predict phase separation of polymer blends more accurately. Diffusion of a single polymer chain in the blend can be approximated using the reptation theory: Dreptation ¼

kB Ta2 3N ξb 2

2

ð7Þ

where a is the step length of a primitive chain, b is the bond length, N is the number of monomers in the chain, and ξ is the frictional coefficient per polymer chain. The reptation theory indicates that D is related to N  2. Since lengths a and b are related by: 2

a2 ¼

4N e b 5

ð8Þ

the reptation model can be expressed for each component i as: ! 4kB T N e;i Di ¼ ð9Þ 15ξi Ni 2 where Ne is average number of monomers existing between each two entanglement points of the polymer chains. It should be noted that reptation behavior of polymer mixture occurs when N4Nc, where Nc is the critical degree of polymerization. The approximate value of N c is 300 monomer units (Brochard et al., 1983; Klein, 1978). The diffusion coefficient could be expressed as:  2  ∂ f ðc Þ D¼M ð10Þ 2 ∂c If mutual diffusion of the binary polymer blend is controlled by the slower moving component, the slow mode theory by deGennes (1980) may be used to determine the mobility (Kramer et al., 1984): 1 1 1 ¼ þ M M1 M2

ð11Þ

where M1 and M2 are the mobility of components 1 and 2, respectively. Combining Eqs. (9)–(11) and considering that Ne,1 ¼ Ne,2 ¼ Ne and ξ1 ¼ ξ2 ¼ ξ, the following expression for the mobility in a

M¼

4υN e cð1  cÞ 15ξ½N 2 c þ N 1 ð1  cÞ

ð12Þ

The linear temperature gradient used in this study is expressed as:   T2 T1 ðx x1 Þ þ T 1 ð13Þ T ðxÞ ¼ x2  x1 where T1 and T2 are temperatures at x1 and x2, respectively. The following scaling relations are used below to non-dimensionalize the governing equation: ! 4υkB N e T c a2 Dimensionless time : t n ¼ t ð14aÞ 15ξL4 kB T c L 2 Dimensionless diffusion coefficient : D ¼   υ kBυT c a2

ð14bÞ

Dimensionless concentration : cn ¼ c

ð14cÞ

Dimensionless temperature T n ¼

T Tc

ð14dÞ

Dimensionless length : xn ¼

x L

ð14eÞ

Dimensionless length : yn ¼

y L

ð14fÞ

where Tc is the critical temperature and L is the sample length which is normally 10–100 mm. The dimensionless diffusion coefficient D is the driving force for phase separation, and is one of the parameters that will be varied in this paper. The following form of the dimensionless interaction parameter is used in the model:

χ ¼ αþ

β

Tn

ð15Þ

where α and β are constants determined experimentally and represent the entropic and enthalpic contributions, respectively. The value of the parameter α for a binary polymer blend without any specific intermolecular interactions is –1 r α r1 and parameter β has an order of magnitude of 10  1. An example of a low noise binary polymer mixture is a perdeuterated polybutadiene (DPB) and protonated polybutadiene (HPB) blend. For this mixture, α and β (after normalization) have values of –5.34  10  4 and 8.44  10  4, respectively (Jinnai et al., 1993). These values of α and β are fitted to the expected linear dependence on T  1 used in the model simulations. Lastly, it should be noted that Flory’s interaction parameter is weakly dependent on concentration in polymer blends and thus not included in Eq. (15). The following partial differential equation is obtained by combining Eqs. (1)–(4), (6), and (12)–(15) to govern the phase separation process in a polymer blend under an externally imposed linear temperature gradient:  
 ∂cn 1 1 ln cn lnð1  cn Þ n ¼ D  þ  þ α 1 2c N1 N2 N2 ∂t n N1 ! n ðN 1  N2 Þc þ N1 ð1  2cn Þ ∇T n ∇cn ðN 2 cn þ N 1 ð1 cn ÞÞ2   c n ð1  c n Þ ∇T n ∇cn  4αD N 2 c n þ N 1 ð1  c n Þ ! ðN 1  N2 Þcn2 þ N 1 ð1 2cn Þ þD ðN 2 cn þ N1 ð1  cn ÞÞ2  
2 1 1  2 þ χ T n ∇cn N 1 c n N 2 ð1  c n Þ

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

 þD

cn ð1 cn Þ N 2 cn þ N 1 ð1 cn Þ



!
2 1 1 þ T n ∇cn N 1 cn2 N 2 ð1  cn Þ2   cn ð1 cn Þ þD N 2 cn þ N 1 ð1 cn Þ ! 1 1 þ  2 χ T n ∇2 c N1 cn2 N2 ð1  cn Þ2   4kB N e a2 1 þ n n 18c ð1  c Þ 15ξ ! ðN 1  N 2 Þcn2 þ N1 ð1  2cn Þ ∇T n ∇cn ∇2 cn ðN 2 cn þ N 1 ð1  cn ÞÞ2   4kB N e a2 c n ð1  c n Þ þ N 2 c n þ N 1 ð1  c n Þ 15ξ ! n 1 2c ∇T n ∇cn ∇2 cn 9ð1 cn Þ2 cn2 ! ðN 1  N2 Þcn2 þ N 1 ð1  2cn Þ þ ðN 2 cn þN 1 ð1  cn ÞÞ2 !
2 1  2cn T n ∇cn ∇2 cn 18ð1  cn Þ2 cn2 ! 3cn2  3cn þ1 n
n 2 2 n  T ∇c ∇ c 9ð1  cn Þ3 cn3 !   c n ð1  c n Þ 1  2cn þ T n ∇2 c n ∇2 c n N 2 c n þ N 1 ð1  c n Þ 18ð1  cn Þ2 cn2 !   1 ðN 1  N2 Þcn2 þ N 1 ð1  2cn Þ ∇T n ∇3 cn þ 18cn ð1  cn Þ ðN 2 cn þN 1 ð1  cn ÞÞ2    1 cn ð1  cn Þ ∇T n ∇3 cn þ 9cn ð1  cn Þ N 2 cn þN 1 ð1  cn Þ    1 c n ð1  c n Þ T n ∇4 c n  n n n n 18c ð1  c Þ N 2 c þ N 1 ð1  c Þ

887

concentration flux is zero and this is a no-flux boundary condition (Chan and Rey, 1995a,b):  jxn ¼ 0 ¼ ∇3 cn ¼ 0



ð19Þ

Therefore, at xn ¼0, the following boundary conditions are used:  ∂cn  h1 þgcn ¼  n ∂x xn ¼ 0 γ

ð20aÞ

∂3 cn ∂ 3 cn þ n n2 ¼ 0 3 n ∂x ∂x ∂y

ð20bÞ

The natural boundary condition is applied in addition to the noflux boundary condition at the surfaces without any preferential attraction to one of the polymer components in the blend (i.e., at x* ¼ 1, at y* ¼0 and at y* ¼ 1). As noted above, the no-flux boundary condition refers to a system in which no mass is exchanged through its boundary with the surrounding (i.e., j¼0). The natural boundary condition is obtained from the variational analysis (Chan and Rey, 1995a,b), and is expressed in generalized form as:
n ∇c U n ¼ 0

ð21Þ

where n is the outward unit normal to a bounding surface. Therefore, for the three surfaces with no preferred attraction for one of the polymer components, the following boundary conditions are applied: ∂cn ¼0 ∂xn

at

xn ¼ 1; and 0 ≤ yn ≤ 1

∂3 cn ∂ 3 cn þ ¼0 3 ∂xn ∂xn ∂yn2 ð16Þ

In Eq. (16), the terms with D are the driving force for phase separation. The infinitesimal concentration fluctuations are always present in the polymer blend initially even in the single-phase region. These infinitesimal concentration fluctuations are sufficient to drive the process of phase separation by the spinodal decomposition mechanism. Therefore, the dimensionless initial concentration can be expressed as:

 cn t n ¼ 0 ¼ cn0 þ δcn t n ¼ 0 ð17Þ where c0* is the dimensionless initial concentration and δc*(t* ¼0) represents the deviation from the average initial concentration c0* (i.e., the infinitesimally small concentration fluctuations that are present in the blend initially). In the short range surface potential case, the model square domain is composed of four surfaces (0 rx* r1, 0 ry* r1), where only one of them (x* ¼0) has a preferred attraction to one of the polymer components. Each surface will have two boundary conditions. Assuming an external surface potential in the system for the domain side with surface attraction, the first boundary condition is:  ∂cn  h1  gcn þ γ n  ¼0 ð18Þ ∂x xn ¼ 0 where h1, g and γ are the three parameters that characterize the surface phase diagram and can be related to the Flory-Huggins lattice theory (Henderson and Clarke, 2005). The second boundary condition is that no penetration of material is possible through the boundary surface. In other words, at this surface, the

∂cn ¼0 ∂yn

at

at

xn ¼ 1; and 0 ≤ yn ≤ 1

0 ≤ xn ≤ 1; and yn ¼ 0 and 1

∂3 cn ∂ 3 cn þ ¼0 ∂yn3 ∂yn ∂xn2

at

0 ≤ xn ≤ 1; and yn ¼ 0 and 1

ð22aÞ

ð22bÞ

ð22cÞ

ð22dÞ

The method of lines is applied to solve partial differential Eq. (16) computationally. The spatial related derivatives are discretized by the finite difference method; in particular, the central difference scheme method is used for the spatial discretization using point-value solution. The model domain is spatially discretized by a N  N mesh, where N ¼256 in this study. All boundary conditions are also spatially discretized. The above boundary conditions will be incorporated in the governing Eq. (16), which will then be reduced to a system of ordinary differential equations (ODEs) where only time derivatives remains. The system of ODEs can be solved by generic solvers. The stiff ODE solver called CVODE (〈https://computation.llnl.gov/casc/sundials〉), based on the Generalized Minimal Residual Iteration (GMRES) method incorporating the Backward Differentiation Formulas (BDF) method, was used to solve the proposed model at different time steps. At every time step, the local error is estimated and required to satisfy tolerance conditions. The step is reiterated with a reduced step size whenever that error test fails. The fixed-leading-coefficient form of the BDF method is also applied to increase the accuracy of the approximations of the derivatives at each time step based on their values at previous times. For all simulations, the absolute tolerance is set to 10  10 and the relative tolerance is set to 10  8. The number of nodal points is kept at 256  256, and each simulation on average took two weeks on a computer with Intel i7 CPU core, due to high nonlinearity and stiffness of the model.

888

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

Fig. 1. Phase diagram for a symmetric polymer blend with polymer degree of polymerization of N1 ¼N2 ¼ 1000, where T* is the dimensionless temperature and c* is the dimensionless concentration. The solid (dotted) curve represents the binodal (spinodal) line. The temperature range for the different quenches studied in this paper for critical (c0* ¼ 0.5) and off-critical (c0* ¼ 0.3, 0.4) concentrations are shown.

3. Results and discussion In this section, numerical results will be presented and discussed for SDPS with short-range surface potential. The focus is on the formation and evolution of the phase separated structures and the different factors that control the morphology to fabricate functional polymeric material using the SDPS method. Hence, only the early and intermediate stages of phase separation are presented and discussed in this study. In the late stage of spinodal decomposition, the morphology coarsens thus destroying the practicality and functionality of the resultant composite material. Quenching quickly to a lower temperature will stop the coarsening process and lock in the desired phase separated morphology. The interaction parameter values for the typical polymer blend DPB/HPB are used; i.e., α ¼  5.34  10  4 and β ¼8.44  10  4 (Jinnai et al., 1993). Fig. 1 shows the phase diagram used in this paper. The solid curve represents the binodal line, whereas the dotted curve is the spinodal line. Fig. 1 also shows the range for the temperature gradients investigated in this paper for the critical (c0*¼0.5) and off-critical (c0*¼0.3, 0.4) quenches. The ranges of the other parameters explored in this paper are: 0.5rh1 r6.0 and 3  105 rDr8  105. These ranges were selected to meet the objectives of this paper mentioned above and were based on values and orders of magnitudes used successfully by Chan et al. (Lee et al., 2002, 2003, 2004; Hong and Chan, 2010; Chan and Rey, 1995a,b) in previous studies under different setups. In these previous studies, Chan et al. (Lee et al., 2002, 2003, 2004; Hong and Chan, 2010; Chan and Rey, 1995a,b) found the general tendency of the dependent variables on the temperature, temperature gradient, initial concentration and diffusion coefficient. The exact functional form of the dependency of the dependent variables on these independent parameters were out of the scope of their manuscript and also in this present one. 3.1. Effect of different quench depths on surface enrichment The morphology formation of the binary polymer blend under different temperature gradients, surface attraction strengths to one polymer component, concentrations, quench depths and diffusion coefficients has been studied for short-range surface potential cases. Typical results are presented and discussed below. For cases that the attracting surface was completely wetted by a component, the rate of the growth of the wetting layer is compared with literature values. It is well known that morphology formation and its evolution following the

thermal-induced spinodal decomposition process depends on the temperature that blends are quenched to (known as quench depth), the diffusion coefficient, the degree of polymerization, and the miscibility of the two components. Furthermore, for the SDPS, the amount of surface attraction to the preferred component could affect the nature of the surface wetting (i.e., complete or partial wetting), and also the rate of the enrichment of the surface (i.e., attraction of preferred polymer to the surface) by the preferred polymer. To study the effect of the above parameters on the formation of the morphology, each parameter was varied sequentially. Phase separation induced by short-range surface attraction potential is initiated in the layer close to the surface and proceeds to include and engage the whole domain due to the intermolecular attraction forces between like polymers within the domain. The temperature at the surface with favorable attraction to one of the component is T1*, while the temperature on the other side of the domain with no surface attraction is T2*. It should be noted that there are two types of phase separation mechanisms studied through all cases of surface potential leading to different morphologies: (1) phase separation at the surface with preferential attraction to one of the polymers, and (2) phase separation within the bulk. Phase separation within the bulk is mainly governed by SD, while the phase separation at the surface is induced by both the attraction of one polymer that is preferred by that surface in addition to the SD mechanism. It has been noticed that there is a competition between the mechanisms of phase separation at the surface and within the bulk. The morphologies obtained in each case are found to be controlled mostly by the mechanism that is more dominant during phase separation. During the early stage, the surface is always prominent in initiating the phase separation prior to any bulk development and interference. During intermediate stage, however, the bulk may govern the morphology development during phase separation. According to the results of this work, domination of the bulk intensifies over time mostly for deep quenches. Fig. 2 shows a typical phase separation for the short-range surface potential case for different off-critical quench depths from the onephase region into the two-phase region with dimensionless initial concentration of c0*¼0.3. In this case, the parameter values are: D ¼4  105, h1 ¼0.5, g¼ 0.5 and T2*¼ 0.20. The temperature gradient is created by setting T1* to: (a) 0.22, (b) 0.24 and (c) 0.26. It can be observed that phase separation at the surface is initiated earlier in the form of partial wetting and also reaches complete wetting earlier in the deeper quench depth case where T1*¼0.22 than the shallow quench depths where T1*¼0.24 and 0.26. In Fig. 2, complete wetting is obtained at dimensionless time t*¼0.254 for T1*¼0.22, at t*¼0.524 for T1*¼0.24, and at t*¼1.624 for T1*¼0.26. The morphology emphasizes that increasing the quench depth will accelerate the transition from the early stage to the intermediate stage. The phase separation mechanism is primarily triggered by the surface during the early stage. During the intermediate stage, the surface still controls the phase separation in its vicinity due to the fact that there is no sign of morphology shift from complete wetting to partial wetting during this stage. The morphologies in Fig. 2 also reveal that for a shallow quench, the layer of the preferred polymer component grows thicker on the surface (i.e., large surface enrichment growth) (Geoghegan et al., 2000). As typical of off-critical quenching conditions within the spinodal region, droplets are formed during the SDPS. The obtained morphologies are consistent with published numerical (Brown and Chakrabarti, 1992; Henderson and Clarke, 2005; Yan and Xie, 2006a,b; Lee et al., 2002, 2003, 2004; Hong and Chan, 2010; Chan and Rey, 1995a,b) and experimental work (Tanaka et al., 1990; Kyu et al., 1992). Spinodal decomposition is usually followed by light scattering. Quantitative information can be obtained from the resolved light scattering intensity profile I(s, t) where s is the scattering wavevector. The numerical equivalence of this profile is the structure factor S(k, t): I ðs; t Þ p Sðk; t Þ ¼ ‖Aðk; t Þ‖2 for s ¼ k

ð23Þ

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

889

Fig. 2. Typical morphology evolution during SDPS with short range surface potential to reach complete wetting resulting from various off-critical (c0* ¼ 0.3) quench depths from the one-phase region into the two-phase region. The dimensionless parameter values are: D¼ 4  105, h1 ¼ 0.5, g ¼  0.5 and T2* ¼ 0.20. The dimensionless temperature T1* varies from (a) 0.22 for a deep quench to (b) 0.24 and (c) 0.26 for shallow quenches. The droplet morphology confirms the off-critical quench conditions.

Consequently, the time evolution of the two-phase structures, such as the ones produced during the phase separation process in the polymer blend systems studied in this work can also be measured using the structure factor. In the early stages of phase separation by SD, the dimensionless structure factor grows exponentially and the wave number k is independent of time. The growth of the concentration fluctuations is weakly nonlinear. In the intermediate stages of phase separation, the structure factor continues to increase but at a slower rate than in the early stage of phase separation. The dimensionless structure factor was calculated using MATLAB (http:// www.mathworks.com/products/matlab/) to determine the Fast Fourier Transform (FFT) at different dimensionless times. The FFT is used to confirm that the numerical work in this study is in agreement with the known evolution of scattering profiles as related to the structure factor found in the literature. Fig. 3 presents the typical evolution of the dimensionless structure factor for the off-critical quench case shown in Fig. 2(b) at different dimensionless times. It shows that in the bulk of the polymer blend undergoing SD after a deep quench into the unstable region, the dimensionless structure factor displays a maximum that grows with time. The wave number remains constant during the phase separation process in the early to the beginning of the intermediate stages, which is ubiquitous in thermal-induced SD. This indicates that the phase separation results are in the early stage and beginning of the intermediate stage of SD (Chan and Rey, 1995b). The obtained evolution of the dimensionless structure factor diagram is in good agreement with published numerical (Chan and Rey, 1995a) and experimental (Jinnai et al., 1993) results. To further investigate the morphology development within the bulk, the dimensionless structure factors were developed in the vertical direction (i.e., along the y*-axis) of the bulk domain. This is due to the existence of a temperature gradient within the horizontal direction (i.e., along the x*-axis) of the bulk. Fig. 4 is a logarithmic plot of the structure factor S*(km*, t*) as a function of dimensionless time t* corresponding to the case shown in Fig. 2(c). S*(km*, t*) is the structure factor S*(k*, t*) evaluated at the

Fig. 3. Typical evolution of the dimensionless structure factor for the case shown in Fig. 2(b) for an off-critical quench from the one-phase region into the two-phase region at different dimensionless times. The parameter values are: c0* ¼0.3, h1 ¼ 0.5, g ¼  0.5, D ¼4  105, T1* ¼0.24 and T2* ¼ 0.20.

dimensionless wavenumber km*, which is where the maximum of S*(k*, t*) is located at time t*. It is noticed that S*(km*, t*) increases exponentially in the early stage but slows down as the phase separation enters the intermediate stage. As well, Fig. 4 shows the transition time from the early to intermediate stages of the SD mechanism within the bulk occurs at a later time for shallower quenches at T* ¼ 0.23, 0.2357 and 0.2525, and at earlier time for deeper quenches at T* ¼0.2075 and 0.2225. The obtained diagram is consistent with published experimental (Krausch et al., 1993) and numerical (Chan and Rey, 1995a) results. Fig. 5 shows the time evolution of S*(km*, t*) on a log–log scale for the case of T1* ¼ 0.2525 in Fig. 4. The dimensionless transition time tnt from the early to intermediate stages of spinodal decomposition is determined from the intersection of the tangent lines drawn over early stage and intermediate stages (Chan and Rey, 1995a), which is tt* ¼0.5480 in this particular case. Typically S*

890

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

(km*, t*) has power-law dependence on time as:
n  α S n km ; t n  t n

ð24Þ

where α is the growth exponent. The Lifshitz–Slyozov (LS) behavior is the only growth law anticipated for the phase separation of polymer systems when there is dominant bulk diffusion. The log– log plot of S*(km*, t*) versus t* in Fig. 5 reveals a growth exponent of 0.31, which is determined from the slope of the tangent line for the intermediate stage of SD in the bulk. The result is consistent with published numerical (Henderson and Clarke, 2005) and experimental (Geoghegan et al., 2000) work, and validates the mathematical model consisting of the Cahn–Hilliard theory. To calculate the thickness layer of the polymer that wets the surface, two considerations are taken into account. Firstly, only the cases in which the polymer that completely wets the surface are considered in the process. For this selection the following procedure is applied. Initially, the concentration of component 1 on the hypothetical line (the line adjacent to the surface that is parallel and close to the wetting surface) is calculated for each node. It was observed that for complete wetting the difference between the highest and the lowest concentration along this line is less than 15%. Therefore, this criterion is used to define the surface as complete wetting along the hypothetical line. Secondly, the

Fig. 4. Maximum structure factor S*(km*, t*) as a function of dimensionless time t* for the simulation shown in Fig. 2. This curve is typical of SD, since there is an exponential growth at first but then slows down. The transition time from the early to intermediate stages of SD occurs at later times for shallower quenches at T* ¼0.23, 0.2357 and 0.2525.

Fig. 5. Typical time evolution of the logarithmic structure factor S*(km*, t*) as a function of dimensionless logarithmic time corresponding to the case presented in Fig. 4. The intersection of the tangent lines drawn over early stage and intermediate stage of SD represents the transition time (t* ¼ 0.5480). The slope of the tangent line for the intermediate stage is calculated to be 0.31 which is consistent with the Lifshitz–Slyozov law.

distance between the wetting surface and the first point in the domain in the direction perpendicular to the surface that had the closest value to climit was measured. This was followed by calculating the average of these values. For calculating the climit, Eq. (25) is used based on the difference between the concentration of component 1 at the surface, which represents the highest value of concentration for the spinodal wave, and its concentration at the first minima of the spinodal wave: climit ¼

cmax þ cmin 2

ð25Þ

The average of these thickness layers in the 256 points along the surface represents the thickness of wetting layer at that specific time which varies with time:
 m z tn  tn ð26Þ where z and m represent the thickness of the wetting layer and growth rate, respectively. The growth rate m can be obtained from the slope of a plot of log z(t*) versus log t*. Fig. 6 shows how the surface enrichment layer growth typically changes with time during the early and intermediate stages of SD at the surface for a shallow quench case with parameters: c0*¼ 0.5, g¼  0.5, h1 ¼0.5, D¼ 4  105, T1*¼0.30 and T2*¼ 0.20. It is observed that in the early stage of SD, the growth rate of surface enrichment layer will increase rapidly (where m¼0.5). As the system enters into the intermediate stage, the growth rate decreases (where m¼0.13). This can be explained by the fact that during phase separation at the surface for a shallow quench, due to the lack of competition between the surface and the bulk, the initially formed droplets adjacent to the surface will grow faster in the absence of bulk phase separation commencement. The preferred polymer would then be attracted to the surface till the phase separation is initiated in the bulk. During the intermediate stage, because of the bulk phase separation, the surface enrichment layer growth becomes slower and tends to break up (i.e., partial wetting starts). This behavior has been observed experimentally where there is a transition from the complete wetting to partial wetting of the polymer at the surface for a shallow quench (Geoghegan et al., 2000). Enrichment of the surface will continue during intermediate stage. Fig. 7 shows the morphological changes over time from complete wetting to partial wetting for variations from shallow quenches to deep quenches of a polymer blend where the phase separation mechanism is governed by off-critical quench conditions. This trend has been observed for different values of c, g, and h1 (Yan and Xie, 2006b). Fig. 7(a) shows partial wetting of the layer attached to the surface at t*¼ 8.48  10  1, which is earlier than at the temperatures in (b) and (c) due to its deeper quench temperature. Higher values of

Fig. 6. Typical surface enrichment growth rate m at early (m¼ 0.5) and intermediate (m¼ 0.13) stages of the phase separation for a shallow quench in a shortrange surface potential case with the following parameter values: c0* ¼0.5, D¼ 8  105, h1 ¼0.5, g ¼  0.5, T1* ¼0.30 and T2* ¼ 0.20.

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

891

Fig. 7. Morphology changes from complete wetting to partial wetting for various off-critical quench depths in a short-range surface potential case. The parameter values are: c0* ¼0.4, D ¼ 4  105, h1 ¼ 0.5, g ¼  0.5, and T2* ¼ 0.20. The temperature T1* varies from (a) 0.20 for a deep quench to (b) 0.22 and (c) 0.24 for shallow quenches.

temperature have a lower quench depth. This transition from the complete wetted surface to the partial one will happen at earlier time if the quench depth is increased due to the acceleration of the phase separation in deeper quenches. As a result, the deeper the quench, the quicker the partial wetting will appear at the surface. Fig. 8 shows this phenomenon for a typical off-critical quench, where only the quench temperature at the surface was varied. It can be seen from this figure that in the range of quench depth ε studied, the transition time decreases two orders of magnitude as quench depth increases by one order of magnitude. The quench depth is determined from:

ε¼

χ  χs χs

where

ð27Þ

χs is the interaction parameter at the spinodal temperature.

3.2. Effect of different diffusion coefficients Fig. 9 presents the morphology formation for the following three dimensionless diffusion coefficient values: D¼3  105, 5  105 and 6  105. The temperature gradient remained constant with T1*¼0.25 and T2*¼0.20, and h1 ¼1 and g¼  0.5. The initial concentration is c0*¼0.4 for these cases but, in general, all different initial concentrations follow the same pattern. It can be observed that the lower diffusion value D¼3  105 has less driving force for phase separation. This then requires more time to reach the intermediate stage of the phase separation process compared to the higher value of the diffusion coefficient D¼6  105. Fig. 9 also shows that the phase separation is already in the intermediate stage of SD at t*¼ 9.75  10  3 for D¼6  105; however, for the lower diffusion coefficient of D¼5  105, the process is just starting to enter the

Fig. 8. The transition time from complete wetting to partial wetting for different deep quench depths at the surface in a short-range surface potential case. The parameter values are c0* ¼0.4, D ¼4  105, h1 ¼ 0.5, g ¼  0.5 and T2* ¼ 0.20. The temperature T1* varies from 0.20, for a deep quench to 0.12 for the deepest quench performed.

intermediate stage, and for D¼3  105, the phase separation has not even started yet. Consequently, the driving force for phase separation increases as D increases in the range of D studied. This trend can be predicted in the model Eq. (16) where the rate of change of c* is proportional to D. It is also observed in the morphology in Fig. 9 that as the diffusion coefficient increases, the surface enrichment rate decreases. This is due to the high rate of phase separation within the bulk, which leads to the starvation of the surface in attracting the preferred component. A higher diffusion coefficient value leads to higher rate of phase separation within the bulk, which creates finer morphology structures

892

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

Fig. 9. Morphology formation for three diffusion coefficient values in a short-range surface potential case where (a) D¼ 3  105, (b) 5  105 and (c) 6  105. The temperature gradient remained constant at T1* ¼ 0.25 and T2* ¼0.20. The parameter values are: c0* ¼ 0.4 (off-critical quench), h1 ¼1 and g ¼  0.5. As D increases, so does the driving force for SD phase separation.

in both off-critical (i.e., droplet-type morphology) and critical (i.e., interconnected structure) cases. Lastly, the diffusion coefficient D controls the rate of phase separation, however, it does not affect the type of morphology formed; i.e., an interconnected structure (droplet-type morphology) forms after a critical (off-critical) quench into the unstable region of the binary phase diagram. The results in Fig. 9 are in agreement with extensive numerical work by Lee et al. (2002, 2003, 2004), Hong and Chan (2010), and Chan and Rey (1995a,b).

separation on the deeper quench side accelerates the phase separation process for the whole domain and make the morphology formation process start at earlier times even for the area with the shallower quench depth near the surface. This will lead to the higher attraction of the preferred polymer component to the surface as T2* decreases, which will result in an increase of the rate of surface enrichment and thickness of the wetting layer. To the best of our knowledge, this effect of temperature gradient on surface enrichment, where T1* is constant and T1*4T2*, has been analyzed for the first time in this article.

3.3. Effect of temperature gradient on surface enrichment 3.4. Effect of surface potential on surface enrichment In this section, the effect of the temperature gradient on the morphology formation is studied. Fig. 10 shows a typical morphology evolution for an off-critical quench; the parameter values for the cases shown are c0*¼ 0.4, D¼4  105, h1 ¼0.5, g¼  0.5 and T1*¼0.25. The dimensionless temperature T2* is varied from 0.10 to 0.20. It is important to note that phase separation through the TIPS method is always induced by a temperature quench whether with or without a temperature gradient within the bulk. In Fig. 10, the deepest quench where T2*¼0.10 allowed more polymer component to be attracted to the surface at any specific time. The thickness of the wetting layer for the deeper quench depth (Fig. 10a) is approximately double that obtained for the shallower quench depth (Fig. 10c). This behavioral trend was observed for the different values of g, h1 and D studied in this paper. Since there is a temperature gradient along the polymer blend, phase separation will start at different times. Lee et al. (2004) discussed how a linear temperature gradient affects the phase separation morphology after a quench into the unstable region. The inclusion of a linear temperature in the model creates a linear differential quench temperature profile along the sample, which in turn will create differential stages of phase separation along the sample. It will start first at the deeper quench end (i.e., at T2*), allowing it to progress through the early stage and possibly the intermediate stage of SD before phase separation starts at the shallow quench end at the surface with preferred attraction for one polymer component (i.e., at T1*). The initiation of the phase

Fig. 11 shows the effect of short-range surface potential on the morphology development and evolution during the phase separation for an off-critical. There is a small temperature gradient along the domain, where T1* and T2* are 0.20 and 0.18, respectively. The following three cases of surface potential h1 strength were studied: (a) 0.5, (b) 2.0 and (c) 6. The higher values of h1 result in delaying transition time from complete wetting to partial wetting. When the surface potential is high, more favored polymer would be attracted to the surface during the early stage of phase separation leading to the faster formation of surface complete wetting. As the system approaches the intermediate stage, where phase separation is governed mostly through the bulk, the partial wetting mechanism of the completely wetted layer would be delayed due to the higher surface attraction strength (i.e., h1). In Fig. 11(a), the phase separation within the bulk domain initiated the surface partial wetting at t*¼2.7  10  1, which is earlier than the cases in (b) and (c) where the surface potential is higher. The domain in Fig. 11(b) is about to start its partial wetting at t*¼4.42  10  1, while in (c) there is no sign of partial wetting yet. This complies well with the fact that higher values of surface attraction potential resist against partial wetting and tends to enrich more of the preferred component to the surface. This causes full wetting of the layer closes to the surface. The elongated domain morphology is the result of the system concentration being close to its critical value, where the phase separated

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

893

Fig. 10. Effect of temperature gradient on surface enrichment growth rate where phase separation is governed by off-critical quench condition in a short-range surface potential case. The parameter values are c0* ¼ 0.4, D¼ 4  105, h1 ¼ 0.5, g¼  0.5, T1* ¼0.25, and (a) T2* ¼0.1, (b) T2* ¼0.15 and (c) T2* ¼ 0.2.

Fig. 11. Typical effect of the surface potential on the surface formation where phase separation is governed by off-critical quench condition in a short-range surface potential case. In this case, c0* ¼ 0.4, g ¼  0.5, D¼ 4  105, T1* ¼0.2 and T2* ¼ 0.18, and the value of h1 varies from (a) 0.5 to (b) 2.0 and (c) 6. Higher h1 values delays the transition from complete wetting to partial wetting.

morphology tend to form an interconnected structure (Cahn, 1965; Chan and Rey, 1995a,b). The strip next to the surface wetting layer in Fig. 11(b) and (c) at t*¼ 7.8  10  2 is a result of the effect of a stronger surface potential creating a SD wave normal to the surface into the bulk during the early stage of SD. As time passes, the phase separation in the bulk overcomes the phase separation initiated by the surface leading to the rupture of the SD wave and destruction of the strip (Jones et al., 1991; Krausch et al., 1993).

Fig. 12(a) presents the typical logarithmic structure factor .S*(km*, t*) as a function of logarithmic dimensionless time corresponding to the phase separated morphology formed at x* ¼ 0.125 in Fig. 11. As the surface potential increased from 0.5 to 6, the transition time increased nonlinearly from approximately 4–6. Higher values of surface potential resulted in faster transition time from the early stage to the intermediate stage of SD within the bulk. This is due to the higher attraction of the

894

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

Fig. 13. Typical logarithmic structure factor S*(km*, t*) as a function of dimensionless time t*, presenting the effect of the surface potential on the transition time where phase separation is governed by a critical quench in a short-range surface potential case.

Fig. 12. (a): Typical logarithmic structure factor S*(k*m, t*) as a function of dimensionless time in the vicinity of the surface (x* ¼0.125) corresponding to the case presented in Fig. 11. Higher h1 values resulted in faster transition time from early stage to intermediate stage within the bulk. (b): Typical logarithmic structure factor S*(km*, t*) as a function of dimensionless time t* in the bulk at x* ¼0.5 corresponding to the case presented in Fig. 11. Different h1 values have no impact on the transition time from early to intermediate stages within the bulk for distances farther away from the surface; this is an evidence of a short-range surface potential case.

preferred polymer to the surface. Fig. 12(b) displays the typical logarithmic structure factor S*(km*, t*) as a function of logarithmic dimensionless time transition time at x* ¼0.5 in Fig. 11. The transition time for all three surface potential is the same, which means that the short range surface potential does not interfere with the phase separation process far away from the surface. Fig. 13 shows the typical logarithmic structure factor S*(km*, t*) as a function of logarithmic dimensionless time for a critical quench with the following parameter values c0*¼0.5, g¼ 0.5, D¼ 8  105, h1 ¼2.0, and T1*¼ T2*¼ 0.2. The morphology development occurs fastest near the vicinity of the surface wall (i.e., at x*¼0.125) like in the off-critical quench case above. Thus, morphology development occurs fastest near the surface wall regardless of initial concentration. This figure shows that the time evolution of S*(km*, t*) nearly overlap each other for distances x*Z0.25, which indicates that the phase separation process is dominated by bulk effect in this region, whereas the surface effect dominates for regions 0rx*r0.125. This finding is also consistent with the result shown in Fig. 12 for an offcritical quench. Lastly, the oscillations in the curves shown in Figs. 12 and 13 are artefacts of the finite mesh points used in the calculation of the structure factor.

4. Conclusions In this study, a model composed of the nonlinear CH theory for phase separation and the FHdG free energy theory were used to study numerically the SDPS phenomena of a binary polymer blend quenched into the unstable region of its symmetric phase diagram. The model was solved within a square geometry with a temperature gradient in the direction normal to the surface with a short range surface potential for the preferred polymer component. The numerical results are consistent with past numerical work (Puri

and Binder, 1992; Henderson and Clarke, 2005; Yan and Xie, 2006a,b; Lee et al., 2002, 2003, 2004; Chan and Rey, 1995a,b). The structure factor analysis showed a faster exponential growth at the early stage of phase separation and a slower growth rate at the intermediate stage with a slope of 0.31 through the bulk, which is consistent with the LS growth law. The investigation of surface enrichment rate at the surface wall demonstrated faster growth rate at the early stage with the slope of 0.5. This rate of growth became slower at the intermediate stage with a slope of 0.13 near the surface, which is consistent with published experimental observations. To examine the role of quench depths in the phase separation morphology, the following results were obtained for deep quenches: (a) faster transition time from complete wetting to partial wetting of the surface, (b) higher rate of morphology development within the bulk that contributed to faster transition time from early to intermediate stage, (c) lower surface enrichment due to losing the competition to the bulk in attracting favorable polymer to the surface, and (d) smaller droplets and finer morphology formation. Higher diffusion coefficients led to the increase of driving force for phase separation and consequently faster morphology development. The influence of various temperature gradient ΔT* values on the surface enrichment rate with constant temperature T1* at the surface and various temperature T2* for the opposite end was studied for the first time. The results showed that the thickness of the wetting layer increased by rise of ΔT* value where the side with T2* temperature goes under deep quench. This feature is due to the phase separation starting at earlier stage at the wetting surface region since the other part of the sample that is in more advanced stage could stimulate the initiation of phase separation earlier in all domains. An analysis of the structure factor evolution showed earlier transition time for higher values of h1 for phase separation within the bulk but close to the surface. However, there was no difference observed for transition time within the bulk at distances farther away from the surface at x* Z0.25. As surface potential increased, spinodal wave became more visible in the bulk and the transition time from complete wetting to partial wetting occurred at a later time on the surface. The numerical solutions and calculated morphologies replicate frequently reported experimental observations and numerical work.

Acknowledgments The authors gratefully acknowledge the financial support from the Natural Sciences and Engineering Research Council of Canada and Ryerson University.

M. Tabatabaieyazdi et al. / Chemical Engineering Science 137 (2015) 884–895

References Brochard, F., Jouffroy, J., Levinson, P., 1983. Macromolecules 16, 1638. Brown, G., Chakrabarti, A., 1992. Phys. Rev. A 46, 4829. Cahn, J.W., 1965. J. Chem. Phys. 42, 93. Chan, P.K., Rey, A.D., 1995a. Macromol. Theory Simul. 4, 873. Chan, P.K., Rey, A.D., 1995b. Comput. Mater. Sci. 3, 377. Dalmoro, A., Barba, A.A., Lamberti, G., Grassi, M., d’Amore, M., 2012. Adv. Polym. Technol. 31, 219. deGennes, P.G., 1980. J. Chem. Phys. 72, 4756. Evans, D.F., Wennerstrom, H., 1999. The Colloidal Domain: Where Physics, Chemistry, Biology and Technology Meet, second ed. Wiley-VCH, New York. Flory, P.J., 1953. Principles of Polymer Chemistry. Cornell University Press, Ithaca. Geoghegan, M., Jones, R.A.L, Clough, A.S., 1995. J. Chem. Phys. 103, 2719. Geoghegan, M., Ermer, H., Jungst, G., Krausch, G., Brenn, R., 2000. Phys. Rev. E 62, 940. Gray, J.J., 2004. Curr. Opin. Struct. Biol. 14, 110. Henderson, I.C., Clarke, N., 2005. Macromol. Theory Simul. 14, 435. Hong, S., Chan, P.K., 2010. Model. Simul. Mater. Sci. Eng. 18, 025013. Jinnai, H., Hasegawa, H., Hashimoto, T., Briber, R.M., Han, C.C., 1993. Macromolecules 26, 182. Jinnai, H., Hasegawa, H., Hashimoto, T., Han, C.C., 1993. J. Chem. Phys. 99, 8154.

895

Jones, R.A.L., Norton, L.J., Kramer, E.J., Bates, F.S., Wiltzius, P., 1991. Phys. Rev. Lett. 66, 1326. Klein, J., 1978. Nature 271, 143. Kramer, E.J., Green, P.F., Palmstorm, C.J., 1984. Polymer 25, 473. Krausch, G., Dai, C.-A., Kramer, E.J., Bates, F.S., 1993. Phys. Rev. Lett. 71, 3669. Krausch, G., Kramer, E.J., Bates, F.S., Marko, J.F., Brown, G., Chakrabarti, A., 1994. Macromolecules 27, 6768. Kyu, T., Mustafa, M., Yang, J.-C., Kim, J.Y., Palffy-Muhoray, P., 1992. Stud. Polym. Sci. 11, 245. Lee, K.-W.D., Chan, P.K., Feng, X., 2002. Macromol. Theory Simul. 11, 996. Lee, K.-W.D., Chan, P.K., Feng, X., 2003. Macromol. Theory and Simul 12, 413. Lee, K.-W.D., Chan, P.K., Feng, X., 2004. Chem. Eng. Sci. 59, 1491. Mishra, K., Hashmi, S.A., Rai, D.K., 2013. J. Solid State Electrochem. 17, 785. Norde, W., 2003. Colloids and Interfaces in Life Sciences. Marcel Dekker, New York. Puri, S., Binder, K., 1992. Phys. Rev. A 46, R4487. Ramanathan, T., Abdala, A.A., Stankovich, S., Dikin, D.A., Herrera-Alonso, M., Piner, R.D., Adamson, D.H., Schniepp, H.C., Chen, X., Ruoff, R.S., 2008. Nat. Nanotechnol. 3, 327. Tanaka, H., Yokokawa, T., Abe, H., Hayashi, T., Nishi, T., 1990. Phys. Rev. Lett. 65, 3136. 〈https://computation.llnl.gov/casc/sundials〉. Yan, L.-T., Xie, X.-M., 2006a. Macromol. Theory Simul. 15, 226. Yan, L.-T., Xie, X.-M., 2006b. Polymer 47, 6472. 〈http://www.mathworks.com/products/matlab/〉. Yaseen, M., Salacinski, H.J., Seifalian, A.M., Lu, J.R., 2008. Biomed. Mater. 3, 034123.