Simulation of spherulite growth using a comprehensive approach to modeling the first-order isotropic/smectic-A mesophase transition
A comprehensive modeling and simulation study of the first-order isotropic/smectic-A transition is presented and applied to phase diagram computation and two-dimensional spherulite growth. An approach based on nonlinear optimization, that incorporates experimental data (from 12CB, dodecyl-cyanobiphenyl), is used to determine physically realistic model parameters. These parameters are then used in conjunction with an optimized phase diagram computation method. Additionally, a time-dependent formulation is presented and applied to the study of two-dimensional smectic-A spherulite growth. These results show the growth kinematics and defect dynamics of nanoscale smectic-A spherulite growth in an isotropic phase with an initially radial layer configuration.
Liquid crystallinity and other forms of self-organization are key phenomena both technologically and in Nature. Liquid crystalline order ranges from liquids that show some degree of orientational order to those that, in addition, show various types of translational order. The myriad of types of material that exhibit this behavior range from traditional low molecular mass molecules, currently used in display technology, to biological membranes composed of phospholipids . To date, the main focus of liquid crystal research has been on the simplest of the class of mesophases where some degree of orientational order is present: the nematics. Research on mesophases that also show a degree of translational order, including smectics and columnar liquid crystals, has been less abundant. Recognizing the increasing importance of these mesophases, particularly in biological systems, the need for practical methods to access the time and length scales at which these phenomena occur becomes more important.
Experimental work in this general field has made great progress in the basic understanding of translationally ordered liquid crystals [2, 3, 4, 5, 6]. Nonetheless, it is currently infeasible to access much of the dynamic phenomena of translational phase-ordering processes. Recent experimental work has begun to address these issues [4, 5], but theoretical approaches are currently the only way to access the length scales (nanometers) and time scales (nanoseconds) at which liquid crystal dynamics occur. The use of high-order models in conjunction with advanced numerical simulation techniques has shown a great deal of promise for theoretical study [7, 8, 9, 10]. Recent computational advances have allowed for the possibility of simulation in greater detail than ever before.
to present a comprehensive approach to modeling and simulation of the first-order isotropic/smectic-A transition.
to determine phenomenological model parameters through incorporation of experimental data (from 12CB, dodecyl-cyanobiphenyl).
to efficiently compute a phase diagram predicted by the model and parameter set.
to study the two-dimensional growth kinematics and defect dynamics of an initially radial smectic-A spherulite in an isotropic matrix.
This approach builds upon previous work , which incorporates experimental data into the phenomenological model and derives equations for phase diagram computation. A time-dependent formulation [9, 10] and the nano-scale growth of an initially radial spherulite are presented. This work is organized as follows: a brief background on relevant types of liquid crystalline order is given (Section 2.1), the model and simulation approach are explained (Sections 2.3-2.6), and simulation results are presented and discussed (Section 3).
2 Background and theory
2.1 Liquid crystalline order
Liquid crystalline phases or mesophases are materials which exhibit partial orientational and/or translational order. They are composed of anisotropic molecules which can be disc-like (discotic) or rod-like (calamitic) in shape. Thermotropic liquid crystals are typically pure-component compounds that exhibit mesophase ordering most greatly in response to temperature changes. Lyotropic liquid crystals are mixtures of mesogens (molecules which exhibit some form of liquid crystallinity), possibly with a solvent, that most greatly exhibit mesophase behavior in response to concentration changes. Effects of pressure and external fields also influence mesophase behavior. This work focuses the study of calamitic thermotropic liquid crystals which exhibit a first-order mesophase transition.
An unordered liquid, where there is neither orientational nor translational order (apart from an average intermolecular separation distance) of the molecules, is referred to as isotropic. Liquid crystalline order involves partial orientational order (nematics) and, additionally, partial translational order (smectics and columnar mesophases). The simplest of the smectics is the smectic-A mesophase, which exhibits one-dimensional translational order in the direction of the preferred molecular orientational axis. It can be thought of as layers of two-dimensional fluids stacked upon each other. Schematic representation of these different types of ordering are shown in Figure 1.
2.2 Order parameters and the phenomenological model
Theoretical characterization of mesophase order is accomplished using order parameters that adequately capture the physics involved. These order parameters typically have an amplitude and phase associated with them. In order to characterize the partial orientational order of the nematic phase, a second order symmetric traceless tensor can be used :
where are the eigenvectors of Q-tensor, which characterize the average molecular orientational axes, and are scalars which represent the extent to which the molecules conform to the average orientational axes [13, 14, 15]. Uniaxial order is characterized by and , which correspond to the maximum eigenvalue(and its corresponding eigenvector) of the Q-tensor, . Biaxial order is characterized by and , which correspond to the lesser eigenvalues and eigenvectors, .
The smectic-A mesophase has one-dimensional translational order in addition to the orientational order found in nematics. Characterizing this mesophase can be accomplished through the use of primary (orientational) and secondary (translational) order parameters together . This is accomplished using the tensor order parameter (1) and the complex order parameter :
where is the phase, is the scalar amplitude of the density modulation. The density wave vector, which describes the average orientation of the smectic-A density modulation, is defined as . The smectic scalar order parameter characterizes the magnitude of the density modulation, and is used in a dimensionless form in this work. In the smectic-A mesophase the preferred orientation of the wave vector is parallel to the average molecular orientational axis, .
A Landau-de Gennes type model for the first order isotropic/smectic-A phase transition is used that was initially presented by Mukherjee, Pleiner, and Brand [11, 7] and later extended by adding nematic elastic terms [17, 18, 8]:
where is the free energy density, is the free energy density of the isotropic phase, terms 1-5 are the bulk contributions to the free energy, terms 6-7 are couplings of nematic and smectic order; both the bulk order and coupling of the nematic director and smectic density-wave vector, respectively. Terms 8-10 are the nematic and smectic elastic contributions to the free energy, respectively. The order parameters are defined in (1-2), is temperature, / are the hypothetical second order transition temperatures for isotropic/nematic and isotropic/smectic-A mesophase transitions (refer to  for more detail), and the remaining constants are phenomenological parameters. Further explanation and justification for the use of this high-order model can be found in .
2.3 Homogeneous free energy and parameter determination
Following past work , in order to compute the phase diagram from the free energy equation (3) a homogeneous uniaxially-ordered volume assumption can be used, resulting in a simplified point-volume free energy density:
A major complication of applying this model is the determination of a suitable set of model parameters. In order to overcome this challenge, a nonlinear programming formulation is derived which then allows to the application of nonlinear solution methods. Utilizing experimental data [19, 20, 21] and minimization criteria of the homogeneous free energy (4) [7, 22], a nonlinear optimization problem is formulated:
where the objective function minimizes the change in the layer spacing between the unknown bulk transition value, , and the known value at some minimum valid temperature for the model, . Each constraint is evaluated at the corresponding temperatures , the bulk transition temperature, and , the minimum valid temperature for the model parameter set.
2.4 Phase diagram determination
Utilizing the homogeneous free energy (4) and model parameters determined from the nonlinear optimization solution in Section 2.3, a phase diagram for the system can be computed. Figure 2 shows a schematic of the phase diagram of a first-order isotropic/smectic-A phase transition. Coexistence regions are present, due to the first-order nature of the transition, where both the isotropic and smectic-A phase are either stable or metastable.
Minimization of the homogeneous free energy (4) is computationally challenging due to the presence of three degrees of freedom. This can be simplified by parameterization of the free energy equation  using minimization invariants and the assumption that the smectic-A phase exists. The resulting free energy equation, parameterized as a function of the nematic scalar order parameter (see (1)), is:
where is the nematic scalar order parameter with the assumption that the smectic-A phase is stable/meta-stable. The minima of (2.4) are easily found in that they are the roots of a polynomial, . Once the minima are determined, the validity of the assumption that the smectic-A phase exists (at the specified temperature) can be tested by i) verifying that the computed and the corresponding values of / are real, and ii) verification of the minimization criteria of the full homogeneous free energy (4). The equations to determine and from are derived from the application of the minimization criteria to the homogeneous free energy equation (4):
2.5 Time-dependent formulation
The Landau-Ginzburg time-dependent formulation is used to capture the kinetics of the phase transition. Due to the higher order derivative term in the free energy functional, a higher order functional derivative must be used. Additionally, in order to utilize standard numerical solution techniques, the complex order parameter (2) is separated into its real and imaginary contributions:
The general form of the time-dependent formulation is as follows :
where / is the rotational/smectic viscosity, is the heterogeneous free energy density (3), and the volume. As previously mentioned, a higher order functional derivative must be used due to the second-derivative term in the free energy (3):
where corresponds to the order parameter.
The use of the full Q-tensor in this time-dependent model does not neglect biaxiality as was assumed in Sections 2.3-2.4. Substituting (8), the free energy (3), and high order functional derivative (11) into the time-dependent formulation (9) yields the closed set of simulated equations:
where the asterisk denotes an nondimensionalized value, the superscript denotes the symmetric/traceless portion of a tensor, and is the ratio of the smectic and rotational viscosities. The nondimensionalized model parameters are as follows:
where is the simulation-specific imposed length scale.
2.6 Simulation conditions and computational approach
A square computational domain with imposed length scales of (approximately 25 layers) and (approximately 50 layers) were used in two separate simulations. Referring to Figure 3a, both symmetry and Neumann boundary conditions were used to simulate bulk conditions. Symmetry conditions for the Q-tensor (1) must take into account vector symmetry, which results in the following boundary conditions :
The initial condition for both simulations was a smectic-A spherulite in an initially radial layer configuration (see Figure 3b). The radius of the spherulite was initially set to , a value on the order of the layer spacing at (Section 2.4), . The initial value used for , , and the layer spacing correspond to the homogeneous values at , determined from the computed phase diagram. The Heaviside step function was used to generate the initial spherulite. The constraint that the spherulite does not impinge on the domain boundaries was verified post-simulation.
A commercial package, Comsol Multiphysics, was used to solve the time-dependent model (2.5). Quadratic Lagrange basis functions were used for the Q-tensor variables and quartic Hermite basis functions used for the complex order parameter components. Standard numerical techniques were utilized to ensure convergence and stability of the solution. This platform does not support adaptive mesh refinement, thus a uniform mesh was used with a density of approximately nodes/ for the 25 layer simulation and nodes/ the 50 layer simulation. Previous simulations using this model and numerical method have shown good agreement with both past experimental and theoretical findings [9, 10]. Additionally, exhaustive past work using this numerical method and the Landau-de Gennes model for the first-order isotropic/nematic phase transition [25, 26, 27, 28] has served to further validate this simulation approach.
3 Results and discussion
3.1 Phase diagram computation
The algorithm for phase diagram computation presented in Section 2.4 was implemented using a high-level numerical programming language . The resulting phase diagram is shown in Figure 4; refer to the phase diagram schematic, Figure 2, for a detailed explanation of the features.
The experimentally observed bulk transition temperature and nematic scalar order parameter  are exactly matched using the parameter determination method presented in Section 2.3. Additionally, the smectic layer spacing on the order of  is also well reproduced below the lower stability limit of the isotropic phase, . Above this temperature the model predicts a layer spacing trend that increases approaching the super-heated stability limit of the smectic-A phase.
3.2 Two-dimensional spherulite growth
Simulation results from the growth of an initially radial spherulite are shown in Figure 5 from the 50 layer simulation (see Section 2.6). Past work on the isotropic/nematic mesophase transition [14, 30, 25, 26, 27, 28, 24] and smectic-A filament growth  provide a great deal of insight into the two general growth processes observed: shape dynamics  and self-similar growth . As will be shown, these simulation results definitely show that the growth of an initially radial textured smectic-A spherulite in an isotropic phase follows a similar topological growth process as observed in its nematic counterpart [25, 26, 27, 28].
Shape dynamic growth
The shape dynamic growth regime is dominated by transient texturing and interfacial forces. When the spherulite radius is on the order of the smectic coherence length:
long-range energy effects (gradients in molecular and layer orientation) are dominant over short-range energy effects (gradients in the bulk order parameters and ). As the spherulite radius increases, interfacial anchoring affects the core of the spherulite less and short-range energy gradients become dominant in this region. This results in increased texturing and defect dynamics [13, 14]. Referring to Figure 6, this process is observed for an initially radial smectic-A spherulite.
The shape dynamic growth regime of this spherulite morphology has been found to have topologically similar dynamics compared to the growth of initially radial nematic spherulites [25, 26, 27, 28]. In agreement similar nematic growth simulations [25, 26, 27, 28], the initially imposed isotropic core is observed to undergo a self-selected transition to a disclination, shown in Figure 7a. This initial morphology minimizes elastic long-range energy through the formation of a high-energy defect. As the spherulite radius exceeds the smectic coherence length (15), long-range elastic gradients become energetically favorable compared to gradients in short-range order. Subsequently, the disclination splits into a pair of disclinations [33, 14]. This is a topologically equivalent texture that minimizes short-range energy, where the energy of a defect is proportional to its the square of its strength . The smectic-A layer configurations corresponding with this splitting process is shown in Figure 6.
The 25 layer simulation, with a more refined mesh, was used to both verify mesh independence of the 50 layer simulation (see Section 3.2.2) and obtain a more refined view of the smectic disclinations observed in this work. In order to identify disclination defects, the degree of biaxiality (see Section 2.2) can be computed as follows [34, 35]:
The core of the disclination is shown in Figure 7a as was determined from the 25 layer simulation. The simulated structure of this disclination is both topologically correct and its biaxial halo/uniaxial center structure well-reproduced, compared to its nematic counterpart [36, 37, 38, 35, 39]. The spherulite length scale of the 25 layer simulation was not large enough to observe the splitting of the disclination into a pair of , thus Figure 7b shows the partially biaxial core structure of the computed in the 50 layer simulation. The core of this defect is also in agreement with its nematic counterpart, with a solely biaxial core .
The radius of the spherulite versus time was determined for both simulations. Figure 8a shows this data and a comparison of the results from the 50 layer and more refined 25 layer simulations show good agreement. Deviation in the early stages of growth result from the more refined 25 layer simulation capturing the formation and core structure of the disclination more accurately. The convergence of the evolution of the two spherulite radii establishes the accuracy of the 50 layer simulation results beyond the radius of convergence.
Following the shape dynamic growth regime, the spherulite forms a scale-independent, or self-similar, shape . Spherulite growth obeys a power law in this growth regime, where is both experimentally observed [40, 41] and theoretically predicted [26, 42] to approach for phase-ordering processes quenched below the lower stability limit of the isotropic phase. A power law fit of the final two points of the spherulite radius evolution (Figure 8, 50 layer simulation) shows that the simulation results in this work correctly predict a self-similar growth regime, where is found. Additionally, the results from Figure 8 show that the transition from the initial shape dynamics regime to the self-similar growth regime occurs at a relatively small length scale compared to the growth of initially radial textured nematic spherulites .
Figure 9 shows the order parameter profiles during the beginning of the self-similar growth regime (Figure 5iii). Good agreement is observed between the bulk values of the order parameters (determined in Section 3.1) and those found in the dynamic simulation. The imposition of layer compression/expansion and curvature from the spherulite morphology and interfacial anchoring result in a decrease of the smectic-A order from that of the ideal homogeneous layer configuration. Additionally, a local increase in smectic-A order is observed in the region immediately outside of the spherulite core where layer compression/expansion is minimized by the presence of the defect structure within the core (See Figure 6iii). This phenomenon is an example of the competition between gradients in long- and short-range order. The existence of such high smectic-A ordering in the presence of substantial layer curvature is in agreement with both theoretical and experimental observations that layer bending is a low energy distortion in the smectic-A mesophase [11, 43].
A comprehensive approach to modeling the first-order isotropic/smectic-A phase transition was presented and applied to phase diagram computation and growth of an initially radial smectic-A spherulite in an isotropic phase. Summarizing the results determined from this work are as follows:
Shape dynamics in the early spherulite growth period were found in agreement with the past theoretical work on a similar system, the first-order isotropic/nematic transition [14, 30, 25, 26, 27, 28, 24], where the initial imposed isotropic core forms a self-selected disclination. As the spherulite radius increases, this disclination splits into two disclinations, which is governed by a competition of long- and short-range ordering in the presence of interfacial anchoring [13, 14] (Figures 5-7).
A self-similar spherulite growth regime was observed following the initial shape dynamic growth regime where a power law fit of the spherulite radius was shown to approach (Figure 8), in agreement with past experimental [40, 41] and theoretical [26, 42] studies of mesophase spherulite growth under a “deep” quench.
The general structure of both and smectic disclinations were found in agreement with studies of their nematic counterparts [33, 36, 37, 38, 35, 39]. Note that due to computational constraints, the mesh density used in this simulation was adequate only to the extent of resolving the general structure of the defect core and not resolution with full detail.
The simulation results presented show that the use of a high-order phenomenological model and experimentally based model parameters results in a substantially more complete reproduction of the physical smectic-A system. Though current computational resources restrict this initial work to two-dimensions, these results show a strong correlation to past experimental [33, 36, 37, 38, 35, 39, 40, 41] and theoretical observations [14, 30, 25, 26, 27, 28, 42, 24]. These promising results show that three-dimensional simulation of this model could be used to study the formation and dynamics of smectic-A spherulites at length and time scales currently inaccessible via experimental study.
This work was supported, in part, by a grant from the Natural Science and Engineering Council of Canada.
- Fisch, M. R. Liquid Crystals, Laptops and Life, volume 23 of World Scientific Series in Contemporary Chemical Physics. World Scientific, (2004).
- Harrison, C., Adamson, D. H., Cheng, Z., Sebastian, J. M., Sethuraman, S., Huse, D. A., Register, R. A., and Chaikin, P. M. Science 290(5496), 1558–1560 (2000).
- Harrison, C., Cheng, Z., Sethuraman, S., Huse, D. A., Chaikin, P. M., Vega, D. A., Sebastian, J. M., Register, R. A., and Adamson, D. H. Phys. Rev. E 66(1), 011706 Jul (2002).
- Ruiz, R., Sandstrom, R., and Black, C. Advanced Materials 19(4), 587–591 (2007).
- Michel, J.-P., Lacaze, E., Goldmann, M., Gailhanou, M., de Boissieu, M., and Alba, M. Physical Review Letters 96(2), 027803 (2006).
- Kleman, M. and Friedel, J. Reviews of Modern Physics 80(1), 61 (2008).
- Mukherjee, P. K., Pleiner, H., and Brand, H. R. European Physical Journal E: Soft Matter 4, 293–297 (2001).
- Biscari, P., Calderer, M., and Terentjev, E. Phys Rev E Stat Nonlin Soft Matter Phys 75(5), 051707 May (2007).
- Abukhdeir, N. and Rey, A. Solid State Phenomena 139, 135–140 (2008).
- Abukhdeir, N. and Rey, A. New Journal of Physics 10(6), 063025 (17pp) (2008).
- de Gennes, P. and Prost, J. The Physics of Liquid Crystals. Oxford University Press, New York, second edition, (1995).
- Abukhdeir, N. and Rey, A. In Modelling and Simulation. International Association of Science and Technology for Development (IASTED), ACTA Press, (2007).
- Rey, A. and Denn, M. Annual Review of Fluid Mechanics 34(1), p233 – (2002).
- Yan, J. and Rey, A. D. Phys. Rev. E 65(3), 031713 Feb (2002).
- Rey, A. D. Soft Matter 3, 1349 – 1368 (2007).
- Toledano, J.-C. and Toledano, P. The Landau Theory of Phase Transitions: Application to Structural, Incommensurate, Magnetic, and Liquid Crystal Systems (World Scientific Lecture Notes in Physics). World Scientific Pub Co Inc, (1987).
- Brand, H. R., Mukherjee, P. K., and Pleiner, H. Physical Review E: Statistical, Nonlinear, and Soft Matter Physics 63, 061708/1–061708/6 (2001).
- Mukherjee, P. K., Pleiner, H., and Brand, H. R. The Journal of Chemical Physics 117(16), 7788–7792 (2002).
- Coles, H. J. and Strazielle, C. Molecular Crystals and Liquid Crystals 55, 237–50 (1979).
- Urban, S., Przedmojski, J., and Czub, J. Liquid Crystals 32, 619–624 (2005).
- DasGupta, S. and Roy, S. Physics Letters A 306, 235–242(8) (2003).
- Prostakov, A. E. Crystallography Reports (Translation of Kristallografiya) 47, 856–858 (2002).
- Barbero, G. and Evangelista, L. R. An Elementary Course on the Continuum Theory for Nematic Liquid Crystals (Series on Liquid Crystals , Vol 3). World Scientific Publishing Company, (2000).
- Abukhdeir, N., Soule, E., and Rey, A. Arxiv preprint arXiv:0805.3539 (2008).
- Wincure, B. and Rey, A. D. The Journal of Chemical Physics 124(24), 244902 (2006).
- Wincure, B. and Rey, A. Continuum Mechanics and Thermodynamics 19(1), 37–58 June (2007).
- Wincure, B. and Rey, A. Nano Letters 7(6), 1474–1479 (2007).
- Wincure, B. and Rey, A. D. Liquid Crystals 34(12), 1397–1413 (2007).
- Eaton, J. Gnu Octave Manual. Network Theory Ltd., (2002).
- Sharma, D. and Rey, A. D. Liquid Crystals 30(3), p377 – (2003).
- Rey, A. D. and Abukhdeir, N. M. Journal of Non-Newtonian Fluid Mechanics 153(2-3), 177–182 August (2008).
- Sethna, J. P. Statistical Mechanics: Entropy, Order Parameters and Complexity. Clarendon Press, (2007).
- Schopohl, N. and Sluckin, T. J. Phys. Rev. Lett. 59(22), 2582–2584 Nov (1987).
- Kaiser, P., Wiese, W., and Hess, S. J. Non-Equilib. Thermodyn 17, 153–169 (1992).
- Kralj, S., Virga, E. G., and Zumer, S. Phys. Rev. E 60(2), 1858–1866 Aug (1999).
- Hudson, S. D. and Larson, R. G. Phys. Rev. Lett. 70(19), 2916–2919 May (1993).
- Sonnet, A., Kilian, A., and Hess, S. Phys. Rev. E 52(1), 718–722 Jul (1995).
- Sigillo, I., Greco, F., and Marrucci, G. Liquid Crystals 24(3), p419 – 425 (1998).
- Andrienko, D. and Allen, M. P. Phys. Rev. E 61(1), 504–510 Jan (2000).
- Dierking, I. Appl. Phys. A: Mater. Sci. Process. 72, 307–310 (2001).
- Dierking, I. and Russell, C. Physica B (Amsterdam, Neth.) 325, 281–286 (2003).
- Huisman, B. A. H. and Fasolino, A. Physical Review E (Statistical, Nonlinear, and Soft Matter Physics) 76(2), 021706 (2007).
- Kleman, M. and Lavrentovich, O. Soft Matter Physics: An Introduction. Springer, New York, first edition, (2003).