Computational Materials Science 38 (2006) 325–339 www.elsevier.com/locate/commatsci
Computational thermodynamics of multiphase polymer–liquid crystal materials Susanta K. Das, Alejandro D. Rey
*
Department of Chemical Engineering, McGill University, 3610 University Street, Montreal, Que., Canada H3A 2B2 Received 28 January 2005; accepted 10 October 2005
Abstract De-mixing in polymer–nematic liquid crystal mixtures is simulated using a previous model that includes phase separation–phase ordering–texturing processes. It is shown that the activation of phase separation–phase ordering–texturing is associated with four internal length scales, whose relative magnitudes gives rise to three regimes of distinct morphological and textural features, as follows: (i) Phase separation and formation of few and large polymer droplets, leads to a defect lattice in the nematic matrix through the prevailing strong anchoring at the polymer–liquid crystal interface; (ii) When the number of polymer droplets increases and their size decreases, a random polymer droplet morphology emerges; due to weak anchoring conditions the nematic matrix is random and weakly oriented; (iii) When the droplet size is of the order of the nematic ordering length scale, phase ordering is frustrated and small droplet morphology emerges. This paper shows that coupling between phase separation–phase ordering–texturing processes lead to complex de-mixing morphologies with varying degrees of positional and orientational order, thus augmenting ways of controlling material architectures beyond the well-known droplet patterns in off-critical quenches of spinodal decomposition. 2006 Elsevier B.V. All rights reserved. PACS: 61.72.y; 61.72.Bb; 64.90.+b; 64.75.+g; 61.82.Rx; 81.30.t Keywords: Multiphase materials; Polymer–liquid crystal; Morphology; Computer simulation; Thermodynamics; Texture
1. Introduction Liquid crystals are orientational ordered viscoelastic, textured, and anistropic materials, used as electro-optical, sensor, and structural applications [1,2]. As in other organic and inorganic material systems, blending and mixing is used to improve performance and reduce costs. Electro-optical applications as light valves for tunable windows [3,4] usually involve blends of liquid crystals with flexible polymers. These materials usually de-mix, forming complex morphologies that reflect the interaction between thermodynamic phase separation due to the presence of polymer and liquid crystal, phase ordering due to the presence of the molecular
*
Corresponding author. Tel.: +1 514 398 4196; fax: +1 514 398 2753. E-mail address:
[email protected] (A.D. Rey).
0927-0256/$ - see front matter 2006 Elsevier B.V. All rights reserved. doi:10.1016/j.commatsci.2005.10.009
order in liquid crystals, and texturing due to the presence of orientation in liquid crystals. Experiments [2–5], theory [1– 4], and simulation [6–9], have found that the development of patterns during de-mixing driven by thermal quenches involves interplay of phase separation and phase ordering. It has been found that phase ordering can drive the demixing through spinodal decomposition driven by phase ordering fluctuations [10–12]. Since liquid crystals are characterized by orientation, in most cases these materials exhibit textures, i.e., defects and orientation gradients due to confinement effects. If we consider morphology of polymer droplets embedded in a liquid crystal matrix, the curvature droplets will exert a distortion that can translate into orientation gradients and defect generation. Thus, de-mixing involving polymers and liquid crystals have to take into account three distinct but interacting processes: (1) phase separation, (2) phase ordering, and (3) texturing. The role
326
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
and impact of texturing on pattern formation is not fully understood, and is the main motivation of this work. The motivation of this work can be made more specific by consideration of a few standard equations. The early stages of spinodal decomposition in isotropic polymer solutions is frequently described using the Flory–Huggins– Cahn–Hilliard (FHCH) model [13–15]. In this model, the free energy functional F reads: Z W ~ Lu 2 F ¼ f H/ ð/Þ þ ðr/ðr; tÞÞ dV ; ð1Þ 2 2 V where W > 0 is the barrier height of the double well dimensionless homogeneous free energy f~ H/ ð/Þ, / is the polymer volume fraction, and Lu > 0 is a gradient energy coefficient. r, t, V represents the position vector, time and volume of the medium respectively. The diffusive interface separating the polymer-rich and solvent-rich phases is characterized by an interface thickness ‘/ of the order of L/ ‘/ ¼ pffiffiffiffiffi . ð2Þ W The interface tension cu is of the order of pffiffiffiffiffi ð3Þ c/ ¼ W ‘/ ¼ L/ W . In the early stages of phase separation and for off-critical quenches of phase separation, the FHCH model in 2D predicts the emergence of a random droplet morphology described by X /ðr; tÞ /ðr; 0Þ ¼ Ak ðk; tÞeikr ; ð4Þ k
where k is wave vector and t is time. In the beginning of the intermediate regime the characteristic droplet radius is R ¼ aðc/ Þm ¼ aW m ‘m/ ¼ aLm/ W m=2 ;
ð5Þ
where m > 0 and a is a real constant; this expression implies that the size of the predicted droplet in the intermediate regime is an increasing function of the interfacial tension. When the solvent is mesogenic, the separating phases can be isotropic or liquid crystalline [16], and the phase separation–phase ordering couplings drives new processes, as shown in this paper. For uniaxial nematic liquid crystalline phases, the order is described by the tensor order parameter Q = S(nn d/3) [17] where S is the scalar order parameter, n the average molecular orientation, and d is the unit tensor. Hence for polymer–nematic liquid crystalline mixtures the energy functional is augmented by the following new terms: Z W ~ L/ 2 F ¼ ðf H/ ðQ; /; T Þ þ f~ HQ ðQ; /; T ÞÞ þ ðr/Þ 2 2 V L/Q LQ 2 þ ðruÞðrQÞ þ ðrQÞ dV . ð6Þ 2 2 In this mixture the interfacial thickness between isotropic and nematic phases become: LQ ‘Q ¼ pffiffiffiffiffi . W
ð7Þ
The coupling energy term describes, among other effects, the ability of concentration gradients $/ to affect the orientation: L/Q ðr/Þ ðr QÞ 2 L/Q ½Sðr/Þ ðr nnÞ þ ðr/Þ ðrSÞ ðnn I=3Þ. ¼ 2
f/Q ¼
ð8Þ The ability of a concentration gradient (D/)/‘/ to orient the director is given by the anchoring energy can which according to Eq. (8) is estimated by can ¼
SðD/ÞL/Q . ‘/
ð9Þ
The internal length scale associated with the anchoring energy can is known as the extrapolation length ‘ext [17] and in the present case it is given by LQ LQ LQ LQ pffiffiffiffiffi . ‘ext ¼ ¼ ‘/ ¼ ð10Þ can SðD/ÞL/Q SðD/ÞL/Q W According to the magnitudes of ‘ext and R, two wellknown anchoring regimes arise [17]: (i) ‘ext R, weak anchoring, where the drop has insignificant effect on the nematic matrix, and (ii) ‘ext R, strong anchoring, where the presence of the drop creates a texture in the nematic matrix. In addition according to the sign of L/Q two interfacial orientation arise [18,19]: (a) homeotropic anchoring L/Q > 0, n ? $/, (b) planar anchoring L/Q > 0, nk$/. Based on the above observations we conclude that the phase separation–phase ordering–texturing process is intrinsically multiscale, and is characterized by the following set of length scales: 9 8 > > > > = < ; X¼ ‘/ ; ‘Q ; ‘ext ; R |fflffl{zfflffl} > > |{z} |{z} > > ; : phase separation phase ordering texturing ð11Þ where the underline text indicates the process associated with each scale. Different length scale orderings will give rise to different regimes. For example if R ‘ext > ‘/ ‘Q ;
ð12Þ
we expect that a few large drop embedded in a textured nematic matrix, and where the defect–defect interaction leads to positional order, as observed experimentally [20– 23]. This paper seeks to establish the role of the four length scales on the morphology of the two-phase isotropic– nematic blend and on the texture of the nematic matrix. The emphasis is on establishing fundamental mechanisms heretofore unidentified in the literature, and not to simulate a specific chemical system, or reproduce experimental data.
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
This paper is an extension of an ongoing research program [16,18,24–26] in computational polymer–liquid crystal thermodynamics and rheology. The particular contribution of this paper is the study of the role of surface tension (i.e., interface thickness) on texturing and pattern formation. The theoretical model used in this paper was formulated and analyzed previously [8,27,28]. In this paper we adopt this previous formulation [8,27,28], due to its computational tractability and accuracy. This paper is organized as follows. The first part presents the dynamical model and computational methods. The second part presents the results and discussions. Finally, conclusions are given. 2. Dynamical model and computational methods In this section we present a brief discussion of the governing equations, as given by [8]. A more lengthy discussion of the equations used in this paper is given in [18]. 2.1. Free energy To describe ordering and orientation of liquid crystal molecules we use the tensor order parameter, defined by [17,18,29] d Q¼ f ðuÞ uu dA; 3 Z2 Z
ð13Þ
where u is a unit vector parallel to the liquid crystal molecules, d is a unit tensor, f(u) is an orientation distribution function, Z2 denotes the unit sphere, and dA is the differential area. In uniaxial state, the tensor order parameter can be written as Q = S(nn I/3), where S, as discussed above, is the molecular alignment along the uniaxial director n, and is given by S ¼ 32 ðn Q nÞ. If there are no external fields, then the free energy functional of the system will consists of the bulk free energy and a non-local free energy that penalizes gradients in composition and tensor order parameter, which according to [8] reads: h
g
f ð/; QÞ ¼ bDFmix =N T ¼ f ð/; QÞ þ f ð/; QÞ;
ð14Þ
f h ð/; QÞ
/ ð1 /Þ ln / þ lnð1 /Þ þ v/ð1 /Þ nI nA 3 C0 ð1 /Þ I0 2 þ ln ð1 /Þ Q : Q ; 4 nA nA 2
¼ ck B T
ð15Þ
g
f ð/; Qij Þ ¼
L/ LQ1 2 2 ðr/Þ þ L/Q ðoi /Þðoj Qij Þ þ ðok Qij Þ 2 2 þ
LQ2 ðoi Qik Þðoj Qjk Þ; 2
ð16Þ
327
where I0 becomes Z 2p Z p 3 d I0 ¼ exp C0 ð1 /ÞQ : rr sin2 hdhdx; 2 3 0 0 ð17Þ C0 ¼ ðva þ 5=4ÞnA ;
b ¼ 1=k B T .
ð18a; bÞ
Here, f is the total free energy density, f h is the homogeneous energy density, and f g is the gradient energy density [17]. T is the absolute temperature, kB is the Boltzman constant, and / denotes the polymer volume fraction. v is the Flory–Huggin’s (FH) interaction parameter [30], the terms va and (5/4)nA in C0 indicate the orientation-dependent attractive (Maier–Saupe; MS) interactions between the mesogens [31] and excluded volume interactions [32] between mesogenic molecules respectively. nI and nA are the number of segments on the isotropic (monomer or polymer) component and the number of segments (axial ratios) on the mesogen, respectively. c denotes per unit volume of a cell or segment. LQ1 and LQ2 are the Frank elastic constants [33]. r is a unit vector in spherical coordinates. If a = va/v, where a represents the relative strength of interactions, we can find the nematic–isotropic transition (NIT) temperature [8,16,18,27]: sNI ð/Þ ¼
anA ð1 /Þ ; 4:55 1:25nA ð1 /Þ
ð19Þ
where the threshold level is set to C0(1 /) = 4.55 [17]. 2.2. Derivation of kinetic equations and numerical simulation Considering the dynamical random phase approximation [19] for the polymer concentration /, the non-equilibrium equation of motion which ensures local conservation of the material can be written as [8,18] o/ df ð/; QÞ of of ¼ c0 r 2 r ¼ c0 r 2 ; ð20Þ ot d/ o/ or/ where the thermodynamic force that drives the flux is given by the gradient of the chemical potential l = df/ d/, and c0 is the mobility, assumed to be constant. While deriving equation of motion for the non-conserved order parameter Q, the local rate of change is assumed to be proportional to the local thermodynamic force (df/dQ). Therefore, the equation of motion for Q can be written as [8,18] ½s ½s oQ df ð/;QÞ of of ¼ n0 r ¼ n0 ; ot dQ oQ orQ
ð21Þ
where n0 is the transport coefficient taken as constant, and [s] denotes symmetric and traceless tensor. Substituting Eqs. (14–18) into Eqs. (20) and (21) yields the following governing system of equations [8,18]:
328
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
h ½s o/ 2 of 4 2 ð22Þ L/ r / L/Q r ðr rQÞ ; ¼ c0 r o/ ot h ½s ½s oQ of 2 L/Q rr/ LQ1 r Q LQ2 rðr QÞ . ¼ n0 oQ ot ð23Þ
In this paper we use LQ1 = LQ2 = LQ, and thus take into account Frank elastic anisotropy [33]. We introduce the following dimensionless quantities as x ~x ¼ ; ‘
y ~y ¼ ; ‘
T Te ¼ ; T
ð24Þ
where ‘ and T* represents system (external) length scale and temperature respectively, the dimensionless kinetic equations (22) and (23), then becomes " #½s ! ~h o/ 2 of 4 2 ~ ~ ~ ~ ~ e ð25Þ ¼ K/ T r r / K/Q r ðr rQÞ ; o~t o/ " #½s oQ of~ h 2 2 ~ ~ ~ ~ e þ K/Q r / þ KQ ðr Q þ rðr QÞÞ ; ¼ N T K/ o~t oQ ð26Þ where the dimensionless numbers are: ~t ¼
n0 L/ c0 t ðdimensionless timeÞ; c0 ‘2
ð27Þ
‘ 2 n0 ðmobility ratioÞ; ð28Þ c0 2 ‘2 ‘ 2 K/ ¼ ¼ ðsystem size/interface thicknessÞ ; ‘/ L/ =ck B T
N¼
ð29Þ L/Q L/ 1 ‘Q ¼ ðnematic interface thickness/ SD/ ‘ext
K/Q ¼
(extrapolation lengthÞSD/Þ;
ð30Þ
LQ KQ ¼ L/ ¼
LQ =ck B T ‘Q ¼ ðnematic interface thickness/ L/ =ck B T ‘/ interface thicknessÞ.
ð31Þ
In this study, the system size (length) scale, ‘, used as a characteristic length scale for non-dimensionalization purpose is a constant and hence the computational domain is a unit square. Depending onpthe ffiffiffiffiffiffi particular emphasis in this paper, we use K/ and 1= K/ ¼ ‘/ =‘ to describe results. The physical significance of changing K/ is a change in interface thickness ‘/ since the system size is considered
to be constant and equal to ‘. A physical and real way to change the interface thickness ‘/ is through addition of a surfactant. The governing equation (25) for the conserved parameter, /, indicates a time-dependent nonlinear parabolic partial differential equation. Whereas, the governing Eq. (26) for Q is a set of five coupled, time-dependent, nonlinear, partial integro-differential equations. For numerical simulation, the governing dimensionless set of partial differential equations (25) and (26) is transformed into a set of ordinary differential equations using second-order centered finite difference approximations for space discretization. Time integration of the system is performed with a stabilized explicit Runge–Kutta–Chebyshev adaptive time integration algorithm [34]. Throughout the numerical computations, mesh and time step independence was established. The computational domain we considered here is the unit square ð0 6 ~x 6 1; 0 6 ~y 6 1Þ and discretized the entire domain into 256 · 256 (65 536) cells. All double integrations were performed using the composite Simson’s rule [35] with sufficient accuracy. The initial conditions were taken in the form of the isotropic uniform mixtures with small random fluctuations and periodic conditions were imposed at the boundary. For simplicity in the results section, due to presentation reasons, we used x (x ¼ ~x) and y (y ¼ ~y ) instead of ~x and ~y respectively. In all calculations we use the following parametric set: N = 1.0, K/Q = 0.2, KQ = 0.1. In this paper, we study incompressible mixtures of flexible polymers made up of nI = 25 monomers and low molar mass rod-like mesogens consisting of nA = 2 monomers, the relative strength of interactions is set to be a = 2.5. 2.3. Characterization of kinetic process Phase ordering, phase separation and texture formation kinetics are analyzed using standard structure factors [9,18]. To characterize the emerging textures, we calculated the temporal evolution of the structure factor for both compositional and orientational fluctuations. The structure factor, W(k, t), is assumed to be proportional to the experimentally observed small angle light, X-ray or neutron scattering intensity I(w, t). W(k, t) is the squared magnitude of the Fourier transform of the concentration fluctuations, jA0(k, t)j2, and orientation fluctuations, jB0(k, t)j2; i.e. 2 NN 1 NN 1 2pi fmk þnk g X X W/ ðk; tÞ ¼ ½/ðr;tÞ /0 e½NN 1 2 ; ð32Þ m¼0 n¼0 2 NN 1 NN 1 2pi fmk þnk g X X WS ðk; tÞ ¼ ½Sðr;tÞe½NN 1 2 ; ð33Þ m¼0 n¼0 where, in our two-dimensional study, k = (k1, k2, 0) and NN is the number of nodes along each coordinate axis. The computed values of k1 and k2 are confined to the range fc 6 ki 6 fc, {ki} = k1, k2, where fc = 1/2D is the Nyquist critical frequency [36] and D is the node spacing of the
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
computational mesh in real space. The maximum structure factors are defined by ¼ W/ ðk /max =2pÞ; Wmax / k /max ðt Þ
Wmax ¼ WS ðk S max =2pÞ; S
ð34a; bÞ
k S max ðt Þ
where and are the time-dependent peak wave numbers for / and S, respectively. 3. Results and discussion Fig. 1 shows the temperature ð Te Þ–concentration (/) equilibrium phase diagram of the system considered in this paper, for nI = 25, nA = 2, and a = 2.5 (Eqs. (14–18)). Details of the procedure for construction and description of the phase diagram for such a system are presented in our previous work [16]. In this paper, we study the morphology following a temperature quench from the isotropic, homogeneous phase into the IN coexistence region below the triple point line. The quench corresponds to the point denoted by A in Fig. 1, lying below the phase separation spinodal curve and above the NIT line, corresponds to a state that is initially unstable with respect to phase ordering but metastable with respect to phase separation. Therefore, it is expected that the initial stage is driven by phase ordering fluctuations, which ignites the phase separation–phase ordering–texturing processes. Next we show and discuss the model predictions by thermal quenches into point A. As in other phase separation studies using the FHCH model [13–15], the evolution is described in the early regime, the intermediate regime, and the late stage. In this paper we focus on the early part of the intermediate regime [13,14]. The intermediate regime, as the name implies, follows the initial linear instability stage, and in this paper the beginning of the intermediate stage was established by calculating the extent of phase separation (see Eq. (47)). All morphological features discussed
329
below (Figs. 3, 5 and 7) are for elapsed times corresponding to the intermediate regime. Further coarsening processes at later times are beyond the scope of this paper. For brevity we only show graphical results that highlight the new aspects of phase separation– phase ordering–texturing processes; results for quenches in other areas of the temperature ð Te Þ–concentration (/) equilibrium phase diagram can be found in our previous work [18]. Fig. 2 shows the equilibrium order parameter S as a function of dimensionless temperature for the binodal concentrations and that of point A. At the temperature corresponding to point A, the value of the order parameter in the equilibrium polymer-rich phase is S = 0, and in the nematic-rich phase is S = 0.653. During phase ordering, the initial order parameter evolves towards the binodal values. 3.1. The three morphological–texture regimes As mentioned above, in this paper we study the role of K/ on the phase ordering–phase separation–texturing processes. By changing the magnitude of K/ we alter the length scale ordering (see Eq. (29)) and hence the morphology and texture of the blend. As mentioned above in this paper we seek to find new structures from thermodynamic instabilities. In practice K/ can be varied using geometry (i.e. change system size), material properties (i.e. change gradient
1.2 left binodal concentration quenching point concentration right binodal concentration
1.0
0.8
2.5 S
Phase separation
I I+ I
2.4
0.6
spinodal
N+I
2.3 Temperature, T
0.4
2.2 2.1
Phase ordering spinodal
0.2
A α = 2.5
2.0
N+I
nI = 25 N+I
NIT line
1.9
0 0
nA = 2
0.2
0.4
0.6
0.8
1.0
1.2
KbT/Ua (χa-1) 1.8 0.0
0.2
0.4
0.6
0.8
1.0
Polymer concentration, φ Fig. 1. Temperature ð Te Þ–concentration (/) equilibrium phase diagram of the system considered in this paper, for nI = 25, nA = 2, and a = 2.5 (Eqs. (14–18)). Adapted from [16] with permission of Elsevier Publishers.
Fig. 2. Equilibrium order parameter S as a function of dimensionless temperature for the binodal concentrations and that of point A. At the temperature corresponding to point A, the value of the order parameter in the equilibrium polymer-rich phase is S = 0, and in the nematic-rich phase is S = 0.653. During phase ordering, the initial order parameter evolves towards the binodal values.
330
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
energy parameters using chemistry or surfactants), and operating conditions (i.e. change temperature). According to the magnitude of K/ the following three regimes and length scale ordering are observed by solving Eqs. (25) and (26). 3.1.1. Colloidal crystal–defect lattice regime: K/ < 103 The main feature of this regime is the presence of few nearly spherical polymer droplets embedded in a nematic matrix containing a defect network (see Fig. 3); the polymer droplets are located in a square lattice, and due to the positional order of the droplets we refer to this regime as colloidal crystal. In this large droplet regime the size of the droplets is of the order of the system size. The size of the polymer droplets as a function of K/ is discussed in Fig. 11. The dimene is given by the following equation sionless droplet radius R (see Fig. 11 below): n e ¼ a ‘/ ; R ‘
ð35Þ
where a 25, n 6/5. The size of the drops is of the system size. Since KQ = 0.1, then using Eq. (31) we find ‘/ = 10‘Q. Using these results in Eq. (30) we find ! e 1=n R 1 ‘/ ‘ext ¼ ¼ 0:2 ! 1=n 10SD/ ‘ext ‘ 10SD/ðaÞ e 1=n ¼ 10SD/ðaÞ1=n ‘ext . !R ð36Þ ‘ Using characteristic values observed in the simulations, S = 0.5, D/ = 0.28, we find that the characteristic length scale ordering is R ‘ext > ‘/ ‘Q
ð37Þ
indicating strong anchoring conditions and hence strong texturing; see Fig. 11 for exact calculation of length scales. Fig. 3a shows the gray scale plot of the concentration profile /ð~x; ~y Þ at ~t ¼ 77:92 1005 for K/ = 103; black (white) represents the isotropic (nematic) phase. At the early stages (not shown here), the system forms a bicontinuous
Fig. 3. (a) Gray scale plot of the concentration profile /ð~x; ~y Þ at ~t ¼ 77:92 1005 for K/ = 103; black (white) represents the isotropic (nematic) phase. At the early stages (not shown here), the system forms a bicontinuous network that rapidly breaks up into polymer circular drops. (b) Corresponding gray scale plot of the scalar order parameter S; black (white) denotes nematic (isotropic) order. (c) Nematic texture; the arrows represent the local nematic director, and defects are marked with small solid circles.
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
ordering is described in concentration-order parameter fields of the form: /ðx; yÞ ¼ /o ðx; yÞ cosðkxÞ cosðkyÞ; Sðx; yÞ ¼ S o ðx; yÞ cosðkxÞ cosðkyÞ;
700
3000 *
t = 77.94
600
*
t = 64
500
*
t = 54
400
*
t = 44
300
*
t = 34
200 *
t = 24
100 0 0
2
4
6
8
a
*
t = 54
1500
*
t = 44
1000
t = 34
500
t = 24
*
*
2
4
6
8
10
*
b
k /2π
k /2π 35
ψ ψ
2000
30
max s
* max
max
Peak wave number
2500
ψsmax/ψφmax
*
t = 64
2000
0 0
3000
φ
1500 1000 500
c
*
t = 77.94
2500
10
*
0 0
ð38Þ
where (/o, So) are the amplitudes, and k is the wave vector of the square cell. Since R > ‘ext, the texturing process is due to the coupling between concentration and tensor order parameter gradient. The texture selection process is coupled the positional ordering process of the polymer drops and by the almost circular shape of the drops. Fig. 3c shows the texture formation; for brevity the initially random textures are not shown here. The figure shows the embedded polymer droplets in a highly textured nematic matrix. The texture is characterized by a defect lattice consisting of singular wedge disclination lines of strength s = 1/2. Due to the emergence of polymer droplets, inducing normal orientation of the director along the polymer–liquid crystal interface, disclination lines of strength s = 1/2 nucleate. The charge of the defects s = 1/2 is due to the normal director anchoring at the droplets; for tangential orientation defects of charge s = +1/2 would nucleate. As elongated fibrils break into droplets, more s = 1/2 defects nucleate, finally forming
Orientational structure factor, ψs
Compositional structure factor, ψφ
network that rapidly breaks up into polymer circular drops (see Fig. 3a). During the phase separation the isotropic droplets grow in the nematic matrix. The growth of polymer droplets in the nematic matrix is governed by the elastic forces from the deformation of the director field of the nematic matrix. Fig. 3a shows that the polymer droplets have positional order, described by fourfold symmetry. Fig. 3b shows the corresponding gray scale plot of the scalar order parameter S; black (white) denotes nematic (isotropic-) order. Fig. 3b shows that eventually the interconnected structure transforms into an isotropic droplet domain everywhere, suggestive of spinodal decomposition (SD) due to phase ordering. Fig. 3b shows the growth of a positionally ordered array of isotropic domains in the nematic matrix. The morphology is that of a filled nematic (FN) crystal, or colloidal crystal [37]. To the best of our knowledge there are no theoretical predictions for formation of crystalline filled nematic (CFNs) materials during polymer–liquid crystal de-mixing. Similar positional ordering driven by liquid crystalline order in filled nematics has been observed [38]. Periodically aligned chains of droplets in a nematic matrix were observed experimentally [23,39] and suspensions of small oil or water droplets in anisotropic solvents revealed strong attractive forces between the droplets induced by the elastic deformation of the anisotropic suspending medium [20–23,40]. The crystal
331
25
ks
20
kφ
* max
15 10 5
10
20
30
40
50
60
Dimensionless time,
t*
70
0
80
d
0
5
10
15
20
Dimensionless time, t*
Fig. 4. (a) Compositional structure factor, (b) orientational structure factor, (c) maximum of structure factors, and (d) peak wave number of structure factors following a quench to point A for, nI = 25, nA = 2, a = 2.5, K/ = 103, N = 1.0, K/Q = 0.2, and KQ = 0.1.
332
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
a stable crystalline FN. Each s = 1/2 disclination is a triple junction. The network has defect–defect and defect– droplet edges. The droplets form a square lattice, and within each lattice there are two s = 1/2 defects. To find the topological charge corresponding to each isotropic droplet we can use the Poincare–Brouwer theorem [41], if M is a surface and n is a vector field, then the total charge ST of M is equal to the Euler characteristic v of M: S T ¼ S v V þ S E E þ S F F ¼ v;
ð39Þ
where V is the number of vertices, E is the number of edges, F is the number of faces, Sv is the charge on each vertex, SE is the charge on each edge, and SF is the charge on each face. Using Eq. (39) we find: 1 ðS drops ÞV þ 2 F ð0ÞE ¼ v ð40Þ 2 and hence the charge associated with each polymer drop is Sdrops = Sv = +1. Hence the stability of the lattice follows from the presence of positionally ordered positive and negative charges. The slight ellipticity of the drops creates
a degree of freedom in the orientation of each pair of defects; some pairs are vertical while others are on the diagonals of the square lattice. Defect lattices are common in liquid crystals and have been discussed by Bouligand [42]. Long-lived defect square lattices composed of +1 and 1/2 wedge disclination lines in uniaxial discotic nematic phases were prepared and studied by Fathollahi and White [43]. A defect square cluster in thermotropic nematic polymers under external fields known as a Lehmann cluster, involves two s = +1/2 and two s = 1/2 disclinations, and was found to be a stable quadrupolar structure [44,45]. Hence, the emergence of a defect lattice with fourfold symmetry in droplet phaseseparated morphologies is neither unexpected nor unusual. The 2N defect lattice in our case contains s = 1/2 disclinations, which it is expected to have lower Frank elastic energy than an N defect lattice of s = 1 disclinations [17]; here N denotes the number of defects. Fig. 4a–c shows the temporal evolution of the compositional structure factor W/, the orientational structure factor WS, and the maximum structure factors Wmax / , 3 Wmax S , for the temperature quench at A for K/ = 10 . The
Fig. 5. (a) Gray scale plot of the concentration profile /ð~x; ~y Þ at ~t ¼ 63:95 1005 for K/ = 104; black (white) represents the isotropic (nematic) phase. At the early stages (not shown here), the system forms a bicontinuous network that rapidly breaks up into polymer circular drops. (b) Corresponding gray scale plot of the scalar order parameter S; black (white) denotes nematic (isotropic) order. (c) Nematic texture; the arrows represent the local nematic director, and defects are marked with small solid circles.
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
time evolution of the structure factor for W/ slightly deviates with the results of the Cahn theory for isotropic SD [9,13– 15]. The orientational structure factor WS profile (Fig. 4b) has also slightly deviates the shape characteristics of that of a classical SD process. The results of time evaluation of both the structure factors (Fig. 4a and b) are consistent with the reported results [8]. The time evolution of maximum max structure factors Wmax (Fig. 4c) shows that Wmax / , WS S evolves much faster and is larger than Wmax / . It indicates that the nematic ordering fluctuations initially induce the SD and the concentration changes within the domains subsequently take place due to the coupling between the two order parameters. Fig. 4d shows the evolution of peak wave numbers k S max ðt Þ and k /max ðt Þ. Both show a weak decrease and saturation. In addition k S max ðt Þ > k /max ðt Þ, indicating that the length scale for ordering is smaller than for concentration, in agreement with the parametric value KQ = 0.1. 3.1.2. Random droplet–random texture regime: 103 < K/ < 105 The main feature of this regime is the presence of a large number of nearly spherical very small polymer droplets embedded in a random nematic matrix lattice (see Fig. 5); the droplets are randomly distributed, and hence the name random droplet.
In this small droplet regime the size of the droplets is much smaller than the system size. The size of the polymer droplets as a function of K/ is discussed in Fig. 11. The e is given by the following dimensionless droplet radius R equation (see Fig. 11 below): n ‘/ a e R¼a ¼ n; ð41Þ K/ ‘ where a 104, 6/5 < n < 17/5. Using these results in Eq. (30) we find ! e 1=n R 1 ‘/ ‘ext < 0:2 ! < 1=n 10SD/ ‘ext ‘ 10SD/ðaÞ e 1=n < 10SD/ðaÞ1=n ‘ext . !R ð42Þ ‘ Using characteristic values observed in the simulations, S = 0.3, D/ = 0.2, we find that the characteristic length scale ordering is ‘ext > R > ‘/ ‘Q
2.5
Compositional structure factor, ψφ
2.0
t* = 52 *
t = 45 15
t = 34 *
t = 24 t = 14
5
t* = 34 t* = 24 t* = 14
0.5
5
10
15
0.0 0
20
k*/2π
a
5
10
15
20
k*/2π
b 35
25
ψ
max
30
k*s max
25
k*φ max
s
Peak wave number
ψ
20 max
t* = 45
1.0
*
0 0
max φ
/ψ φ
15
max
t* = 52
1.5
*
10
t*= 64
Orientational structure factor,ψs
t* = 64 20
ψs
10
5
0 0
ð43Þ
indicating weak anchoring conditions and hence random texturing due to the lack of orienting internal surfaces; see Fig. 11 for exact calculation of length scales. We recall that phase ordering in the absence of anchoring effects leads to a
25
c
333
20 15 10 5
10
20
30
t*
40
50
60
0 0
70
d
2
4
t*
6
8
10
Fig. 6. (a) Compositional structure factor, (b) orientational structure factor, (c) maximum of structure factors, and (d) peak wave number of structure factors following a quench to point A for, nI = 25, nA = 2, a = 2.5, K/ = 104, N = 1.0, K/Q = 0.2, and KQ = 0.1.
334
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
random textures by activation of the Kibble mechanism [42,46]. We note that as K/ increases more and more smaller droplets arise, and the ratio ‘ext/R increases. Fig. 5a shows the gray scale plot of the concentration profile /ð~x; ~y Þ at ~t ¼ 63:95 1005 for K/ = 104; black (white) represents the polymer-rich (mesogen-rich) phase. The morphology corresponds to a random droplet. The polymer-rich drops are nearly circular and display positional disorder. Fig. 5b shows the corresponding gray scale plot of the scalar order parameter S; black (white) denotes nematic (isotropic) order. From these figures we can see clearly that the isotropic droplets form initially in an LCrich matrix. The random droplet has no positional order. Fig. 5c shows the corresponding nematic orientation texture. There are a large number of defects, as in common phase ordering processes. The anchoring on the polymer droplets is weak, and the random location and reduced size of the droplets leads to a random texture with many wedge disclinations of strength s = 1 and s = 1/2. The CFN observed in the previous regime is now absent and instead a randomly filled nematic (RFN) structure emerges.
Fig. 6a and b shows the temporal evolution of the compositional structure factor, W/, and the orientational structure factor, WS, for the temperature quench at A for K/ = 104. The time evolution of the structure factor for W/ (Fig. 6a) deviates from the results of the Cahn theory for isotropic SD [13–15] and consistent with the reported results [8], in particular showing a wide range of wave numbers and a significant tail. The orientational structure factor WS profile (Fig. 6b) shows a chimney shape. The funnel width of the structure factor WS becomes very small at large time step. It indicates that the phase ordering plays a significant role to form the dynamical textures in this case. The max time evolution of maximum structure factors Wmax is / ; WS max shown in Fig. 6c. Fig. 6c shows that W/ evolves much faster than Wmax S . It indicates that the concentration changes within the domains take place initially induced by the classical phase separation SD and the nematic ordering subsequently take place due to the coupling between the two order parameters. Fig. 6d shows the evolution of peak wave numbers k S max ðt Þ and k /max ðt Þ. Both show a decrease and saturation. In addition k S max ðt Þ > k /max ðt Þ, indicating that
Fig. 7. (a) Gray scale plot of the concentration profile /ð~x; ~y Þ at ~t ¼ 61:98 1005 for K/ = 106; black (white) represents the isotropic (nematic) phase. At the early stages (not shown here), the system forms a bicontinuous network that rapidly breaks up into polymer circular drops. (b) Corresponding gray scale plot of the scalar order parameter S; black (white) denotes nematic (isotropic) order. (c) Director field corresponding to figure 7a and 7b; black corresponds to isotropic polymer and white corresponds to pure liquid crystals (LCs): the arrows represent the local nematic director.
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
the length scale for ordering is smaller than for concentration, in agreement with the parametric value KQ = 0.1. 3.1.3. Small droplets regime: K/ > 105 The main feature of this regime is the absence of nematic order, and the presence of a small droplets microstructure (see Fig. 7). For sufficiently thin interfaces, the number of drops increases at a significantly higher rate. In this case the length scale ordering disappears and the smallest morphological features become of the same order as the interface thickness and nematic correlation length, R ‘/ ‘Q . ð44Þ Hence phase ordering is frustrated and no nematic ordering arises. The scalar order parameter is nearly zero. A cascade of tiny droplets microstructure arises. The characteristic thickness of the droplets is n ‘/ e R a ; ð45Þ ‘ where a = 104, n = 17/5. See Fig. 11 for exact calculation of length scales. Considering Eq. (26) we see that phase ordering frustration through small-scale morphology emerges from a balance between ! of~ h KQ/ ~ 2 ð46Þ r / 0 ! Q ’ 0; S ’ 0 oQ Te K/
335
and the small-scale concentration gradients inhibit ordering. Fig. 7a shows the gray scale plot of patterns of concentration profiles at ~t ¼ 61:98 1005 for K/ = 106. The formation of small droplets with low polymer concentration is shown. Fig. 7b shows the gray scale plot of patterns of the scalar order parameter S at the same time step. The characteristic value of S is slightly larger than zero. The compositional of the small droplet morphology is reflected in the weak scalar order parameter which fluctuates with very small amplitude around S = 0. The droplets are randomly distributed and their length fluctuates, as can be seen from Fig. 7b. Fig. 7c shows the director field nð~x; ~y Þ corresponding to Fig. 7a and b. In tandem with the scalar order parameter, which is slightly greater than zero, the director field has no significant long-range correlations, and hence phase ordering is absent. Fig. 8a and b shows the temporal evolution of the compositional structure factor, W/, and the orientational structure factor, WS, for the temperature quench at A for K/ = 106. The time evolution of the structure factor for W/ (Fig. 8a) shows secondary maximum. The double peaks indicate the presence of two dominant length scales: the droplet radius and the thickness of droplet interface. The orientational structure factor, WS, profile (Fig. 8b) is almost zero at all the time steps (here we show the final time step only). It indicates phase ordering frustration.
*
0.0004
t = 62
10
Orientational structure factor,ψ s
Compositional structure factor, ψφ
12
*
t = 54
0.0002
*
8
t = 49 *
t = 44 6
*
t = 34 *
t = 24
4
-0.0004
0
5
10
15
0
20
k*/2π
a
Peak wave number
max
15
20
25
30
35
k*/2π
30
max
ψφ
8 φ
10
35
10
6 4 2 0 0
5
b
12
ψ
0
-0.0002
2 0
c
t* = 62
25 20 15 k*s max
10 5
10
20
30
40 *
t
50
60
0 0
70
d
5
10
15
20
25
30
35
40
t*
Fig. 8. (a) Compositional structure factor, (b) orientational structure factor, (c) maximum of structure factors, and (d) peak wave number of structure factors following a quench to point A for, nI = 25, nA = 2, a = 2.5, K/ = 106, N = 1.0, K/Q = 0.2, and KQ = 0.1.
336
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
The time evolution of the maximum structure factors ;1 Wmax , W/max ;2 corresponding to the two peaks is shown / in Fig. 8c. It indicates that the small droplets structure is formed due to the classical phase separation SD and in the absence of nematic ordering (S = 0). Fig. 8d shows the evolution of peak wave numbers k /max ;1 ðt Þ, k /max ;2 ðt Þ. It shows an exponential decrease and saturation [9,47].
as a function of the x-coordinate, at a fixed y = 0.49, ~t ¼ 63:95 1005 , corresponding to the beginning of the intermediate region. The upper figure (Fig. 9a) shows that as K/ increases the wavelength and amplitude of the concentration pattern decreases. This indicates that as K/ increases the intensity of de-mixing decreases and the scale 1.0
= 10
V
= 10
φ
V
φ
φ
4
= 10
6
0.6
E*
In this section we present a comparative evaluation of the main distinguishing features of the morphologies, textures, and order in the three regimes discussed above. Fig. 9a and b shows the concentration (/), and scalar order parameter (S) corresponding to the three regimes,
0.8
3
V
3.2. Comparative regime characterization
0.4
1.0 y = 0.49 V
o
0.8
V V
= 10
3
0.2 4
φ
= 10
φ
= 10
6
0 0
φ
0.6
φ
10
20
30
40
50
60
70
80
t* 0.4
Fig. 10. Dimensionless extent of phase separation component, E*, as a function of dimensionless time, t ¼ ~t 105 .
0.2
0 0
0.2
0.4
0.6
a
0.8
1.0
1.2
x 1.0 y = 0.49
0.8
V
= 10
φ
= 10
S
V
V
0.4
3
φ
0.6
4
6
φ
= 10
0.2
0
b
0
0.2
0.4
0.6
0.8
1.0
1.2
x
Fig. 9. Evolution of (a) concentration field, /, and (b) orientational scalar order parameter, S, along x-direction at a fixed y for different values of K/ following a quench to point A for, nI = 25, nA = 2, a = 2.5, N = 1.0, K/Q = 0.2, and KQ = 0.1.
e ¼ 2R=‘ dimensionless droplet Fig. 11. Dimensionless length scales ð D diameter, ‘~ext ¼ ‘ext =‘ dimensionless extrapolation length, ‘~Q ¼ ‘Q =‘ dimensionless nematic interface thickness) as a function of dimensionless concentration interfacial thickness ‘~/ ¼ ‘/ =‘. The slopes are indicated on each curve. The graphical insets shows the computed morphology (i.e. gray-scale visualization of /ð~x; ~y Þ). The full line of slope +1 shows how the other length scales compare with the interfacial thickness. The left vertical dash-line indicate the tiny droplets cascade-random droplet transition, while the right vertical dash-line indicate the random droplet-colloidal crystal transition.
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
of de-mixing decreases from colloidal droplets to small droplets morphology. The lower figure (Fig. 9b) shows that as K/ increases the order in the nematic phase decreases and eventually phase ordering does not occur (small droplets regime). Fig. 10 represents the profiles of extent of phase separation as a function of dimensionless time for all the three regimes discussed above. The extent of phase separation is a useful measure to characterize the kinetics of phase-
337
separated structures developed by Ariyapadi and Nauman [48], and defined as E ¼
2 N X N ð/X;i;j /X;0 Þ 1 X ; N 2 j¼1 i¼1 ð/X;0 /X;l Þð/X;u /X;0 Þ
ð47Þ
where /X,i,j is the concentration of component X at the (i, j)th nodal point, /X,0 is the initial concentration, /X,l and /X,u are respectively the lower and upper binodal
Fig. 12. Fourier transform (spectrum) of the emerging textures for different values of K/, at a late time step following a quench to point A for, nI = 25, nA = 2, a = 2.5, N = 1.0, K/Q = 0.2, and KQ = 0.1.
338
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
concentration of component X, and N is the total number of nodal points in the computational domain considered in the computation. If there is no phase separation at all, e ¼ 0. Other/X,i,j = /X,0 for all nodal points, and thus E wise, for a perfectly phase-separated structure, /X,i,j approaches either one of the two binodal points, /X,l or e ! 1. Fig. 10 shows that no phase /X,u, and therefore E separation occurs up to t* 20 for all of the regimes. Significant phase separation only occurs at 20 < t* < 60 when E* increases exponentially. At t* > 60 phase separation saturates. The extent of phase separation decreases with decreasing droplet size. The present model predicts that as the characteristic length scale of the droplets decreases, de-mixing strength also decreases. e ¼ 2R=‘ Fig. 11 shows the dimensionless length scales ( D ~ dimensionless droplet diameter, ‘ext ¼ ‘ext =‘ dimensionless extrapolation length, ‘~Q ¼ ‘Q =‘ dimensionless nematic interface thickness) as a function of dimensionless concentration interfacial thickness ‘~/ ¼ ‘/ =‘. The quantity e ¼ 2R=‘ is extracted from the computed solution vector D (Q and /ð~x; ~y Þ). The quantity ‘~ext ¼ ‘ext =‘ is calculated using the computed solution vector (Q and /ð~x; ~y Þ) and Eq. (30). The quantity ‘~Q ¼ ‘Q =‘ is calculated using Eq. (31): ‘~Q ¼ KQ ‘~/ ¼ 0:1‘~/ . The slopes are indicated on each curve. The graphical insets shows the computed morphology (i.e. gray-scale visualization of /ð~x; ~y Þ). The full line of slope +1 shows how the other length scales compare with the interfacial thickness. The left vertical dash-line indicate the small random droplet transition, while the right vertical dash-line indicate the random droplet-colloidal crystal transition. The average equivalent dimensionless droplet diameter is a function of the dimensionless interface thickness ‘/; for small droplet morphologies we represent the equivalent diameter by the smallest length scale. The figure clearly shows the presence of two power lawscaling regimes:
Fig. 12 shows the Fourier transform of the emerging concentration morphologies for increasing K/. In the colloidal crystal regime (upper left) the spectrum has the fourfold symmetry corresponding to the square lattice arrangement. Increasing K/ first destroys the positional order and leads to the circular ring, typical of SD [7,10] and reflecting a random droplet phase. Further increase in K/ leads to the colloidal droplets-to-small droplets transformation and the emergence of the small droplets morphology. The model predicts that at the intermediate stage, the ring radius and thickness increase with K/.
e / ð‘/ Þn ; R
In the colloidal crystal regime the polymer drops organize in a square lattice, and due to strong anchoring, the director orients such that wedge disclination defect lines of strength s = 1/2 nucleate and form a defect lattice. As droplet become smaller and smaller, the random droplet regime takes over. The morphology is that of random droplets, and the nematic matrix displays many randomly positioned defects. The size and lack of positional order produce little effect on the matrix, which in the absence of orienting surfaces leads to many topological defects. Finally when the concentration length scale is of the order of the phase ordering scale, phase ordering is frustrated and a small droplet morphology arises. The small droplets are randomly oriented but the nematic order parameter is close to zero. This work identifies the conditions and mechanisms that lead to morphological transformations in mixtures of polymers and liquid crystals. It provides new mechanisms in how to create patterns vastly different among the droplet morphologies characteristic from
ð48Þ
where n = 17/5 corresponds to the tiny droplets cascade regime, and n = 6/5 to the colloidal crystal regime. The figure shows that the computed extrapolation length ‘ext is a decreasing function of the interface thickness ‘/ (see Eq. (30) for ‘ext); the reason being that as ‘/ increases the denominator S(D/) increases. The figure shows that as the random droplet-colloidal crystal transition occurs e dewhen R ‘exp. As ‘/ decreases the droplet size D creases to values comparable to the interface thicknesses (‘/, ‘Q), in which case ordering is frustrated and the small droplets morphology arises. The physical interpretation that can be extracted from Fig. 11 is that by adding a surfactant, phase separation produces smaller droplets and when the droplets are sufficiently small, weak anchoring sets in and hence the presence of drops does not produce a texture since the director does not sense the presence of curvature.
4. Conclusions Scaling analysis of coupled phase separation–phase ordering–texturing processes indicates the possibility of colloidal crystal ordering as well as phase ordering frustration. To test the nature of these regimes a previously formulated model obtained by merging the Flory– Huggins–Chan–Hilliard model of phase separation, with the Onsager–Maier–Saupe phase ordering model, and with the Landau–Gennes tensor order parameter model, was analyzed numerically. The three coupled phase separation–phase ordering–texturing processes are associated with four length scales: the droplet size R (texturing), extrapolation length ‘ext (texturing), concentration interface thickness ‘/ (phase separation), and nematic order interface thickness ‘Q. Length scale ordering gives rise to the following regimes: (i) R > ‘ext > ‘/ ‘Q: colloidal crystal–defect lattice texture regime, (ii) ‘ext > R > ‘/ ‘Q: random droplet–random texture regime, (iii) R ‘/ ‘Q: small droplet regime–frustrated phase ordering regime.
S.K. Das, A.D. Rey / Computational Materials Science 38 (2006) 325–339
off-critical quenches in the classical spinodal decomposition process. It was shown that including phase ordering and texturing processes provides new opportunities to create complex multi scale architectures. Acknowledgement This work was supported in part by the ERC Program of the National Science Foundation under Award Number EEC-9731680, and by the Natural Science and Engineering Research Council of Canada. References [1] P.S. Drzaic, Liquid Crystal Dispersions, World Scientific, Singapore, 1995. [2] National Research Council Report, Liquid Crystalline Polymers, National Academy Press, Washington, DC, 1990. [3] J.W. Doane, in: B. Bahadur (Ed.), Liquid Crystals: Applications and Uses, World Scientific, New Jersy, 1990, p. 362 (Chapter 14). [4] G.P. Crawford, S. Zumer, Liquid Crystals in Complex Geometries Formed by Polymer and Porous Networks, Taylor and Francis, London, 1996. [5] M. Graca, S.A. Wieczorek, R. Holyst, Macromolecules 36 (2003) 6903. [6] A.M. Lapena, S.C. Glotzer, S.A. Langer, A.J. Liu, Phys. Rev. E 60 (1) (1999) R29. [7] T. Kyu, H.W. Chiu, Polymer 42 (2001) 9173. [8] A. Matsuyama, R.M.L. Evans, M.E. Cates, Phys. Rev. E 61 (3) (2000) 2977. [9] J.R. Dorgan, J. Chem. Phys. 98 (1993) 9094. [10] Z. Lin, H. Zhang, Y. Yang, Macromol. Theory Simulat. 6 (1997) 1153. [11] H.W. Chiu, T. Kyu, J. Chem. Phys. 110 (12) (1999) 5998. [12] C.B. Muratov, E. Weinan, J. Chem. Phys. 116 (11) (2002) 4723. [13] P.K. Chan, A.D. Rey, Macromol. Theor. Simulat. 4 (1995) 873. [14] P.K. Chan, A.D. Rey, Macromolecules 29 (1996) 8934; P.K. Chan, A.D. Rey, Macromolecules 30 (1997) 2135. [15] J.W. Cahn, Trans. Metall. Soc. AIME 242 (1968) 166. [16] S.K. Das, A.D. Rey, J. Comput. Mater. Sci. 29 (2004) 152. [17] P.G. de Gennes, J. Prost, The Physics of Liquid Crystals, Oxford University Press, New York, 1993. [18] S.K. Das, A.D. Rey, J. Chem. Phys. 121 (9) (2004) 9733. [19] A.J. Liu, G.H. Fredrickson, Macromolecules 26 (1983) 2817; A.J. Liu, G.H. Fredrickson, Macromolecules 29 (1996) 800. [20] Y.Ch. Loudet, P. Barois, P. Poulin, Nature (London) 407 (2000) 611.
339
[21] J.C. Loudet, P. Poulin, Phys. Rev. Lett. 87 (2001) 165503. [22] P. Poulin, H. Stark, T.C. Lubensky, D.A. Weitz, Science 275 (1997) 1770. [23] I.I. Smalyukh, S. Chernyshuk, B.I. Lev, A.B. Nych, U. Ognysta, V.G. Nazarenko, O.D. Lavrentovich, Phys. Rev. Lett. 93 (11) (2004) 117801. [24] S.K. Das, A.D. Rey, Mol. Simulat. 31 (2–3) (2005) 201. [25] A.D. Rey, D. Grecov, S.K. Das, J. Ind. Eng. Chem. Res. 43 (23) (2004) 7343. [26] S.K. Das, A.D. Rey, in: Proceedings of NSTI-Nanotechnology Conference, Boston, Massachusetts, USA, March 7–11, vol. 3, 2004, p. 91. [27] A. Matsuyama, T. Kato, J. Chem. Phys. 105 (1996) 1654. [28] A. Matsuyama, T. Kato, J. Chem. Phys. 108 (1998) 2067. [29] A.D. Rey, T. Tsuji, Macromol. Theory Simul. 7 (1998) 623. [30] P.J. Flory, Principles of Polymer Chemistry, Cornell University, Ithaca, 1953. [31] W. Maier, A. Saupe, Z. Naturforsch. A 14a (1959) 882. [32] L. Onsager, Ann. N.Y. Acad. Sci. 51 (1949) 627. [33] F.C. Frank, Discuss. Faraday Soc. 25 (1958) 1. [34] P.J. Van der Houwen, Appl. Numer. Math. 20 (1996) 261. [35] A.W. Al-Khafaji, J.R. Tooley, Numerical Methods in Engineering Practice, The Dryden Press, Saunders College Publishing, London, 1986. [36] W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flannery, Numerical recipes in FORTRAN: the art of scientific computing, second ed., Cambridge University press, New York, 1992 (Chapters 12–13). [37] B.I. Lev, K.M. Aoki, P.M. Tomchuk, H. Yokoyama, Condens. Matter Phys. 6 (2003) 169. [38] V.G. Nazarenko, A.B. Nych, B.I. Lev, Phys. Rev. Lett. 87 (7) (2001) 075504-1. [39] P.R. Ten Wolde, D. Frenkel, Science 277 (1997) 1975. [40] P.D. Olmsted, W.C.K. Poon, T.C.B. Mcleish, N.J. Terrill, A.J. Ryan, Phys. Rev. Lett. 81 (1998) 373. [41] R.S. Millman, G.D. Parker, Elements of Differential Geometry, Prentice Hall, New Jersey, 1977. [42] Y. Bouligand, Defect and Textures, in: D. Demus, J. Goodby, G.W. Gray, H.-W. Spiess, V. Vill (Eds.), Physical Properties of Liquid Crystals, Wiley-VCH, Weinham, Germany, 1999. [43] B. Fathollahi, J.L. White, J. Rheol. 38 (1994) 1591. [44] S. Hudson, E.L. Thomas, Phys. Rev. E 44 (1991) 8128. [45] J.E. Zimmer, J.L. White, Carbon 21 (3) (1983) 323. [46] A.D. Rey, J. Chem. Phys. 120 (5) (2004) 2010. [47] C. Sagui, A.M. Somoza, R.C. Desai, Phys. Rev. E 50 (1994) 4865. [48] M.V. Ariyapadi, E.B. Nauman, in: R.J. Roe (Ed.), Computer Simulation of Polymers, Prentice Hall, 1992.