Theoretical and Numerical Studies of Entropy-Driven Self-Assembled Smectic-A Monolayer Membranes by Hao Tu B.Sc., University of Science and Technology of China; Hefei, Anhui, P. R. China, 2008 A dissertation submitted in partial fulfillment of the requirements for the degree of Doctor of Philosophy in Department of Physics at Brown University PROVIDENCE, RHODE ISLAND May 2013 c Copyright 2013 by Hao Tu This dissertation by Hao Tu is accepted in its present form by Department of Physics as satisfying the dissertation requirement for the degree of Doctor of Philosophy. Date Robert Pelcovits, Ph.D., Advisor Recommended to the Graduate Council Date Jay Tang, Ph.D., Reader Date Thomas Powers, Ph.D., Reader Approved by the Graduate Council Date Peter Weber, Dean of the Graduate School iii Curriculum Vitae Hao Tu is born on June 4th, 1986 in Hefei, Anhui, P. R. China. He attended University of Science and Technology of China in 2004 and graduated with a B. Sc. in applied physics in 2008. He started his Ph.D. program in physics at Brown University in September, 2008. He joined the research group of Prof. Pelcovits in 2009. His major research interest is the theory of free energy of Smectic-A monolayer self assembled from solution. He received the 2012-2013 Physics Merit Dissertation Fellowship of Brown University. Publications C. Nadir Kaplan, Hao Tu, Robert A. Pelcovits and Robert B. Meyer. ”Theory of depletion-induced phase transition from chiral smectic-A twisted ribbons to semi- infinite flat membranes”, Physical Review E, 2010, 021701. iv Acknowledgements During my pursuit of a Ph.D degree of physics at Brown University, I have been fortunate to have received help from a number of people. Without of their help, I would not have been able to accomplish the work presented here. First, I thank my advisor Prof. Robert Pelcovits. Under him, I have learned much about conducting research, including the basis of independent thinking, how to overcome bottlenecks, and so on. He has taught me a lot about liquid crystalline phases from the perspective of a scholar and also offered me help in my non-academic life as a friend. Second, I would like to thank my committee members, Prof. Jay Tang and Prof. Thomas Powers. They have given me valuable suggestions during my preliminary exam, thesis writing and dissertation defense. I appreciate the help given by schol- ars at Brandeis University: Prof. Zvonomir Dogic, Prof. Robert Meyer, C. Nadir Kaplan, Prerna Sharma, Mark Zakhary, Edward Barry and Thomas Gibaud. Collab- oration with them has been enjoyable and helpful and their experimental work has provided plenty of ideas for our theoretical model. The combination of theory and experimentation has made this work more meaningful. Prof. Michael Kosterlitz from Brown University also helped me a lot during a study group about quantum phase transitions. And I appreciate all the guidance given by the staff in the department. I also thank those friends who offered me help in both academic and non-academic v fields. My officemates, especially Chenjie Wang, have always been available to discuss problems raised in research and coursework with me. My roommates have helped me in many aspects in life. I shall always cherish the friendship of Mingming Jiang, Dongfang Li and Shu Wang, together with my friends since childhood, Jiwei Liu, Chuan Qin and Heming Zhen, just to name a few. I cannot list all the names of those who have supported me throughout this experience here but many thanks regardless. This work is dedicated to my parents. vi Contents Curriculum Vitae iv Acknowledgments v 1 Introduction 1 1.1 Nematic and Smectic Phases of Liquid Crystal . . . . . . . . . . . . . 2 1.2 Director Field, Order Parameter . . . . . . . . . . . . . . . . . . . . . 3 1.3 Landau Theory for the Nematic-Isotropic Phase Transition . . . . . . 4 1.4 Elasticity Theory of Liquid Crystals . . . . . . . . . . . . . . . . . . . 6 1.5 Topological Defects in Liquid Crystals . . . . . . . . . . . . . . . . . 8 1.6 Order Parameter and Elastic Theory of Smectic Phase . . . . . . . . 10 1.7 Landau Theory for Nematic-Smectic-A Phase Transition and Its Sim- ilarity with Superconductors . . . . . . . . . . . . . . . . . . . . . . . 11 1.8 Coherence Length, Penetration Length and TGB Phase . . . . . . . . 13 1.9 Birefringence of Liquid Crystals . . . . . . . . . . . . . . . . . . . . . 14 1.10 Excluded Volume Effect and Depletion Force . . . . . . . . . . . . . . 16 1.11 An Introduction to Membrane Systems . . . . . . . . . . . . . . . . . 19 1.12 Monte Carlo Simulation . . . . . . . . . . . . . . . . . . . . . . . . . 20 2 Basic Experimental Methods and Our Core Theoretical Model 22 2.1 Experimental Methods and Facts . . . . . . . . . . . . . . . . . . . . 23 2.2 Core Theoretical Model . . . . . . . . . . . . . . . . . . . . . . . . . . 24 2.3 Application of the Model to Flat Disks . . . . . . . . . . . . . . . . . 28 3 Application of the Model to Twisted Ribbons 33 3.1 Derivation of the Free Energy for Membranes with Arbitrary Geometry 34 3.2 Derivation of Euler-Lagrange Equations for Twisted Ribbons . . . . . 37 3.3 Results on Twisted Ribbons . . . . . . . . . . . . . . . . . . . . . . . 40 3.4 Conclusions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 4 Instability from Disks to Twisted Ribbons 53 vii 4.1 Experimental Facts . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54 4.2 Model for Rippled Disks and the MC Method . . . . . . . . . . . . . 56 4.3 Results on Rippled Disks . . . . . . . . . . . . . . . . . . . . . . . . . 62 4.4 Instability from a Disk to a TR upon Stretching . . . . . . . . . . . . 64 4.5 Summary . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 73 5 Smectic-A Disks Formed by Mixture of Molecules with Opposite Handedness 75 5.1 Experimental Facts . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76 5.2 Theoretical Model and Numerical Method . . . . . . . . . . . . . . . 78 5.3 Results on Mixture Membranes . . . . . . . . . . . . . . . . . . . . . 81 5.4 A Future MC Model for the Membranes . . . . . . . . . . . . . . . . 91 6 Conclusion 94 viii List of Tables 4.1 Dependence of geometrical parameters of rippled disks on the line tension. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 ix List of Figures 1.1 Landau free energy density for a first order phase transition. . . . . . 5 1.2 A demonstration of splay twist and bend deformations in a nematic liquid crystal. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 1.3 A sketch of the directors near 2d disclinations in nematics. . . . . . . 9 1.4 An illustration of excluded volume effects. . . . . . . . . . . . . . . . 17 1.5 A demonstration of the depletion force. . . . . . . . . . . . . . . . . . 18 2.1 Experimental images of flat disks and twisted ribbons and the phase diagram. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 24 2.2 Experimental images of flat disks utilizing different techniques. . . . . 25 2.3 Electron micrograph image of the edge of a flat disk. . . . . . . . . . 25 2.4 Comparison of the theoretical and experimental retardance of a flat semi-infinite disk. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 2.5 Solutions of the director field on finite disks with different radii. . . . 31 2.6 Optimal radius of the finite disks as a function of the bare line tension. 32 3.1 The local coordinate system on a patch of the monolayer. . . . . . . . 35 3.2 Visualization of the director field on a twisted ribbon. . . . . . . . . . 40 3.3 Tilt angle θ as a function of ρ on twisted ribbons with different widths and pitches. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 3.4 Phase diagrams exhibiting disks and twisted ribbons. . . . . . . . . . 43 3.5 Dependence of width and pitch of the twisted ribbons on γ with dif- ¯ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ferent k. 44 3.6 Edge-to-area ratio of TRs as a function of the line tension. . . . . . . 46 3.7 Maximum tilt angle at the edge of the TRs on the phase boundary be- tween semi-infinite disks and TRs as a function of k, ¯ with comparison to the maximum tilt angle at the edge of semi-infinite disks. . . . . . 46 3.8 Width and pitch at the phase boundary between semi-infinite disks and TRs as a function of k. ¯ . . . . . . . . . . . . . . . . . . . . . . . 47 3.9 Phase diagram for stronger chirality including semi-infinite flat disks, finite disks and twisted ribbons. . . . . . . . . . . . . . . . . . . . . . 48 3.10 Plot of different terms in the free energy density for TRs across the phase boundary between flat membranes and TRs. . . . . . . . . . . 49 4.1 Several snapshots from experiments during the phase transition pro- cess from disks to TRs. . . . . . . . . . . . . . . . . . . . . . . . . . . 54 4.2 Stretching of a flat membrane utilizing two optical traps to make a TR. 55 x 4.3 The force felt by the optical trap upon stretching a flat membrane to make a TR. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 4.4 Shape of a rippled disk. . . . . . . . . . . . . . . . . . . . . . . . . . . 58 4.5 Solutions of the director field on rippled disks and the corresponding birefringence images. . . . . . . . . . . . . . . . . . . . . . . . . . . . 61 4.6 Theoretical model for the shape of the stretched membrane. . . . . . 67 4.7 An example of a ribbon-like shape perturbation on a flat semi-infinite disk with a Gaussian-bumped edge. . . . . . . . . . . . . . . . . . . . 67 4.8 The variation of the different terms in the free energy perturbation as a function of the extent of stretching for a protrusion with shape of a Gaussian bump. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 4.9 The total variation in the free energy for the ribbon-like perturbation as a function of the extent of stretching with different choices of the mean curvature modulus and the line tension with the Gaussian bump protrusion. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 4.10 The variation of the different terms in the free energy perturbation as a function of the extent of stretching for a protrusion with shape of a equilateral triangle. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72 4.11 The total variation in the free energy for the ribbon-like perturbation as a function of the extent of stretching with different choices of the mean curvature modulus and the line tension with the equilateral triangle protrusion. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 72 5.1 Experimental images of the mixture membranes and a schematic of the cusps. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76 5.2 Fluorescence image of mixture membranes. . . . . . . . . . . . . . . . 77 5.3 An example of the model shape for the mixture membrane. . . . . . . 78 5.4 The phase separation of the theoretical model when the chirality is allowed to vary locally. . . . . . . . . . . . . . . . . . . . . . . . . . . 81 5.5 An example of solution of the director field on a mixture membrane and the corresponding birefringence image. . . . . . . . . . . . . . . . 82 5.6 Free energy per unit area of a ’cookie disk’ as a function of height of the cusps, size of the radial bulges and the spacing between the cusps. 83 5.7 Contribution from different terms in the free energy for the mixture membranes as functions of the height of the cusps and the size of the radial bulges. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85 5.8 Phase diagram of the mixture membrane and flat disks. . . . . . . . . 88 xi Chapter One Introduction 2 1.1 Nematic and Smectic Phases of Liquid Crystal The commonly known phases of matter include solid, liquid and gas. These phases are distinguished by their different physical properties such as symmetry groups, densities and chemical compositions. Among these physical properties, symmetry change is related to many phase transitions such as solid-liquid phase transitions. Besides these well-known phases, there still exist intermediate phases whose symme- try property is in between two of these common phases. Liquid crystalline phases are a class of these intermediate phases usually consist of highly anisotropic molecules which are rod-like or disk-like. Like crystals, liquid crystals have broken orientation- al symmetry and in certain phases such as in smectics, there is long range order in one of the spacial dimensions as well. Meanwhile, like liquids, liquid crystals have no long range translational order in at least one space direction. According to the details of the symmetry groups, liquid crystals can be further categorized into nematic, smectic-A (Sm-A), smectic-C (Sm-C ), discotic, cholesteric and so on. Nematic liquid crystals are characterized by orientational order of the constituent molecules in which at least one of the principal axes of the molecules are roughly oriented in the same direction and the absence of translational long range order in all three spacial dimensions. On the other hand, smectic phases not only exhibit long range orientational order but also have long range positional order in one spacial dimension. Specifically, the molecules are organized in a layered structure. The difference between Sm-A and Sm-C phases is that in the former the molecules in each layer are nearly parallel to the layer normal while in Sm-C phases the molecules are tilted to a preferred direction with a finite angle with respect to the layer normal. The cholesteric phase is characterized by the chiral property of the composing molecules which usually results in a helical structure in the orientational 3 ordering of the molecules. The period of the helical structure is defined as the pitch. For a complete review of liquid crystals, we refer the readers to Ref. [1, 2]. 1.2 Director Field, Order Parameter The director field is a concept that is introduced to describe the local preferred direction of the anisotropic molecules composing liquid crystals. It is a vector field with unit length and will be referred to as n through out this work. In general in the nematic and smectic phases the ordering is such that the heads and tails of the molecules are not distinguished and thus there is no difference between n and −n. The order parameter is a physical quantity introduced to measure the magnitude of order in a phase usually defined near a phase transition. It is a measure of the extent to which symmetry has been broken within the ordered phase. The order parameter for a uniaxial nematic phase is a second rank traceless tensor in the form[2] 1 Sαβ = S(nα nβ − δαβ ) (1.1) 3 where S is the magnitude of the order parameter tensor which measures the extent to which the molecules are aligned, nα is the α-th component of the director field and δαβ is the Kronecker delta. 4 1.3 Landau Theory for the Nematic-Isotropic Phase Transition In the Landau theory of phase transitions, the free energy is written as a power expansion of the order parameter with phenomenological coefficients near a phase transition[3]. The form of the free energy must agree with the symmetry the matter assumes. A general form of the Landau free energy density is fL = C1 φ2 − C2 φ3 + C4 φ4 + . . . + D1 (∇φ)2 (1.2) where φ is an order parameter. Landau theory has been successful in qualitatively describing phenomenon near second and some first order phase transitions. In the specific example of nematic-isotropic phase transition, the free energy must be a scalar function of the order parameter tensor. The general form of the free energy density can be written as[2] 1 1 1 1 f = ASαβ Sαβ − BSαβ Sβγ Sγα + CSαβ Sαβ Sγδ Sγδ − χa Hα Sαβ Hβ (1.3) 2 3 4 2 where χa = χk − χ⊥ is the magnetic susceptibility anisotropy and H is external magnetic field. Usually the three coefficients A, B and C are functions of temperature and pressure. Landau’s assumption is that B and C are constants independent of temperature while A = a(T − T ∗ ) where a is a constant independent of temperature. Later we shall see that T ∗ is related to the critical temperature where the phase transition occurs. In the uniaxial case without external magnetic field, inserting the form of the 5 Figure 1.1: Landau free energy density for a first order phase transition at different temperatures. Tc is the temperature where the phase transition occurs, T ∗ and T ∗∗ are the critical supercooling and superheating temperatures respectively. Reproduced from Principles of condensed matter physics, page 170. order parameter tensor (1.1) into (1.3), we obtain 1 2 1 f = AS 2 − BS 3 + CS 4 (1.4) 3 27 9 The equilibrium value of the scalar order parameter S should minimize (1.4). As a function of S, the free energy density is plotted in Fig. 1.1. From the figure, we can see that when the temperature is above T ∗∗ , there is only one minimum at S = 0, indicating an isotropic phase. After this temperature is reached from above, there are two minima of the curve at S = 0 and S = (B/4C)[1 + (1 − 24AC/B 2 )1/2 ] of which the latter corresponds to the nematic phase with the isotropic phase being the absolute minimum. This temperature T ∗∗ is the critical temperature a nematic phase can be superheated to. At the temperature 1 B2 Tc = T ∗ + (1.5) 27 aC the two minima have equal free energy density and below it, the nematic phase 6 becomes the equilibrium phase and the isotropic phase becomes metastable. After the temperature T ∗ is reached, the minimum at S = 0 vanishes and nematic phase is the only phase existing. So this temperature is the lowest temperature that the isotropic phase can be supercooled to. Note that the order parameter jumps from 0 to a finite value at the phase transition indicating that the phase transition between nematic and smectic phases is first order. 1.4 Elasticity Theory of Liquid Crystals It costs energy to distort the director field in a nematic phase. By symmetry argu- ment, the free energy should be even in n and proportional to (∇n)2 to the lowest order. When the distortion ∇n is small on the molecular scale, the free energy density of this distortion is given by[3] 1 1 1 fd = K1 (∇ · n)2 + K2 (n · ∇ × n)2 + K3 (n × ∇ × n)2 (1.6) 2 2 2 where Ki is the elastic constants. The three terms above correspond to three kinds of deformation in the director field, namely splay, twist and bend respectively. See Fig. 1.2 for an illustration of these deformations. When the composing molecules are chiral, which result in a cholesteric phase, the second term should be modified to 1 K2 (n · ∇ × n − q)2 (1.7) 2 where q defines the sign and strength of the chirality. It is easy to show that the pitch of the cholesteric phase is 2π/q. Eq. 1.6 is the basic formula for the elastic theory of liquid crystals and the free energy is also called the Frank free energy. 7 Figure 1.2: A demonstration of (a) splay, (b) bend and (c) twist deformations in a nematic liquid crystal. Reproduced from Ref. [3], page 300. In theoretical discussions of the distortion free energy of liquid crystals, one often takes the one-constant approximation where all three elastic constants are taken to be equal[1] K = K1 = K2 = K 3 (1.8) which results in 1 1 fd = K[(∇ · n)2 + (∇ × n)2 ] = K(∇n)2 (1.9) 2 2 Although many experimental data have shown significant differences in values of the three elastic constants, this approximation is widely used in theoretical work and has indeed produced many useful results. A well known example of the application of this elasticity theory is the widely studied Frederiks transition[4, 5, 6, 7, 8, 9]. In the geometry of the Frederiks transi- tion, a nematic liquid crystal cell is placed between two parallel boundary surfaces. These surfaces are rubbed so that the directors at the surfaces are anchored to be 8 parallel or perpendicular to them. Then an electric field or magnetic field is applied to the nematic cell which favors a direction of the director field that is incompatible with the anchoring direction. When the strength of the field is increased, a transition from an undistorted state to a distorted state characterized by a finite rotation of the directors at the center of the cell and the corresponding relaxation in the director field toward the anchoring state at the boundaries can be observed at a critical field strength. The Frederiks transition is a perfect example of the interplay between the elasticity and the response to an external field of a nematic liquid crystal. 1.5 Topological Defects in Liquid Crystals Topological defects are one type of defects characterized by a core region where the definition of order parameter is ambiguous and a far field region where a slow deformation in the order parameter is assumed. Topological defects exist in systems where usually a continuous symmetry is broken (kinks in Ising models are examples where a discrete symmetry is broken) and are given different names according to the type of the broken symmetry. Topological defects cannot be eliminated simply by a continuous transformation of the order parameter. The sign and strength of one topological defect can be defined by looking at the variation of the order parameter on a curve encompassing the core region. Topological defects with opposite sign and equal magnitude of strength can annihilate each other. The free energy to create a topological defect includes the core energy and the distortion energy in the order parameter. For the detail of topological defects, we refer the readers to Ref. [10]. Examples of topological defects in liquid crystals include integer disclinations at the surface of nematics, ±1/2 disclinations inside 3-d nematics and dislocations in 9 Figure 1.3: A sketch of the directors whose local direction is represented by the tangential direction of the lines near a 2d disclination with different strengths in nematics and the corresponding Schlieren textures. Reproduced from Current Applied Physics, 12, 1387 (2012). smectic phase. Shown in Fig. 1.3 are sketches of the director field near disclinations with different strengths[11]. Note that on a closed contour encircling the core of the disclinations, the director field rotates by ±2sπ where s is the strength of the disclination. Also shown in this figure are two Schlieren textures[1, 12, 13, 14] cor- responding to the director field near the disclinations. The physical reason of the formation of such textures is deferred to a later section. 10 1.6 Order Parameter and Elastic Theory of Smec- tic Phase The smectic liquid crystalline phase is characterized by its layered structure. The correlation function is liquid like inside each layer while it is crystal like in the direction perpendicular to the layer. Mathematically, the mass density wave is used to describe the layered structure[3] which is given by ρ(x) = ρ0 + [hψ1 ieiq0 ·x + c.c.] (1.10) where q0 = (2π/d)n0 is the wave vector related to the period d of the layered struc- ture and the abbreviation c.c. stands for complex conjugate. The higher harmonics have been neglected in this expression. hψ1 i is a complex order parameter which is 0 for the nematic phase and has a finite amplitude and finite phase angle for the smectic phase which is written as hψ1 i = |hψ1 i|e−iq0u (1.11) The position of the n-th layer is given by q0 · x − q0 u = 2πn (1.12) The distortion of the layers in a smectic phase costs elastic free energy. At first sight, the free energy should be proportional to (∇u)2 , which is similar to the nematic distortion free energy proportional to (∇n)2 . But by more careful argument, one kind of distortion that is rigid rotation of the layers and the director fields together should not cost any energy. So the lowest order in the elastic free energy density should 11 read 1 fs = B(∇k u)2 + D(∇⊥ u + δn)2 + fd (1.13) 2 where the first term is related to the compression or expansion of the layer spac- ing and fd is the Frank free energy given in Eq. (1.6). Among the three kinds of distortions, namely, splay, twist and bend deformations, the splay deformation cor- responds to the lower energy excitation in smectics because only splay deformation can result in an unchanged layer spacing while it is impossible to keep the layer spacing constant in a bend or twist deformation. Thus, bend and twist deformations are expelled from deep inside the smectic phase[15]. 1.7 Landau Theory for Nematic-Smectic-A Phase Transition and Its Similarity with Supercon- ductors It has been mentioned previously that the main order parameter that distinguishes the smectic-A phase from the nematic phase is the first harmonic of the mass density wave ψ1 . Keeping in mind that a rigid rotation of the smectic layers and the directors introduces a phase factor into the order parameter, we should replace the space derivative term ∇⊥ φ in Eq. (1.2) by (∇⊥ − iq0 δn)ψ1 while keeping the ∇k ψ1 term unchanged. Thus the Landau-Ginzburg free energy density for the nematic-smectic- A reads[3] 1 f = r|ψ1 |2 + ck |∇k ψ1 |2 + c⊥ |(∇⊥ − iq0 δn)ψ1 |2 + g|ψ1 |4 + fd (1.14) 2 where fd is again the Frank free energy density. 12 As pointed out by de-Gennes, the form of free energy density in Eq. (1.14) describing the nematic-smectic-A transition is similar to the Ginzburg-Landau free energy density describing the superconductor-normal phase transition[15, 16] which is given by  2 β 4 ~2  2 2eA + 1 (∇ × A)2 − 1 H0 · (∇ × A) fs = fn + α|ψ| + |ψ| + −i∇ − ψ 2 4m ~c 8π 4π (1.15) where fs and fn is the superconductor and normal state free energy density respec- tively, ψ is the wave function of the coherent Cooper pair wave function or in other words, the order parameter of the superconductor phase, A is the magnetic vector potential, H0 is the external magnetic field and cgs units are used here. Note that if we replace A in this formula by δn and make corresponding replacements to the coefficients, the formula Eq. (1.14) is recovered. It should be pointed out that there is no term corresponding to the splay deformation in the superconductor free energy and there is no difference between the elastic constant related to bend and twist deformation. The expulsion of twist and bend deformations from bulk smectic phase is the analog of the well-studied Meissner effect in superconductors[17]. The term contributed by the external magnetic field is similar to the contribution of chirality in the twist term in Eq. (1.7), so we can argue that the effect of chirality is similar to an external field. If the chirality is too strong, the smectic-A phase is destroyed similar to the destruction of superconductor order when the external magnetic field is beyond a certain critical value. 13 1.8 Coherence Length, Penetration Length and TGB Phase Seeing the similarity between the smectic-A phase and the superconductor phase, we can define two characteristic lengths in the smectic-A phase similar to the corre- sponding lengths in superconductors, namely, the coherence length and penetration length. The coherence length is the length scale over which the order parameter is correlated. Above the transition temperature Tc , it describes the size of the smectic domain inside the nematic phase while below Tc , it is the characteristic length scale over which a perturbation of order parameter relaxes. The size of a dislocation core inside the smectic-A phase is comparable to the coherence length. Mathematically, it is given by[16] q ξk, ⊥ = ck, ⊥ /|r| (1.16) It is well known that the London length in superconductors describes how deep the external magnetic field can penetrate into the bulk[18]. Similarly, two penetration lengths can be defined in smectic-A phase to describe the length scales on which a twist or bend deformation can penetrate into the bulk. The mathematical definition of the penetration lengths are given by[16] s k, ⊥ 1 ck, ⊥ K2,3 g λ2,3 = (1.17) q0 2|r| Depending on the ratio between the coherence length and penetration length κ = λ/ξ named the Ginzburg-Landau parameter, superconductors can be categorized √ into two kinds. When κ < 1/ 2, it is called type I superconductor and when √ κ > 1/ 2, it is type II. When the external magnetic field reaches a critical value 14 hc1 (T ) from below in a type II superconductor, vortex lines are formed inside the bulk which can form a vortex lattice[16]. Inside the core regions of these vortices, it is normal conducting phase while outside the core region far from the vortices, it relaxed back to the superconductor phase. This phase can only exist between two critical fields hc1 (T ) and hc2 (T ) and is named Abrikosov phase[19, 20, 21]. There is no such coexistence phase for type I superconductor. By the similarity between superconductors and smectic-A phases, it is expected √ that when κ > 1/ 2 and the counterpart of the external magnetic field, the chirality, is tuned between the lower and upper critical value, there should exist a phase where dislocations can coexist with smectic-A phase similar to the Abrisokov phase. This resulting twist grain boundary phase (TGB) firstly predicted by Renn and Lubensky[22] theoretically consists of smectic slabs separated by grain boundaries composed by arrays of screw dislocations. The TGB phase is a widely studied example of the interplay between smectic order and chirality[23, 24, 25, 26, 27]. 1.9 Birefringence of Liquid Crystals When the refractive index of a substance depends on the directions of the polarization and propagation of the light, the substance is said to be birefringent. Birefringence is commonly seen in anisotropic materials such as crystals and liquid crystals. The simplest example of birefringence is that inside uniaxial materials. The axis is called the optical axis. Generally speaking, when light is incident with an arbitrary angle into the material, the light will decompose into two rays whose polarizations are parallel and perpendicular to this optical axis respectively. These two rays are called extraordinary light and ordinary light respectively and the corresponding refractive 15 indices are denoted by ne and no . The strength of the birefringence is measured by the difference between these two indices ∆n = ne − no . Inside a uniaxial liquid crystal phase like uniaxial nematics, the local optical axis coincides with the director field and the o-wave is always polarized perpendicularly to the director. The well known method to observe the birefringence inside the liquid crystals is to put a slab of the substance between two crossed polarizers. When the linear polarized light enters the medium, it will decompose into extraordinary light and ordinary light with different refractive indices. These two waves acquire a phase difference upon exiting the slab and rejoin. Thus, the local polarization of the outgoing light from the slab depends on this local phase difference, or in other words, the local optical axis orientation and the thickness of the slab. After transmitting the second polarizer, the intensity of the light is no longer uniform in space and is determined by this local polarization. In a simplified theoretical calculation, the slab and the two polarizers are put parallel to the x-y plane and the light is incident from the z-axis. It is assumed that the local directors lie in the x-y plane and their directions do not vary with the z-coordinate. The outgoing light intensity is shown to be[16]   2 2 πd I = I0 sin (2β) sin (ne − no ) (1.18) λ0 where I0 is the incident light intensity, β(x, y) is the angle made by the director to the polarization of the incident light, d is the thickness of the slab and λ0 is the wave- length in the vacuum of the incident light. Thus, the intensity of the outgoing light depends on the local directions of the director field and the famous Schlieren texture is formed. The Schlieren texture is a powerful tool to visualize the dislocations inside nematics. 16 1.10 Excluded Volume Effect and Depletion Force We all know from thermodynamics that the free energy function F = E − T S controls the phase transition behavior where E is internal energy and S is entropy. In many solution systems with added macromolecules, the added particles floating around interact with each other only when they are close enough. In theoretical models, this kind of interaction is usually modeled by hard core repulsions. Thus, the internal energy E is approximately constant in these systems and the factor that dominates the free energy at constant temperature is the entropy. This fact implies that many entropy-driven effects should determine the behavior of these systems. Excluded volume effects and depletion force are two examples of the entropy effects in colloidal systems. First let us consider a solution with spherical colloidal particles with radius r interacting by hard core repulsion. Imagine the colloidal particles are added one by one. When the first colloidal particle is present, a finite volume is occupied where the second particle cannot fit. The possible volume that the center of mass of the second particle can reside in is thus reduced by an amount v1 = 4π(2r)3/3, so that the positional entropy of the second particle is reduced. Likewise, when the third particle is added, the number of possible positions is further reduced by the combination of the excluded volumes of the two previously added particles. If the excluded volumes of the previously two added particles overlap, then the total excluded volume is reduced and the entropy to add the third particle is increased and is favored. This makes the state in which the colloidal particles closely packed favored from the entropy perspective. When the shape of the colloidal particles is rod-like instead of sphere, the excluded 17 Figure 1.4: An illustration of excluded volume effects for hard spheres on the left (Reproduced from http://www.zeably.com/Excluded volume) and hard rods on the right (Reproduced from Bax- ter et. al. Composite Sci. Tech. 71, 1273 (2011)). volume becomes orientation-dependent. It is obvious that the excluded volume is reduced and thus the entropy is increased when the rod-like molecules align with each other. Onsager utilized this idea in his theory of nematic order[28]. Depletion force is another kind of entropy effect whose basic idea is demonstrated in Fig. 1.5. When there are not only colloidal particles but also some smaller particles called depletant present in the solution, a force that tends to assemble the large particles appears[29]. Usually the depletant unit is macromolecules with flexible tails. By the accumulation of the colloidal particles, the excluded volume is reduced and the positional entropy of the depletant is increased. More importantly, when a larger volume can be occupied by these macromolecules, their flexible tails can assume much more configurations which enhances the configuration entropy by a large amount[30]. This explains the appearance of the depletion force. Depletion force is often utilized to make self-assemblies[31, 32, 33, 34, 35]. The systems considered in this work are also formed with the presence of depletion force. 18 Figure 1.5: A demonstration of the depletion force. Reproduced from http://soft- matter.seas.harvard.edu/index.php/Repulsion - Steric(entropic). 19 1.11 An Introduction to Membrane Systems Self assembled membranes have attracted much attention both experimentally and theoretically in physical, chemical and biological studies. For example, in nanotech- nology, the membranes self assembled on interfaces can dramatically change the inter- facial properties[36]. In biology, certain kind of peptide can spontaneously form mem- branes when salt is added to the solution which might find applications in the field of biomaterials[37]. Among these membrane systems, a lipid bilayer is one that has been widely studied due to its novel physical properties, its ability to bear proteins or other structures on it and its close relation to biological applications[38, 39, 40]. The constituent molecules of lipid membranes are amphiphilic with different head and tail structures. These membranes can form closed structures like vesicles with the size on the micrometer scale when placed in aqueous solution which is an important structure in many biological processes[41, 42]. The hydrophilic heads of the com- posing molecules in a vesicle are in contact with the aqueous solution outside while the enclosed region can hold different substances. Surfactant membranes dividing the oil and water region in microemulsions are another kind of interesting structure that can form many complicated 3-d topologies[43]. There are different theoretical models and numerical studies constructed to ex- plore the behavior of the membranes and the structures attached to these membranes[46]. For example, fluid mosaic model are proposed to model the interplay between a lipid membrane and the proteins that are embedded in it[47]. On the micrometer scale, elasticity model for membranes described as smooth surfaces in the language of dif- ferential geometry have been formulated[69, 75, 48, 49, 50]. We should mention here the most studied Helfrich model which will be central to this work. In the Helfrich 20 model, the free energy of a membrane is related to its curvatures and is written as[69] 1 1 Z Z Z Z FH = k (c1 +c2 −c0 ) dA+ k¯ 2 c1 c2 dA = k (2H −c0 ) dA+ k¯ 2 KG dA (1.19) 2 2 where c1 , c2 are the two eigenvalues of the curvature tensor, c0 is the spontaneous curvature, H = (c1 + c2 )/2 is the mean curvature and KG = c1 c2 is the Gaussian curvature, k is the mean curvature modulus or the bending rigidity and k¯ is the Gaussian curvature modulus. The Helfrich model has been largely successful to explain phenomena related to vesicles[51, 52, 53]. When there is certain kind of field such as director field present on membranes, the behavior of the system can be more complicated. Membranes with nematic order on them are examples of the interplay between the field and the underlying geometry[54, 55, 56]. 1.12 Monte Carlo Simulation Monte Carlo (MC) simulation is a numerical method first introduced after computers were invented and found its first application in the invention of atom bomb. The systems solved by MC simulation usually are not fully deterministic in microscopic scale, or in other words, stochastic in some perspectives. These stochastic processes can be realized numerically by generating a sequence of random numbers. Although different sequences of the random number usually result in different microscopic final states of the system, the physical properties, such as pressure, order parameter and chemical composition of the system usually fall in a region within statistical error. Thus much information like the equilibrium properties, the fluctuations and so on can be retrieved for the system[57]. 21 There are many different algorithms to decide whether a MC step is accepted or discarded during the simulation. The most common used criterion is the Metropolis method. The basic process is as this: the change of the energy ∆U resulted from the test move is computed; if ∆U < 0, i.e. the energy of the system is lowered, then the step is accepted; otherwise a random number r which is uniformly distributed in [0, 1) is generated, and if r < exp(−β∆U) where β = 1/kB T is the inverse temperature, then the step is accepted; otherwise the test move is discarded and the system remains unchanged. MC simulation has been more and more widely used in almost all fields of scientif- ic research with the computers becoming more and more powerful. It has been largely successful in the field of physics for solving problems beyond analytical reach. Well known examples of the application of MC simulation include percolation[58, 59, 60], diffusion limited aggregation[61, 62], properties of magnetic systems[63], and so on. Chapter Two Basic Experimental Methods and Our Core Theoretical Model 23 2.1 Experimental Methods and Facts The system that this work is concentrated on is composed of fd viruses. A single fd virus assumes a 880 nm long filamentous shape with diameter of 6.6 nm[64] whose persistence is several times longer than the macromolecule itself[65]. The chirali- ty of fd viruses can be tuned by varying the temperature[33]. They are chiral at room temperature and the chirality is weakened upon heating until total disappear- ance at about 60◦ C. The interaction between the molecules is purely repulsive in aqueous solution. By increasing the concentration, they can form isotropic, nematic (cholesteric) and smectic phases sequentially[66]. Upon adding polymers into the solution, the interaction becomes attractive under the mechanism of depletion force and the molecules can self-assemble into colloidal phases and sediment to the bot- tom due to the larger molecular weights[67]. Among the various phases observed, the smectic-A single layer membranes whose thicknesses are about one molecular length have attracted much attention. The viruses are almost perfectly straight in these membranes so that they will be modeled as hard rods in this work. The smectic-A monolayer membranes self-assembled assume different shapes like round flat disks, twisted ribbons and disks with π-walls when the concentration of Dextran which provides the depletion force is varied. Examples of the experimental images of flat disks and twisted ribbons are shown in Fig. 2.1. Specifically, when the concentration of Dextran is lowered from above, a first order phase transition from round disks to twisted ribbons is observed. To further understand the structure in the interior of these membranes, several techniques are used especially on the round disks which are the easiest to make quan- titative observations on. For example, face-on phase contrast images are obtained to distinguish the regions with different refractive indices, fluorescence images can 24 Figure 2.1: Experimental images of flat disks (left) and twisted ribbons (center) and the phase diagram (right). As for the phase diagram, note that on y-axis, when temperature is lowered, the strength of the chirality is increased while on the x-axis, the concentration of the Dextran corresponds to the line tension in the theoretical model. Phase diagram reproduced from Ref. [33]. show the location, direction and movement of several single viruses in the assembly, retardance images calculated from LC-PolScope images show the slow axis orienta- tions at different spots. These images contain much information about the director field which are important guides to the theoretical research. Generally, these images show that the directors stay parallel to the layer normal near the center of the disks and lie down tangential to the edge at the periphery. We show several images in Fig. 2.2. Analysis of the images concludes that the penetration depth of the twist deformation mentioned in Sec. 1.8 is about 0.5 µm[66]. The effective interaction of the membranes with Dextran is reflected in the quantity named line tension γef f and is measured by analyzing the fluctuations of the edge of a flat disk. Besides these coarse grained techniques, individual viruses can be visualized in an electron micrograph as shown in Fig. 2.3 and the shape of the edge of a flat disk is observed. 2.2 Core Theoretical Model In order to understand the behavior of the monolayer smectic-A membranes, the free energy needs to be constructed for the monolayer. Theoretically, we consider a 25 Figure 2.2: Experimental images reproduced from Ref. [66]. (a) Phase contrast image. (b) Fluorescence image. One labeled molecule (pointed by an arrow) at the edge is seen lying down while the molecules in the interior stay perpendicular to the layer. (c) Retardance image. (d) Orientation of the slow axis. All scale bars are 5 µm. Figure 2.3: Electron micrograph image visualizing individual viruses at the edge of a flat disk. The scale bar at the bottom is 100 nm. Note that the molecules are tangential to the edge at the periphery and are parallel to the layer normal at the center. Reproduced with permission by Dogic and coworkers. 26 free energy density constituted by several parts[68]. The first part is the curvature elastic energy given by Helfrich[69, 70] related to the deformation of the monolayer from a flat one and is purely geometrical. The Helfrich free energy density is given by 1 ¯ G fH = k(2H − c)2 + kK (2.1) 2 where H and KG are the mean and Gaussian curvatures of the monolayer, k is the bending rigidity and k¯ is the Gaussian curvature modulus and c is the spontaneous curvature reflecting the up-down asymmetry of the membrane. Given the up-down symmetric nature of the fd virus, c is always taken to be 0 in this work. The second part involves the smectic order and the elastic energy of deformation of director field which is written in Eq. (1.14) with the chiral term in the Frank free energy density Eq. (1.7). For clarity, we rewrite this de Gennes free energy density explicitly here 1 fN A = r|ψ1 |2 + ck |∇k ψ1 |2 + c⊥ |(∇⊥ − iq0 δn)ψ1 |2 + g|ψ1 |4 2 (2.2) 1 1 1 + K1 (∇ · n)2 + K2 (n · ∇ × n − q)2 + K3 (n × ∇ × n)2 2 2 2 There are several assumptions made to simplify this free energy density stated as following. Firstly, we assume a uniform and non-zero smectic order parameter ψ1 through out the membrane. Though this is not exactly correct based on the ex- perimental facts that the viruses lie down at the edge so that the smectic order is transformed to nematic order there, this approximation reduces the complexity of the model by a large amount and still provides good numerical results shown in the following chapters in this work. Under this assumption, the first and second term in Eq. (2.2) are simply constant and can be absorbed into the definition of the free energy density, and the third term becomes a pure tilt energy density describing 27 the energy cost resulted from the deviation of the director from the layer normal. Secondly, the one constant approximation is taken. Thus f becomes 1 1 fN A = K[(∇ · n)2 − 2qn · (∇ × n) + (∇ × n)2 + q 2 ] + C sin2 θ (2.3) 2 2 where K = K1 = K2 = K3 is the single Frank constant, C is a constant proportional to c⊥ in Eq. (2.2) and θ is the tilt angle of the director to the layer normal. We also incorporate a third term which effectively models the interaction of the membrane with the surrounding environment including the depletion force. For simplicity, we take this term to be proportional to the length of the edge of the monolayer as following I fedge = γ dl (2.4) where γ is the theoretical bare line tension. The experimental line tension is given by the sum of this quantity and the free energy gain by tilting the directors at the edge[33]. The total free energy of the monolayer membrane is thus the summation of the three terms as mentioned above integrated over the area of the membrane Z I F = (fH + fN A )dA + γ dl (2.5) Eq. (2.5), (2.3) and (2.1) are the core formulae in our theoretical model for the free energy of the smectic-A monolayer membrane. 28 2.3 Application of the Model to Flat Disks We should postpone the author’s work to the following chapters and for the com- pleteness of this work, review the previous theoretical results in Ref. [66, 71] where the theoretical model is applied to flat disks. A flat disk is the simplest structure of the smectic-A membranes observed in the experiments. The Helfrich free energy vanishes in this case. First, we consider a flat semi-infinite disk filling the x ≥ 0 half x-y plane. By symmetry, it is seen that the orientation of the director is purely in the y-z plane. The de Gennes free energy per unit length reduces to[66] ∞   1 1 Z 2 2 F = K(dθ/dx − q) + C sin θ dx (2.6) 0 2 2 Using a variational method, we derive the Euler-Lagrange (EL) equation of the director field to be d2 θ 1 λ2t − sin 2θ = 0 (2.7) dx2 2 where λt = (K/C)1/2 is the penetration depth. The boundary conditions(BC) for this differential equation are as following. Far from the edge, i.e. at x = ∞, the directors are parallel to the layer normal meaning θ = 0. At the edge, according to the experimental observations, the viruses are lying down into the plane of the membranes tangential to the edge suggesting a fixed BC with a 90 degree tilt at the edge. But purely from a theoretical point of view, at least in the model proposed by far, there is no argument requiring such a fixed BC and the free BC is also a reasonable choice. Actually for the previous works reviewed here, only the free BC at the edge is chosen and good results are reproduced compared to experiments. The free BC at the edge requires that the directors there are free to relax so that the 29 torque on the directors should vanish, meaning ∂f =0 (2.8) ∂(dθ/dx) which results in dθ λt = ± sin θ (2.9) dx The solution to the EL equation is θ0 x θ(x) = 2 arctan[tan exp(− )] (2.10) 2 λt with θ0 = − arcsin(qλt ) the maximum tilt angle at the edge. The tilt angle assumes an exponential decay into the interior of the disk as expected. The theoretical retardance at one particular point is proportional to the product of sin2 θ and the thickness of the membrane at that point. The thickness of the monolayer is taken to be constant. This is not true at the edge where the thickness reduces as the directors are tilting over, however the assumption produces good fit of the theoretical retardance to the experimental data as shown in Fig. 2.4. Next, we consider the more realistic flat disks with finite radius R. In this case, the tilt is purely in the polar angle direction and the angle is only a function of the radial coordinate r. The free energy reads[71] " 2 # R 1 sin4 θ 1  1 sin 2θ dθ Z F = K + −q + K 2 + C sin2 θ 2πrdr + 2πγR (2.11) 0 2 2r dr 2 r 2 The corresponding EL equation is given by d2 θ dθ 1 r2 r2 + r − 2qr sin 2 θ − sin 2θ cos 2θ − sin 2 θ sin 2θ − sin 2θ = 0 (2.12) dr 2 dr 2 2λ2t 30 Figure 2.4: Comparison of the theoretical (dashed line) retardance of a semi-infinite flat disk which is modeled by Eq. 2.6 and experimental (points) retardance of a disk with large radius as a function of the distance from the edge. The theoretical values have been convolved with the resolution function so that the peak of the retardance is shifted to a negative position on x-axis. Reproduced from Ref. [66]. and the boundary conditions are θ = 0, r = 0 (2.13) dθ sin 2θ + − q = 0, r = R dr 2r The solution to this differential equation is no longer analytically achievable. Thus a differential equation solver in MATLAB is used to obtain the numerical solution. From now on, we use dimensionless units where length is measured in the unit of penetration depth and energy is measured in the unit of the Frank constant. Several solutions with different radii of the disk with qλt = 0.71, a value derived from experiments, are shown in Fig. 2.5. Also shown in this figure are the solutions to the semi-infinite disks in solid lines. It is found that when the radius of the finite disk R ≥ 5, the tilt angle profile is well fitted with the semi-infinite profile. We will use this as a criterion to distinguish a semi-infinite disk to a finite disk. 31 Figure 2.5: Solutions of the tilt angle (points) on finite disks with radii (a) R = 2, (b) R = 3, (c) R = 4 and (d) R = 5. The solid lines are the solutions on semi-infinite disks and the vertical dashed lines indicate the center of the finite disks. Note that when the radius reaches R = 5, the tilt profile of a semi-infinite disk is almost indistinguishable from the profile of the finite disk. Reproduced from Liq. Cryst., 36, 1157(2009). In order to obtain the equilibrium radius of the disks at a given bare line tension, the free energy per unit area is calculated with different radii and then the minimum is looked for. It is found that when γ is increased, the optimal radius increases until a critical value γc is reached where the equilibrium radius becomes infinity as shown in Fig. 2.6. At this point, a second order phase transition from finite disks to semi-infinite disks occurs. 32 Figure 2.6: Optimal radius of the finite disks as a function of the bare line tension γ. The vertical dashed line indicates the location where the optimal radius becomes infinity. Reproduced from Liq. Cryst., 36, 1157(2009). Chapter Three Application of the Model to Twisted Ribbons 34 3.1 Derivation of the Free Energy for Membranes with Arbitrary Geometry The major results of this chapter are published in Ref. [68]. Let us consider the terms related to director field in the de Gennes free energy density in Eq. (2.3). These terms are easily calculated for the disks as shown in the previous chapter due to the simple geometry of the membranes, but their calculation becomes extremely tedious when the geometry becomes more complicated. A general calculable form of these terms is needed for arbitrary geometry of the monolayers. Differential geometry is the tool adopted by us to achieve this goal. We refer the readers to Ref. [72, 73, 74, 75, 76] for details of the underlying mathematical concepts which are also briefly outlined below. In the framework of differential geometry, a membrane embedded in Euclidean space can be modeled by a position vector with two coordinates as follows Y = Y(u1, u2 ) (3.1) The following geometrical quantities are needed in our calculation Yi = ∂i Y, Yij = ∂i ∂j Y = Γkij Yk + Lij N, ∂i N = ∂N Yi = −Lij g jk Yk , gij = Yi · Yj , g ij = (gij )−1 , g = det gij , (3.2) Lij = Yij · N, Lij = (Lij )−1 , L = det Lij , Γkij = g km Yij · N 35 Figure 3.1: The local orthogonal coordinate system formed by Y1 , Y2 and the local normal N on a patch of the surface. The director n makes a finite tilt angle θ with respect to the normal. where indices i, j, k = 1, 2 represent the two coordinates, ∂i = ∂ui is the partial derivative with respect to the i-th coordinate, ∂N is the partial derivative in the direction of the layer normal N, gij , Lij and Γkij are the first fundamental form, second fundamental form and the Christoffel symbol of the surface respectively; repeated indices are summed over unless otherwise indicated in the following of this work. The lower and upper indices of g and L indicate the covariant and contravariant components respectively. Note that the area of the local surface element is given by √ dA = gdu1 du2 . The local surface normal with unit length is given by Y1 × Y2 N= √ (3.3) g The set of the vectors {Y1 , Y2 , N} forms a local right-hand orthogonal coordinate system in which the director can be written as n = ni Yi + cos θN (3.4) Fig. 3.1 shows the local coordinate system on a patch of the surface. The fact that the director is a unit vector imposes a constraint on the director n · n − 1 = gij ni nj − sin2 θ = 0 (3.5) 36 and this should be incorporated into the free energy as a Lagrange multiplier. The gradient operator on the surface is given by ∇ = g ij Yi ∂j + N∂N (3.6) Applying this operator to the director field and after some calculation, we obtain the divergence and curl of the director field as follows ∇ · n = ni,i + Γjij ni − cos θLij g ji , ǫ3ji (3.7) ∇ × n = √ [gik (nk,j + Γkjl nl )N + (θ,j sin θ − 2Ljk nk )Yi ] g where the comma in the subscripts stands for a derivative with respect to the indices following it, for example, ni,j = ∂j ni and ǫijk is the totally antisymmetric Levi-Civita tensor. The twist term is calculated to be ǫ3ji n · (∇ × n) = √ [gil nl θ,j sin θ − 2Ljk gil nk nl + cos θgik (nk,j + Γkjl nl )] (3.8) g In the form of these geometrical quantities, the mean curvature and Gaussian cur- vature are given respectively by 1 L H = g ij Lji , KG = (3.9) 2 g Collecting the previous calculation results, we arrive at the final form of the free 37 energy density for a membrane with arbitrary shape apart from the edge energy 1 L K fM = k(g ij Lij )2 + k¯ + [gik gip nk,j np,j − gik gjp nk,j np,i + 2(gik gip Γpjl − gik gjp Γpil )nl nk,j 2 g 2g + (giq gip Γqjl Γpjk − giq gjp Γqjl Γpik + 4gii Ljk Ljl − 4gij Ljk Lil )nk nl + gii θ,j θ,j sin2 θ − gij θ,i θ,j sin2 θ + 4(gij Ljk − gjj Lik )nk θ,i sin θ] ǫ3ji C − Kq √ [gil nl θ,j sin θ − 2Ljk gil nk nl + cos θgik (nk,j + Γkjl nl )] + sin2 θ g 2 (3.10) We point out here that the twist term proportional to q here is identical to the Helfrich-Prost term[77] first introduced to model chirality in lipid membranes. The total free energy of the membrane is given by √ Z I 2 F = (fM + λ(gij ni nj − sin θ)) gdu1 du2 + γ dl (3.11) where λ is a Lagrange multiplier. The formulae (3.10) and (3.11) are the basic formulae of the work that follows. 3.2 Derivation of Euler-Lagrange Equations for Twisted Ribbons Mathematically, a twisted ribbon (TR) is modeled as a surface parameterized by two coordinates, the first of which is the radial distance from the axis ρ and the second is the revolution angle φ. The surface is written as Y = {ρ cos φ, ρ sin φ, bφ} (3.12) 38 in which |ρ| ≤ R where R is the half-width of the ribbon, 0 < φ < 2nπ where n is the winding number and pr = 2π|b| is the pitch. The sign of b determines the handedness of the ribbon; when b > 0, the ribbon is right-handed while b < 0 represents a left-handed ribbon. We focus on a right-handed ribbon without loss of generality. Defining (ρ, φ) as (u1 , u2), the tangential and normal vectors of the surface are given by Y1 = ∂ρ Y = {cos φ, sin φ, 0}, Y2 =∂φ Y = {−ρ sin φ, ρ cos φ, b}, (3.13) 1 N = √ {b sin φ, −b cos φ, ρ} g The corresponding geometrical quantities defined in Eq. (3.2) are calculated to be g11 = 1, g22 = ρ2 + b2 , g12 = g21 = 0, g = ρ2 + b2 , −b b2 L11 =L22 = 0, L12 = L21 = p , L=− , (3.14) ρ2 + b2 ρ2 + b2 ρ Γ212 = Γ221 = , Γ122 = −ρ, 0 otherwise ρ2 + b2 The Gaussian curvature of a twisted ribbon is b2 KG = − (3.15) (ρ2 + b2 )2 and the mean curvature is 0 indicating that the twisted ribbon is a minimal surface. Note that the Gaussian curvature of a twisted ribbon is always negative. With these results, we can calculate from Eq. (3.11) the corresponding EL equa- 39 tions for a twisted ribbon to be 4Kb2 (ρ2 + b2 )−3/2 n1 + 4Kqb(ρ2 + b2 )−1/2 n1 + 2λ(ρ2 + b2 )1/2 n1 = 0, K[−(ρ2 + b2 )3/2 n2,11 −3ρ(ρ2 + b2 )1/2 n2,1 + 2(ρ2 + b2 )1/2 n2 − 2ρ2 (ρ2 + b2 )−1/2 n2 +2bθ,1 sin θ] + Kq[−4b(ρ2 + b2 )1/2 n2 − 2(ρ2 + b2 )θ,1 sin θ] + 2λ(ρ2 + b2 )3/2 n2 = 0, K[−(ρ2 + b2 )1/2 θ,11 sin2 θ − (ρ2 + b2 )1/2 θ,1 2 sin θ cos θ − ρ(ρ2 + b2 )−1/2 θ,1 ) sin2 θ −2bn2,1 sin θ] + Kq[2(ρ2 + b2 )n2,1 sin θ + 4ρn2 sin θ] − 2λ(ρ2 + b2 )1/2 sin θ cos θ + C(ρ2 + b2 )1/2 sin θ cos θ = 0, n21 + (ρ2 + b2 )n22 − sin2 θ = 0 (3.16) It is easily seen from symmetry that the directors do not tilt along the direction of the radial direction Y1 , indicating n1 = 0. Plug this into Eq. (3.16) and after some mathematical effort, the EL equation of the tilt angle θ is shown to be ρ2 + 3b2 2bρ K[−(ρ2 + b2 )θ,11 − ρθ,1 + 2 2 sin θ cos θ + 2 sin2 θ] ρ +b ρ + b2 (3.17) 2 2 2 +Kq[2ρ sin θ − 4b sin θ cos θ] + C(ρ + b ) sin θ cos θ = 0 and the form of the free energy is reduced to sin2 θ    1 ρ 4b 2 Z 2 F = K(θ,1 − q) + K sin 2θ + sin θ (θ,1 − q) + K(g + 3b2 ) 2 2 g g g b2 √ √  Z +C sin2 θ − 2k¯ 2 gdρdφ + γ R2 + b2 dφ g (3.18) The corresponding boundary condition for Eq. (3.17) at ρ = 0 is θ = 0 by symmetry. The boundary condition at the edge ρ = R is the free boundary condition given by 40 Figure 3.2: Visualization of the director field on a twisted ribbon. Note that the arrows are just guide to the eyes because directors preserve head-tail symmetry. ∂f /∂θ,1 = 0 which results in 2b R θ,1 (R) + sin2 θ(R) + 2 sin θ(R) cos θ(R) − q = 0 (3.19) R2 +b 2 R + b2 This latter boundary condition can also be derived by the physical argument that the torque should vanish at the edge or equivalently the twist of the directors is compatible with the spontaneous twist q, meaning n · (∇ × n) = q. Up to this point, we have obtained the EL equation and the corresponding boundary conditions needed to solve for the director field on a twisted ribbon. The solutions of this equation with different parameters are presented in the next section and further discussion made. 3.3 Results on Twisted Ribbons The differential equation for the tilt angle θ Eq. (3.17) is solved with a differential equation solver in MATLAB for TRs with varying parameters R and b. We visualize the solution of the director field on a twisted ribbon in Fig. 3.2. It is seen that the directors stay parallel to the layer normal on the axis of the TR while they tilt 41 0.8 0.8 (a) (b) 0.7 0.7 0.6 0.6 0.5 0.5 0.4 θ 0.4 θ 0.3 0.3 0.2 0.2 0.1 0.1 0 0 0 0.5 1 1.5 2 2.5 3 3.5 0 0.5 1 1.5 2 2.5 3 3.5 ρ ρ Figure 3.3: Tilt angle θ as a function of ρ on twisted ribbons with different widths and pitches. (a) Solutions for twisted ribbons with R = 3.66 and b =1.85(red), 2.35(green), 2.85(blue), 3.35(gray), 3.85(orange) from bottom to top. (b) Solutions for twisted ribbons with b = 3.85 and R =3.66(red), 3.16(yellow), 2.66(gray), 2.16(blue) and 1.66(green) from right to left. by a finite angle θmax from the layer normal to the tangential direction of the edge at the boundary. The tilt directions at opposite edges are opposite. Note that the handedness in the director field agrees with the handedness of the TR. In Fig. 3.3, we show the tilt angle θ as a function of ρ across half-width of TRs with different sets of the widths and pitches. In Fig. 3.3(a) the width R is fixed at 3.66 and the pitch b is varied. The solutions with different pitches almost coincide with each other only with a slight increase in θmax with increasing b. In Fig. 3.3(b) b is fixed at 3.85 and the width is varied. Generally, θ assumes an exponential-decay-like shape in the range of the parameters selected but it is seen from Fig. 3.3(b) that when the width becomes small, the tilt angle profile becomes more linear-like. This is because when the width is comparable to the penetration depth, the directors cannot relax as in the case with a larger width and the linear profile becomes the best choice to lower the free energy density. It is also observed that θmax is non-monotonic with R: it increases first until approximately R = 3.16 and then falls. Having obtained the solution to the director field, the total free energy of a TR with one pitch length can be calculated from Eq. (3.18). Dividing it by the area of 42 the TR calculated through  s  R  2 √ R R R Z Z p Ar = gdρdφ = 2π ρ2 + b2 dρ = 2πb2  1+ + arcsin  −R b b b (3.20) we arrive at the free energy per unit area fT R of a TR, the minimum of which is searched for by varying R and b at fixed k¯ and γ. Since the TR is a structure with ¯ 2 negative Gaussian curvature, a second order Gaussian curvature free energy kK G is added to the total free energy to stabilize the theoretical model[78, 79, 90, 91]. We choose k¯ to be as small as 0.01 in our dimensionless unit which is enough to eliminate the possible divergence in the free energy and meanwhile maintain the results obtained near the phase boundary between flat disks and TRs, where this term can be neglected. The values of R and b minimizing the free energy per unit area correspond physically to the equilibrium size of the TR when it is the stable phase where fT R is negative and lower than the free energy density for disks on the phase diagram spanned by k¯ and γ. When the Gaussian curvature modulus k¯ is chosen to be negative, it is found that the optimal pitch of the TR is several tens times of the penetration depth which is too large compared to the experimental results. So we focus on the positive k¯ side of the phase diagram. This is a natural choice because a positive k¯ favors structures with negative Gaussian curvature such as TRs. We should point out here that typically measured values of the Gaussian cur- vature modulus are negative for lipid monolayer or bilayer membranes, leading to the formation of structures with positive Gaussian curvature like the widely studied vesicles[80]. However, the fd system may be different from those lipid systems usu- ally composed of amphiphilic molecules. Given the fact that never has a vesicle-like shape been observed in such fd systems and that there are still more structures (we 43 Semi-infinite Disks 0.28 0.26 γ Twisted Ribbons 0.24 0.22 0 0.05 0.1 _ 0.15 0.2 0.25 k Figure 3.4: Phase diagram spanned by k¯ and γ including the semi-infinite disks and twisted ribbons. The red squares denote the phase boundary between semi-infinite disks and TRs and the line connecting them are guide to eyes. The dashed line at γ = 0.254 denotes the phase transition from semi-infinite disks to finite disks which is preempted by the transition from semi-infinite disks to TRs. will see examples in later chapters) seeming to exhibit negative Gaussian curvatures appearing in experiments, although k¯ cannot been directly measured at present, we believe that this positive k¯ conclusion is very likely true and we will give a brief explanation for this at the end of this section. The phase diagram containing semi-infinite flat disks, finite disks and TRs is one of the main results of our theoretical model which is presented in Fig. 3.4. The red squares denote the phase boundary between semi-infinite disks at the top and TRs at the bottom. The horizontal dashed line labels the phase transition from semi-infinite disks to finite disks at γ = 0.254. Our computation shows that the free energy per unit area for TRs is always lower than that of finite disks below the dashed line so that this latter transition is preempted by the phase transition from semi- infinite disks to twisted ribbons. This result suggests that in experiments when the 44 4.5 20 (a) (b) 4 15 3.5 3 10 2.5 R b 2 5 1.5 1 0 0.5 0 -5 0.05 0.1 0.15 0.2 0.25 0.05 0.1 0.15 0.2 0.25 γ γ Figure 3.5: Dependence of (a)R and (b)b on γ with k¯ = 0.05(red squares), k¯ = 0.10(green triangles) and k¯ = 0.15(blue circles). The lines are guides to eyes. line tension is lowered from above, large disks should transform into TRs rather than breaking into small disks, which agrees with the experimental observations. The TRs are still metastable near the phase boundary in the disk region and vice versa which indicates that the phase transition is of first order. By looking at the optimal R and b values, we speculate that the value of k¯ should be around 0.15 where R = 3.67 and b = 3.92 to have a good agreement with the experimental data. When converting to physical units, this indicates that k¯ is of one order of magnitude lower than the Frank constant which is a reasonable value. However, the critical line tension γc predicted by our theory is one order of magnitude lower than the experimental value which is approximately 100kT /µm. We attribute this discrepancy to the over-simplified edge model where all the effects of the appearance of an edge are included in a constant modulus. A more realistic edge model is beyond the scope of this work[82]. It is interesting to further look at the behavior of the TRs where they are the stable structure. In Fig. 3.5 we show the dependence of (a)R and (b)b on γ at k¯ = 0.05, 0.10 and 0.15 represented by red squares, green triangles and blue circles respectively. The width only depends weakly on k¯ near the phase boundary and is a decreasing function of k¯ when γ is small. It is seen from Fig. 3.5(a) that 45 when the line tension is reduced, the width of the TRs also reduces monotonically with the rate of decrease being faster at the beginning near the phase boundary. Different features are seen when k¯ is varied. We observe an abrupt drop of the width for k¯ = 0.10 and 0.15 at certain γ values which might be the indication of another transition where the TR suddenly shrank, and the value for this transition to happen shifts to higher γ value with k¯ increasing. As compared with the width, the pitch depends much more strongly on k¯ at all γ values. The higher k¯ results in TRs with shorter pitch because they have a more negative Gaussian curvature. When the line tension is reduced, we observe a non-monotonic behavior for k¯ = 0.05 and 0.10 and a monotonic decrease of the pitch for k¯ = 0.15. Experimental data have shown that the pitch is monotonically decreasing with a decrease of the line tension which further verifies the choice of k¯ = 0.15 as a good estimate in our theoretical model. As with the width, the pitch also drops suddenly for k¯ = 0.10 at the same value of γ when the width drops which supports our assumption that this is another transition where TRs suddenly shrink. For k¯ = 0.15, the pitch value is already small, which may make this transition invisible but it may be observable if the width is measured. We should point out that this sudden shrinking of TRs is not observed in experiments either because the line tension for this to happen is so low that the membranes have already transformed into other phases or that this transition is happening at too low a γ value that is beyond the reach of experiments. The complex behavior of the width and pitch with varying γ is better understood by looking at the edge-to-area ratio as shown in Fig. 3.6. When the line tension is reduced, the edge-to-area ratio increases monotonically whatever value k¯ is. This is because lower line tension favors structures with longer edge length. Also clearly seen in this figure is the shrinking transition of TRs mentioned above for k¯ = 0.10 and 0.15 given by the finite jump in the value of the edge-to-area ratio. 46 5 4.5 4 3.5 Edge-to-Area Ratio 3 2.5 2 1.5 1 0.5 0 0.05 0.1 0.15 0.2 0.25 γ Figure 3.6: Edge-to-area ratio of TRs as a function of the line tension γ. k¯ = 0.05, 0.10 and 0.15 are represented by red squares, green triangles and blue circles respectively. The lines are guide to eyes. 0.85 0.8 0.75 θmax 0.7 0.65 0.6 0.55 0 0.05 0.1 _ 0.15 0.2 0.25 k Figure 3.7: Maximum tilt angle at the edge of the TRs (red circles) on the phase boundary ¯ with comparison to the maximum tilt angle between semi-infinite disks and TRs as a function of k, at the edge of semi-infinite disks (black horizontal line). It is clearly seen that the maximum tilt angle decreases as the phase transition from the semi-infinite disks to TRs happens. 47 4.2 16 4 14 3.8 12 3.6 10 3.4 R 3.2 8 b 3 6 2.8 4 2.6 2 2.4 2.2 0 0 0.05 0.1 _ 0.15 0.2 0.25 k Figure 3.8: R (red squares with y-axis on left) and b (green triangles with y-axis on right) at the phase boundary between semi-infinite disks and TRs as a function of k.¯ Another experimental observation is that the maximum tilt angle at the edge θmax decreases across the phase transition from semi-infinite disks to TRs which our theoretical model agrees with. Shown in Fig. 3.7 is θmax at the phase boundary as a function of k¯ which is always smaller than the maximum tilt angle at the edge of semi-infinite disks whatever the value k¯ is in our theoretical model. Fig. 3.8 shows ¯ We see that generally the the dependence of R and b at the phase boundary on k. values of R, b and θmax at the phase boundary all decrease with k¯ but for R there is an oscillating region within the range from k¯ ∼ 0.05 to k¯ ∼ 0.12. To understand the role of chirality in our theoretical model, we study the phase diagram with stronger chirality q = 1.00 shown in Fig. 3.9. Note that the semi- infinite-disk-to-TR transition denoted by the red symbols still preempts the semi- infinite-finite-disk transition represented by the dashed horizontal line. The phase diagram is similar to the one with q = 0.71, however, the two phase transitions 48 0.6 Semi-Infinite Disks 0.55 γ 0.5 Twisted Ribbons 0.45 0 0.1 0.2 0.3 0.4 0.5 _ k Figure 3.9: Phase diagram with q = 1.00 spanned by k¯ and γ exhibiting semi-infinite flat disks, finite disks and TRs. Red symbols denote the first order phase transition from semi-infinite disks to TRs when γ is reduced and the line is guide to eyes. The dashed line at γ = 0.506 is the transition between semi-infinite and finite disks preempted by the TR transition. are shifted to higher γ values for larger chirality. This is easily understood because stronger chirality favors structures with more twist in the membranes so that the transition occurs at higher values of γ. Actually, we can define an effective line tension combining the effect of the pure edge and the reduction in the free energy by introducing an edge in the director field, i.e. combining the effect of γ and q. The semi-infinite disk model is used to derive the expression for this γef f , and the result is " # ∞  2 K dθ C K Z γef f = γ + −q + sin2 θ − q 2 dx (3.21) 0 2 dx 2 2 where the integral represents the free energy stored in edge with unit length. Plug- ging in the solution for the tilt angle θ in Eq. (2.10) and integrating, we arrive at p γef f = γ + K(1 − 1 − q 2 ) − Kq arcsin q (3.22) 49 Frank Twist 0.15 Tilt Gaussian Edge 0.1 0.05 0 -0.05 f -0.1 -0.15 -0.2 -0.25 0.27 0.275 0.28 γ Figure 3.10: Plot of total Frank elastic energy density excluding the twist term (red squares), the twist term proportional to q (green triangles), tilt free energy density proportional to sin2 θ (blue right triangles), Gaussian curvature free energy density (grey diamonds) and edge energy per unit area (orange circles) in the total free energy density for stable or meta-stable TRs across the phase boundary between flat membranes and TRs (vertical dashed line) as functions of γ at k¯ = 0.15. Plugging in q = 0.71 and 1.00 results in γef f = γ − 0.2647 and γef f = γ − 0.5708 respectively. By converting to γef f , we see that the critical value of γef f , which is 0.012 for q = 0.71 and −0.026 for q = 1.00, does not change too much when q is varied and this value may be constant for such systems. In order to look for the major driving force of the phase transition between flat membranes and TRs, we plot the contributions from different terms to the total free energy density for a stable or meta-stable TR across the phase boundary as a function of γ at k¯ = 0.15 in Fig. 3.10. From this figure, it is observed that TRs become unfavorable with increasing line tension because of the increase in the twist free energy related to chirality which dominates over other terms in the total free energy density, i.e. the chirality drives the phase transition from flat membranes to TRs. It is expected that TRs will never be favored when the membrane is achiral 50 and experimentally twisted ribbons have not been observed in achiral systems. We shall present here a simple-minded model to explain the unusual positive Gaussian curvature modulus as stated previously. There are many theoretical models concerning the Gaussian curvature modulus in membrane systems[83, 84, 85, 86, 87], however most of them predict negative Gaussian curvature modulus. Our starting point is that the major force in fd systems is the entropy-driven depletion force whose effect is to minimize the volume occupied by the viruses. We assume that the viruses are rigid rods whose center-of-masses lie in a smooth plane which is defined as the mid-plane of the monolayer. Given the up-down symmetry of the viruses, it is natural to further assume that the density of the molecules on the mid-plane of the membrane is conserved when the physical parameters, such as the Gaussian curvature modulus and the line tension, remain unchanged. Under these assumptions, the membrane can be viewed as a stack of thin layers each with an infinitesimal thickness. Each layer is defined by the displacement δ from the mid- plane as following when this displacement is small compared to the principal radii of curvatures of the mid-plane[81] Y ′ = Y(u1, u2 ) + Nδ (3.23) where Y defines the mid-plane. The area element of the displaced surface is calcu- lated to be dA′ = dA(1 + 2Hδ + KG δ 2 + O(δ 3 )) (3.24) with dA the area element on the mid-plane. Thus the volume element occupied by the viruses whose centers-of-mass are inside a surface element dA in the mid-plane 51 is given by h/2 1 1 Z dV = dδ(1 + 2Hδ + KG δ 2 )dA = dA(h + Hh2 + KG h3 ) (3.25) −h/2 2 12 where h is the length of the viruses. We see from this equation that a negative Gaussian curvature KG reduces the volume element while the effect of mean curva- ture is always to increase the volume. In other words, it costs energy to introduce a positive Gaussian curvature to the membrane resulting in a positive Gaussian cur- vature modulus. This very simple prediction seems contradictory to the commonly measured negative k¯ in many lipid systems but we should point out that usually in those lipid systems the constituent molecules have different microscopic structures and the behavior may be very different for the fd systems. 3.4 Conclusions In this chapter, the basic theoretical model brought forward in the previous chapter is applied to the geometric shape of TRs with the aid of knowledge from differential geometry. The E-L equation for the tilt angle is derived and solved with a differential equation solver in MATLAB. By fitting the geometrical parameters of the TRs in our model to experimental values, we conclude that the Gaussian curvature modulus k¯ is positive which seems unusual and a brief explanation is given to this based on the depletion interaction. Comparison of the free energy density for flat membranes and TRs results in a phase diagram where the phase transition from semi-infinite disks to TRs is first-order when the line tension is lowered from above, preempting the phase transition to finite disks. According to the results in chapter 2, disks with radius larger than 5 penetration depths can be modeled as semi-infinite membranes. 52 In experiments, the phase transition is indeed observed between such large disks and TRs without appearances of any small disks which our theory confirms. We find that this phase transition is driven by chirality by looking at separate terms in the free energy density. Besides this major point, there are also many other points of agreement between our results and the experimental facts listed as follows: i) the width and pitch of the TRs near the phase boundary predicted by our model are close to experimental values; ii) the maximum tilt angle at the edge of the membranes decreases when disks transform to TRs; and iii) the pitch of TRs decreases with decreasing line tension near the phase boundary. One notable discrepancy, however, is that the critical line tension obtained in our theory is one order of magnitude lower than the experimental value which may be attributed to the over-simplified edge model used. Chapter Four Instability from Disks to Twisted Ribbons 54 Figure 4.1: Several snapshots from experiments during the phase transition process from disks to TRs from left to right. Reproduced from Ref. [33]. 4.1 Experimental Facts In the last chapter, we saw that our theoretical model succeeded in predicting the first-order phase transition from disks to TRs when the line tension is lowered. But the theoretical model did not propose a kind of instability for this transition to take place between two such distinctive shapes. This is the main topic of this chapter. Before proposing theories, we introduce below the experimental discoveries about the process of this phase transition. In the experiments, the phase transition from disks to TRs can be induced by several different ways. Decreasing the temperature which effectively increases the strength of the chirality of the fd viruses, and lowering the concentration of Dextran (effectively the line tension), are the two ways to make stable TRs. It is easier to observe the dynamical process by tuning temperature than by lowering the line tension in experiments while in our theoretical model, the line tension which only appears in an add-on term is easier to tune. So we focus on a model where γ is tuned and we will see that the kind of instability produced by our theory is similar to the process observed by experiments when the temperature is lowered. It is found in experiments that upon lowering the temperature approaching the phase boundary, ripples start to appear on the edge of disks replacing the round shape of the edge. Then these ripples grow in size and protrude radially from the disks with a clearly- seen internal twist structure like that in the TRs. As temperature drops further, 55 Figure 4.2: Stretching of a flat membrane utilizing two optical traps to make a TR. The position of the right hand side of the disk is fixed by one optical trap while the other end is stretched. these protrusions become TRs connected at the center of the disks making a shape resembling spiders and finally they consume all the material in the disks and then break up into separate TRs. Several snapshots of this process are shown in Fig. 4.1[33]. The ripple-edged disks are not stable phases and they will relax back to round disks if the sample is kept at constant temperature for some time. However, we believe that these ripples are the key structures enabling the phase transition between disks and TRs to happen. Thus, a theoretical model based on such kind of ripples attached to the edge of a disk is proposed in section 4.2 and 4.3. Another method to make a TR is to use two optical traps to stretch a membrane. One of the traps holds the position at one edge of the membrane while the other one pulls the other edge. At the initial stage of stretching, the membrane remains flat with a protrusion produced by the force exerted on the edge. With more stretching, twist structure establishes on the protrusion and abruptly the disk transforms to a TR. Several snapshots during this process are shown in Fig. 4.2. This whole process happens in the region of the phase diagram where disks are the stable structure, so it is naturally expected that when the optical traps vanish, the TR should relax back to the disk and this phenomenon is indeed observed in experiments. The force felt by the optical traps during this process is measured in experiments and is shown in Fig. 4.3. It is seen that, at first, the force is almost linearly increasing with the stretching, indicating that no twist has formed. At a critical stretching, the force abruptly drops and this is where the transition happens. We show a brief calculation 56 Figure 4.3: The force felt by the optical trap which is stretching a flat membrane to make a TR. Note that the end of the linear region is the point where the TR forms. Reproduced by permission from M. Zakhary. based on this experimental observation in section 4.4. 4.2 Model for Rippled Disks and the MC Method We propose a simple geometrical shape mimicking the disks with ripples attached to the edge. We assume these ripples can be described mathematically by sinusoidal waves both in the out-of-plane and in-plane directions with the out-of-plane height decaying exponentially into the interior. We define the plane of a flat disk without ripples as the x-y plane and put the center of the disk at the origin. We have applied our theory on two slightly different model shapes of a rippled disk. Shape I is written as[89] hI (r, φ) = Az sin(aφ) exp[−b(R0 − r)] (4.1) 57 and Shape II is written as hII (r, φ) = Az sin(aφ) exp[−b(R0 − (r + Ar cos(aφ + ψ)))] (4.2) with the projection of the edge on the x-y plane defined as R(φ) = R0 + Ar cos(aφ + ψ) (4.3) where h(r, φ) is the out-of-plane height of the membrane at the point given by the polar coordinate (r, φ), Az and Ar are the amplitudes of the ripples in the z-direction and radial direction respectively, R0 is the radius of the disk without ripples attached, a is the number of the ripples along the edge, b is the decay length of the ripples in z-direction into the interior of the disk and ψ is the phase difference between the in-plane ripples and out-of-plane ripples. We should point out that these two model shapes are similar. And we will see that the results on both shapes are also similar, indicating that the qualitative behavior of the rippled disks does not depend on the exact mathematical form chosen for the shape. So we argue that the behavior of the rippled disk is rather a general phenomenon which captures the key feature of the underlying physical process. Fig. 4.4 shows an example of Shape I of the rippled disk and the similar Shape II is not shown. It is noted that there is a twist structure residing in the region of each protrusion which acts as the seed for a TR to grow from. Although it is not clear from experiments whether these rippled disks are stable structures, we assume these model shapes are equilibrium structures and solve for the director field on them. Defining r and φ as the surface coordinates u1 and u2 respectively (see Sec. 3.1 for the definition of u1 and u2 ), we can follow the steps in Chapter 3 as we did for TRs to derive the E-L equations whose solution gives 58 (a) (b) Figure 4.4: Shape I of a rippled disk with Az = 0.5, Ar = 0.5, ψ = 0.0, a = 5, b = 1.0 and R0 = 5.0. (a) Viewed from above. (b) Viewed from a tilted angle. We use large values of Az and Ar for clarity of the image. the director field. However, this is an impossible task because of the complexity of the shape which results in a very complicated form of the E-L equations. Thus we switch our method to obtain the director field from directly solving the differential equations to a Monte Carlo method whose starting point is the free energy itself. It should be pointed out that effectively we are still solving the E-L equations by keeping the MC temperature very low (we will see the value of the temperature in the following part). We discretize the model shapes using a square lattice in the x-y plane with lattice spacing 0.05. We define the height function at each lattice site through Eq. (4.1) by converting Cartesian coordinates to polar coordinates. Note that the shape is prescribed and not allowed to vary during the computation. We do not include the lattice sites which lie outside the periphery given by Eq. (4.3) except those sites which have nearest sites inside the periphery. These serve as ”ghosts” and allow us to impose a BC on the director field. A unit vector is put on each lattice site including the ghost sites representing the director field. These directors interact via the free energy density given in Eq. (3.10). To calculate the geometrical quantities appearing in Eq. (3.10), x and y are defined as u1 and u2 and a finite difference 59 algorithm is used. The derivatives of the director field are also calculated through a finite difference algorithm. Then the free energy can be calculated at each lattice site and the total free energy is obtained by taking the numerical integration over the surface. The free energy density is calculated by adding the total free energy to the edge energy and then dividing it by the total area. In each MC cycle, the director at each lattice site is allowed to move once with the magnitude of the test move selected so that the total acceptance ratio is around 50%, and then the resulting free energy difference ∆F is calculated after each move. The rejection or acceptance of one move is then determined by the usual Metropolis algorithm described in Chapter 1 where the temperature is kept at a low value such as 0.0001. We start the configuration where all the directors are pointing in the z direction and run 50000 MC cycles to let the director field relax and then an additional 50000 MC cycles to collect data. Convergence is carefully checked. We have also checked that using random initial conditions for the director field produces identical results. Unlike the usual method in a MC simulation where the free energy at each cycle is recorded and averaged, we collect the orientations of directors at each cycle and obtain the averaged orientations for each director, and calculate the free energy for this ”averaged configuration”. By doing this, we effectively reduce the effects of thermal noise and the configuration obtained is close to the zero-temperature solution, i.e., the solution to the E-L equations. The validity of this method was checked by comparing the results of this method for a round flat disk to the results presented in Chapter 2 where the EL equation was solved directly. The boundary condition for the director field at the edge is either free or fixed at a π/2 tilt following the tangential direction of the edge for Shape I and only free BC is imposed on Shape II for the reasons to be discussed in the following section. There is a slight modification to the shape defined in Eq. (4.1) and (4.2) used 60 in the computation because these two shapes are problematic in that the Gaussian curvature diverges at the center of the disks. To fix this problem, we take the region inside a certain radius to be perfectly flat with h = 0. This limiting radius is determined as follows. We calculate the Gaussian curvature KG (r, ∆r) for a ring of inner radius r and outer radius r + ∆r with ∆r a small quantity compared to r. Physically, KG (r, ∆r) is expected to monotonically decrease going into the interior of the disk. The limiting radius is thus taken at the value where this quantity starts to increase instead of decrease. Although the introduction of this limiting radius produces a discontinuity in the shape, the discontinuity is tiny and does not affect the physical results. By varying ψ in Eq. (4.1), we see that the minimum of the free energy always appears at ψ = 0 for right-handed membranes and at ψ = π for left-handed mem- branes. Thus we shall neglect the effect of this term in the following and only study the right-handed case without loss of generality. In order to study the instability from disks to TRs, we study the variation of the amplitude of the ripples with the line tension γ reduced from the disk region in the phase diagram presented in Chapter 3. So we mainly focus on the two parameters Az and Ar and fix b at 1.0. Physically, the length b is determined by the interplay between the deformation in the director field and the bending of the membrane so that one penetration depth is a reasonable choice. The number of ripples a can also be varied in our computation but we believe that after the ripples first appear at the edge, this number should not change any more. 61 (a) (b) (c) (d) Figure 4.5: Solutions of director field on rippled disks with Az = 0.5, Ar = 0.5, a = 5 and q = 0.71, other parameters are as stated in the main text. Also presented are the corresponding birefringence images. (a)(b)(e) Free boundary condition. (c)(d)(f) Fixed boundary condition where the directors at the edge lie down in the tangential direction of the edge. 62 Table 4.1: Dependence of the out-of-plane height Az and in-plane bulge size Ar on the line tension of rippled disks with free BC. Left columns are results from Shape I while the right columns are for Shape II. Note the Ar values in the last line are the largest value computed respectively. γ Az Ar γ Az Ar 0.325− 0.0 0.0 0.230- 0.0 0.0 0.295 − 0.325 0.10 0.2 0.198-0.230 0.02 0.1 0.290 − 0.295 0.10 0.3 0.195-0.198 0.05 0.5 0.280 − 0.290 0.10 0.5 0.194-0.195 0.05 0.7 0.245 − 0.280 0.12 0.5 0.191-0.194 0.05 0.8 0.200 − 0.245 0.10 ≥ 0.7 0.180-0.191 0.05 ≥ 1.1 4.3 Results on Rippled Disks Fig. 4.5 (a)-(d) show solutions of the director field on Shape I rippled disks with Az = 0.5, Ar = 0.5, a = 5 and q = 0.71, other parameters are as stated in the previous section. Note that the directors are mostly perpendicular to the x-y plane at the center of the disks as expected. The director field does show some differences when the boundary condition is switched from free BC shown in (a) and (b) to fixed BC shown in (c) and (d) where the directors at the edge lie down and follow the tangential direction of the edge. However the resulting birefringence images are similar as shown in (e) for free BC and (f) for fixed BC. The birefringence images both have darker protrusions and brighter regions between the adjacent protrusions. For the fixed BC case, the brighter region is larger and the contrast is stronger because the directors tilt more than in the free BC case. Experimental observations tend to favor the fixed BC but for the completeness of our theory and consistency with the theory for disks and TRs, both results are discussed below. The corresponding graphs for Shape II with free BC are similar and not shown here. As before, we look for the optimal Az and Ar values minimizing the free energy per unit area for fixed k¯ = 0.15 and γ. The results are shown in Table 4.1 with free BC for Shape I on the left with three ripples (a = 3) and Shape II on the right 63 with five ripples (a = 5). We have tested that these specific choices of the value of a do not influence the qualitative results shown above. Due to the limitation of our computing resources, we could not explore the landscape of the free energy as a function of Az and Ar in full detail but we believe the results shown here are very representative. The Ar values in the last line are the largest test values computed respectively. From Table 4.1, it is seen that when the line tension is high above γc = 0.276, which is the value at the phase boundary between disks and TRs, both vertical height and radial size of the ripples vanish and the disks retain perfect flat round shapes. Close to the critical line tension, the height and radial size both grow. For Az , the rate of grow as a function of γ is slow and it maintains a small value compared to the penetration depth for all γ values computed, which means the disks do not deviate from their planes very much. Meanwhile, the radial size Ar grows rapidly inside a narrow range of γ to a large value. Specifically, Ar changes from 0.2 to a value no less than 0.7 in the range between γ ∼ 0.30 and γ ∼ 0.24 for Shape I, and even more dramatically for Shape II, from 0.1 to not smaller than 1.1 in the range between γ ∼ 0.20 and γ ∼ 0.19. This phenomenon is in accord with the experimental observations mentioned in Sec. 4.1, where the growth of ripples in the radial direction is much more obviously seen than that in the z-direction before those protrusions transform to TRs. The fast growth of Ar inside a narrow range of γ is the indication of an instability from disks to other shapes like TRs. Given the fact that these ripples also have negative Gaussian curvature, they are very likely to transform to TRs as minimal surfaces (zero mean curvature). Thus, we identify this range of γ as the most probable region for the phase transition to happen. Although there are slight numerical differences in the value for the two shapes where it is around 0.27 for Shape I and around 0.20 for Shape II, they are both close to the critical line tension γc = 0.276 found in the previous chapter and the general behaviors of these two shapes are very similar. These results are in accord with the observation that 64 these ripples act as seeds for the formation of TRs. As both shapes yield similar results, we believe that they are not our numerical artifacts of choosing a particular shape of the ripples. However, for fixed BC where the directors are lying down at the edge following the tangential direction, such instability is not seen for Shape I and to save computing time we did not impose this BC on Shape II believing their behaviors should also be similar in this case. Given the fact that our theoretical model for disks and TRs are constructed with free BC in Chapter 2 and 3, we argue that the free BC is still a reasonable choice here. 4.4 Instability from a Disk to a TR upon Stretch- ing In this section we model the instability toward a TR induced by stretching a flat membrane utilizing optical traps as mentioned in Sec. 4.1. Note that this instability happens totally inside the region in the phase diagram where disks are the equilibrium shape. We model this problem as a perturbation problem where the ribbon-like structure is superposed to a flat semi-infinite membrane in the edge region. And for simplicity, only the Helfrich free energy is considered in this section. The idea is to decide whether the stretched membrane should remain flat by looking at the free energy perturbation due to the ribbon-like structure. As before, the unperturbed membrane is defined by Y(u1, u2 ). The perturbation is taken locally in the direction of the surface normal resulting in a perturbed surface 65 given by the following form[75] Y ′ = Y + Ψ(u1 , u2 )N (4.4) where Ψ(u1 , u2) is the amplitude of the perturbation which is taken to be small. From this definition, we can then calculate all the perturbed geometrical quantities defined in Eq. (3.2). We refer the readers to Ref. [75] for detailed formulae which will not be listed here. The resulting variation in the Helfrich free energy is ¯ G )√gdu1 du2 + ¯ G )δ √gdu1 du2 Z Z 2 δFH = (2kδ(H ) + kδK (2kH 2 + kK (4.5) where δ denotes the perturbation to the quantity that follows. The initial unper- turbed shape is a flat semi-infinite membrane where H = 0 and KG = 0 so that the second term drops out and the formula is reduced to ¯ G )√gdu1 du2 Z δFH = (2kδ(H 2 ) + kδK (4.6) Now we construct the coordinate system on the surface. The membrane plane is chosen to be the x-z plane where the semi-infinite surface occupies the z < 0 half- plane. In this specific case, the variations of the geometrical quantities in Eq. (3.2) are simplified to δgij = Ψi Ψj , δg = Ψ21 + Ψ22 , δLij = Ψij , δL = Ψ11 Ψ22 − Ψ212 , (4.7) 1 δH = (Ψ11 + Ψ22 ), δKG = Ψ11 Ψ22 − Ψ212 2 where the subscripts on Ψ denote the derivative taken with respect to the corre- sponding coordinate. Thus the perturbation in the Helfrich free energy is given 66 by k Z Z δFH = (Ψ11 + Ψ22 ) dxdz + k¯ 2 ¯ (Ψ11 Ψ22 − Ψ212 )dxdz ≡ 2kA + kB (4.8) 2 Now let us write down our exact form of the ribbon-like perturbation on the surface. The unperturbed surface is a semi-infinite flat membrane in the x-z plane with a certain shape of a protrusion attached to the edge at z = 0, with the magnitude of the protrusion defined as z0 measuring the extent of the stretching. Recall that a TR is parameterized as (ρ cos φ, ρ sin φ, bφ). Now in the coordinate system, x = ρ cos φ and z = bφ, so y = x tan(z/b), where 2πb is the pitch length. Because the perturbation only happens around (0, z0 ) where the optical trap functions, an exponential decay of the shape is assumed both in x-direction and into the interior of the flat membrane resulting in "  # 2 ze(z−z0 )/λ1   x Ψ(x, z) = x tan exp − (4.9) b λ2 where λ1 and λ2 are two parameters defining the rate of the decay in the z-direction and in the x-direction respectively. b is kept large to ensure that the ribbon-like perturbation is small. As for the protrusions, we used two different shapes which are shown in Fig. 4.6. The first is a Gaussian bump whose edge is given by z = z0 exp[−(x/λ3 z0 )2 ] where λ3 is a proportionality constant measuring the ratio between the half width of the Gaussian bump and z0 . Fig. 4.7 shows an example of the perturbation on this shape. The perturbed structure is obviously a twisted structure like the TRs. The second shape for the protrusion is an equilateral triangle with height z0 . Plugging Eq. (4.9) into Eq. (4.8) gives the resulting perturbation in the Helfrich 67 Figure 4.6: Theoretical model for the shape of the stretched membrane. The colored region is the membrane. Note that the ribbon-like perturbation is not shown here so that the membranes remain flat. The Gaussian bump model is shown on the left while the equilateral triangle model is on the right. Figure 4.7: An example of a ribbon-like shape perturbation on a flat semi-infinite disk with a Gaussian-bumped edge. The parameters are λ1 = λ2 = λ3 = 1.0, b = 104 and z0 = 1.778. Note that the scale on z-direction is small ensuring the ribbon-like structure is a perturbation. 68 free energy. There is an extra variation in the free energy due to the variation in the length of the edge which is proportional to the line tension γ. The total variation in the free energy is given by the sum of these two contributions. We use MATHEMATICA to calculate the corresponding numerical integrals. It is natural to assume that the two characteristic lengths λ1 , λ2 are on the order of penetration depth. This is because if they are too small compared to the penetration depth, then the surface is curved so much that the director field on the membrane cannot respond to the membrane shape and will result in a high energy for the ribbon-like perturbation. They also cannot be large because the perturbation structure should only exist in the region around the point being stretched. So both λ1 and λ2 are kept at 1.0 through out the calculation and we have verified that a slight change of their values does not influence our qualitative results. As for λ3 , it is also natural to expect that the half width of the Gaussian bump induced by the stretching is proportional to the extent of stretching z0 , and this proportionality constant is kept at 1.0 and we have also verified that a slight change of the value does not influence our results. And as mentioned above, without loss of generality, b is kept at 103 which is a large value compared to the other length scales in the problem. We have verified that identical results are produced with b varying from 102 to 104 . Let us present the result for the Gaussian bump protrusion at first. (δH)2dxdz R R Shown in Fig. 4.8 are the dependence of B ≡ (δKG )dxdz, A ≡ and the variation of the edge length δL on the extent of stretching z0 labeled as red squares, green triangles and blue circles respectively. As the stretching becomes larger, both the length of the edge and the Gaussian curvature increase monotonically in absolute value. But the mean curvature term assumes a parabolic shape with a valley at around 1 penetration depth. 69 1.5E-05 1E-05 5E-06 0 A, B, δL -5E-06 -1E-05 -1.5E-05 -2E-05 -2.5E-05 -3E-05 1 2 3 4 5 z0 Figure 4.8: The variation of the integrated Gaussian curvature (B in Eq. (4.8), shown as red squares), the integrated square of mean curvature (A in Eq. (4.8), shown as green triangles) and the variation of the edge length δL (shown as blue circles) as functions of the extent of stretching z0 for a protrusion with the shape of a Gaussian bump. 6E-07 (a) 4E-06 (b) 3.5E-06 4E-07 γ=0.3,k=0.00 γ=0.3,k=0.10 γ=0.3,k=0.05 3E-06 γ=0.4,k=0.10 γ=0.3,k=0.10 γ=0.5,k=0.10 2E-07 2.5E-06 γ=0.6,k=0.10 2E-06 δF δF 0 1.5E-06 -2E-07 1E-06 5E-07 -4E-07 0 -6E-07 -5E-07 1 2 3 4 5 1 2 3 4 5 z0 z0 Figure 4.9: The total variation in the free energy δF for the ribbon-like perturbation as a function of the extent of stretching z0 with different choices of the mean curvature modulus and the line tension with the Gaussian bump protrusion. (a) γ = 0.3, k = 0.0, 0.05 and 0.10 from bottom to top. (b) k = 0.10, γ = 0.3, 0.4, 0.5 and 0.6 from bottom to top. The dashed horizontal line marks the line at δF = 0 below which the perturbation is favored. 70 To be consistent with our previous results on TRs, k¯ is kept at 0.15. The total variation of the free energy is plotted in Fig. 4.9 as a function of the extent of stretching z0 for different choices of the mean curvature modulus k and the line tension γ. In Fig. 4.9(a), γ = 0.3, which is in the region of the disks in the phase diagram computed in Chapter 3, and k = 0.0, 0.05 and 0.10 from bottom to top. In Fig. 4.9(b), k is kept at 0.10 and γ = 0.3, 0.4, 0.5 and 0.6 from bottom to top. The dashed horizontal lines in both figures mark the δF = 0 lines below which the perturbation is favored. It is seen that generally, δF assumes the shape of a parabola, and the perturbation is not favored with δF > 0 at the beginning of the stretching which indicates that the membrane should remain flat at the initial stage of stretching. When the critical length of stretching is reached, the perturbation starts to be favored and this small perturbation may act as a seed where the whole membrane transforms to a TR. This general behavior agrees with the experimental observations mentioned in Sec. 4.1. The critical value of the stretching is z0c = 0.569 ± 0.001 for k = 0.1 and γ = 0.3 . We should point out that this critical value is one order of magnitude lower than the experimental observation where the TR starts to form after the stretching reaches several penetration depths. We attribute this difference to the absence of the director field in our simple model. Adding the director field back should increase the theoretical critical stretching. This is because the stretching process happens totally inside the disk region where the twist structure in the director field in a TR is not favored compared to that in a flat membrane; thus it should shift the transition to larger values of stretching. It is observed that z0c increases with k and γ. We should point out that after δF first changes sign from positive to negative upon stretching, the phase transition to TR happens so the points to the right in the figures where δF changes sign again from negative to positive might be irrelevant to the problem considered here. 71 This general behavior breaks down when k is close to 0.0 where the perturbation is favored at the very beginning of the stretching. Thus we suggest that, when only geometrical factors are considered as in the Helfrich free energy and edge energy, to understand the effect of the mean curvature modulus in the following way. Although the disks and TRs are both minimal surfaces with zero mean curvature where k seems to be irrelevant, the effect of k is to place a free energy barrier between these two phases preventing the phase transition from happening. This conclusion is supported by the data shown in Fig. 4.8 where the perturbation in the mean curvature term remains positive when the stretching is small while the other terms are approaching zero, so that the ribbon-like perturbation is not favored at the beginning of the stretching purely from the mean curvature perspective. Also, the increase of γ to a large value would not favor the ribbon-like perturbation either. This is easy to understand because a TR is a structure with a longer edge which is not favored by a large γ. Now we present the result for the equilateral triangular protrusion. Shown in Fig. 4.10 are the dependence of B, A and δL on the extent of stretching z0 labeled as red square, green triangles and blue circles respectively, similar to Fig. 4.8 for the Gaussian bump protrusion. The trend of change of the terms are identical to those in the Gaussian bump model. It is more obvious that the mean curvature term remains positive when the stretching is small for the triangular protrusion, which supports our assumption that the effect of mean curvature is to place an energy barrier between flat membranes and TRs. Similar to Fig. 4.9, the total δF is plotted with different choices of k and γ in Fig. 4.11 for the equilateral triangle protrusion. The major difference between these two models is that the free energy variation does not become positive again for the triangle model in the range of z0 calculated. But as pointed out in the 72 2E-06 1E-06 0 -1E-06 -2E-06 A,B,δL -3E-06 -4E-06 -5E-06 -6E-06 -7E-06 -8E-06 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 z0 Figure 4.10: The variation of the integrated Gaussian curvature (B in Eq. (4.8), shown as red squares), the integrated square of mean curvature (A in Eq. (4.8), shown as green triangles) and the variation of the edge length δL (shown as blue circles) as functions of the extent of stretching z0 for a protrusion with shape of a equilateral triangle. (a) (b) 2E-07 2E-07 0 0 -2E-07 δF δF -2E-07 -4E-07 γ=0.3,k=0.00 γ=0.3,k=0.10 γ=0.3,k=0.05 γ=0.4,k=0.10 γ=0.3,k=0.10 -4E-07 γ=0.5,k=0.10 γ=0.6,k=0.10 -6E-07 -6E-07 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5 z0 z0 Figure 4.11: The total variation in the free energy δF for the ribbon-like perturbation as a function of the extent of stretching z0 with different choices of the mean curvature modulus and the line tension with the equilateral triangle protrusion. (a) γ = 0.3, k = 0.0, 0.05 and 0.10 from bottom to top. (b) k = 0.10, γ = 0.3, 0.4, 0.5 and 0.6 from bottom to top. The dashed horizontal line marks the line at δF = 0 below which the perturbation is favored. 73 previous context, this might be irrelevant to the problem considered here. As for the critical point where the perturbation starts to be favored, i.e. the point where δF first becomes negative, this shape predicts that z0c = 0.724 ± 0.001 for k = 0.1 and γ = 0.3 which is close to the value given by the Gaussian bump model and they are consistent with each other. Comparing the results given by these two models, we conclude that theoretical- ly, upon stretching, a flat membrane will transform to a TR when the stretching reaches a critical amount on the order of about half a penetration depth, which is qualitatively in accord with the experimental observations. Although quantitative- ly, our prediction for the critical stretching is smaller than the experimental values, we should point out again that only Helfrich free energy is considered in our model and that the addition of the director field should increase this critical value and improve the agreement with experiments. By choosing two different shapes for the prolonged membrane, we believe that our model captures the principles governing this instability and does not depend on the exact mathematical form of the shape chosen. 4.5 Summary In this chapter, two kinds of instability from disks toward TRs are discussed. The first is related to the phase transition spontaneously happening when the line tension is lowered starting from the disks region in the phase diagram. We study a theoretical model where ripples both in the out-of-plane and in-plane directions are attached to the edge of a flat disk acting as seeds for the TRs to grow from. It is discovered that the size of the ripples remains very small deep inside the disks region while it 74 abruptly becomes large in a narrow window of γ when γ is lowered. The second instability is related to the stretching of the membrane using optical traps where the ribbon-like structure is imposed as a perturbation to the flat membrane. The free energy perturbation is calculated by considering the Helfrich free energy and edge energy. It is found that the perturbation is unfavorable at the beginning of the stretching and only after the stretching reaches a critical extent then the perturbation becomes favored and would act as the seed to become a TR. Both theoretical results match the experimental observation qualitatively and improve our understanding of the process of the phase transition from disks to TRs. Chapter Five Smectic-A Disks Formed by Mixture of Molecules with Opposite Handedness 76 Figure 5.1: Experimental images of the mixture membranes. Top left: a membrane with ’cook- ie’ patterns on the edge. Top center: birefringence image of a mixture membrane. Top right: 3D-PolScope image of one cusp on the mixture membrane where the change of the directions of the directors across the cusp are observed. Bottom: schematic of the cusps formed by the rods. Reproduced by permission from Dogic and coworkers at Brandeis University. 5.1 Experimental Facts The systems discussed in previous chapters are composed of wild type fd -wt viruses. A new kind of structure was recently discovered in experiments when the right- handed fd -Y21M viruses and the left-handed fd -wt viruses were mixed together in an aqueous solution with proportionalities selected so that the whole system is close to the achiral limit. The mixture forms monolayer disks with ’cookie’ shapes which are flat disks with out-of-plane cusps located at the edge and in-plane radial bulges between two adjacent cusps as shown in Fig. 5.1. The directors following the edge change directions going from one side of a cusp to the other which can be observed in 3D-LC-PolScope images, an example of which is shown in the top right image of Fig. 5.1, indicating these cusps are defects in the director field. The two kinds of viruses with opposite handedness can be labeled by fluorescence with different colors. It is observed from the resulting fluorescence images shown in Fig. 5.2 that 77 Figure 5.2: Fluorescence image of mixture membranes where molecules with opposite handedness are labeled with red and green respectively. It is seen that the mixture is uniform so that it is close to the achiral limit. Reproduced from permission by Dogic and coworkers. they mix uniformly in the disks without any phase separation which makes the entire membranes effectively achiral except at the edge. The manipulation of the positions of the cusps at the edge by optical traps shows that these cusps are structures at least metastable and the spacings between nearest cusps tend to stay around equilibrium values. The force measured when the cusps are dragged from their equilibrium positions is on the order of 10−2 pN. Recently, there is evidence in experiments that ’cookie’ disks can be equilibrium structures for pure achiral membranes without the need of a mixture[88]. In the following part of this chapter, we follow our previous work to write down a theory for these achiral ’cookie disks’. A Monte Carlo method is used to numerically solve the model. And we determine the physical conditions for these structures to be stable or at least metastable. 78 Z Y X Figure 5.3: An example of the model shape used for the mixture membrane with radius R = 5.0, Az = 1.0, b = 0.5 (see Eq. (5.1)), and Ar = 0.5 (see Eq. (5.3)) in the unit of penetration depth and the 4 defects are distributed along the edge with equal spacing between adjacent cusps. 5.2 Theoretical Model and Numerical Method We now consider the free energy F of a model of the ’cookie disks’ observed in experiments. We still use the general form of the free energy derived in Sec. 3.1. A ’cookie disk’ is an essentially flat, nearly circular Sm-A monolayer lying in the x–y plane with an array of cusps at the monolayer edge rising out of the plane and in–plane radial bulges between adjacent cusps. Alternating regions of left– and right– handed director twist meet at the cusps and the presence of the cusps reduces the director mismatch between the alternating regions of twist. Given the alternating regions of twist and the tendency of the cusps to reduce the Frank energy, neighboring cusps must have heights h(x, y) (in the z direction) of opposite signs. Likewise the alternation of the twist at the edge requires that cusps must appear in pairs. Fig. 5.3 shows an example of the shape. 79 We model the cusps as local peaks or valleys with an exponential decay of the height |h| from the peak of the cusp[89]. Specifically, for a single cusp centered at ri = (xi , yi ), the height is given by   i |r − ri | hi (r − ri ) = (−1) Az exp − (5.1) b where Az is the maximum height at the peak (Note the factor (−1)i ensuring cusps lying alternatively above and below the plane of the disk), b is a characteristic length governing the decay of the cusps into the flat central portion of the disk, and r is the position vector in the x–y plane. For a disk with n cusps (n even) the height at any point r is given by the superposition: n X h(r) = hi (r − ri ) (5.2) i=1 We choose the origin of r at the center of the disk and the cusps are located at ri = R(cos φi , sin φi ) where R is the radius of the disk (excluding the bulges) and φi is the polar angle location of cusp i. We model the in-plane radial bulges between neighboring cusps as follows:   π(φ − φi ) Ri (φ) = R + Ar sin (5.3) φi+1 − φi where Ri specifies the radial coordinate of the disk edge and Ar is the magnitude of the protrusion. The edge of the disk is then given as a piecewise function of these protrusions between neighboring cusps. With an analytic form of the shape specified, we follow our method used in Sec. 4.2 to compute the free energy of the membrane. we discretize the underlying x–y plane using a square lattice of grid size 0.05 and minimize the free energy by means 80 of a Monte Carlo (MC) simulation at low temperature (10−4 in dimensionless energy units with kB = 1) varying the geometrical parameters of the shape. We also carried out a multi-grid computation where the grid is finer at the edge where the directors are twisted and similar results are obtained. Experimentally [33] the director field at the edge of the disk is observed to be tangent to the edge. We impose this boundary condition in our model by introducing ghost directors along the edge of the cookie disks. These directors are fixed tangent to the edge and interact with neighboring directors in the interior of the disk via the Frank free energy. We use a central difference algorithm to compute the derivatives of the director field. The total free energy is computed using numerical integration over all the lattice sites inside the cookie disk. Defects in the director field can Kdef θz2 where θz is the P appear at the cusps and we include a defect core energy angle made by the director with respect to the z-axis and the summation runs over all the lattice sites within one lattice constant of the cusps. The energy Kdef is kept equal to the Frank constant. Experimental observation also suggests that cookie disks are observed to maintain constant area after their formation and thus as we vary the geometric parameters of the shape, we adjust the radius R so that the total area of the disk is kept constant. For each set of geometric parameters we carried out 100000 MC cycles where each director is allowed to move once during each cycle and data are collected over the last 50000 cycles. Equilibrium is carefully checked by verifying that the free energy changes little (< 0.1%) after another 50000 cycles. To reduce the numerical errors induced by the discreteness of the lattice around the edge and cusps, the disks are rotated several times so that their relative position to the lattice is varied and the data are averaged over these rotated configurations. As mentioned in Sec.4.2, since we are effectively solving the EL equations whose solution corresponds to the zero-temperature ground states of the system,we collect 81 Figure 5.4: An example showing the phase separation of the viruses of the theoretical model when the chirality is allowed to vary locally. The coarse grained magnitude of chirality at each lattice site is color coded. the directions of each director in each cycle and compute the averaged directions for each of them and then use this averaged configuration to compute the free energy instead of the usual method that collects the free energy data and takes the average on that. We have checked the validity of this method by considering a flat round disk and obtaining good agreement with the results of our MC simulations and the results obtained analytically in Ref. [71]. 5.3 Results on Mixture Membranes We consider a cookie disk with a fixed number of cusps, specifically, four on a disk with total area approximately 100 for computational simplicity. The decay length b is chosen to be 0.2, of the order of magnitude observed in experiments. If the 82 (a) Z Y X Figure 5.5: (a)An example of the solution to the director field on a ’cookie’ disk with R = 5.41, Az = 0.6, Ar = 0.4, b = 0.2 in unit of penetration depth and the 4 defects are distributed uniformly along the edge.(b)The corresponding birefringence image. local chirality is allowed to vary which physically corresponds to a local assembly of viruses with unique handedness, an obvious phase separation in the protrusion is observed from our computation as shown in Fig. 5.4. But according to the experimental observations, no phase separation is seen. Given this fact, we mainly focus on the achiral limit where q is kept at 0.0 at each lattice site and not allowed to vary. An example of the solution to the director field on a ’cookie disk’ structure is shown in Fig. 5.5(a). Near the center of the disk the directors are approximately perpendicular to the x-y plane as expected because of the tilt energy term. However when approaching the edge, the directors start to tilt over and at last lie down in the x-y plane at the edge. This is because of the boundary condition that is inconsistent with the all-perpendicular configuration which has to be relaxed over a length scale of the penetration depth. It is seen from the figure that the director field establishes defects at the cusps where the relative tilt direction switches going from one side of each cusp to the other side of the same cusp. Fig. 5.5(b) shows the theoretical birefringence image where the gray scale is proportional to sin2 θ. The bulges are the bright regions where the directors tilt more while the cusps are darker because the directors approximately point to the z-direction there. This is qualitatively in agreement with experimental observations. 83 Spacing 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 1.06 1.04 0.92 0.93 0.918 0.925 0.916 1.02 0.92 0.914 0.912 0.915 f f 1 0.91 0.91 0.908 0.905 F/Area 0.98 0.906 0.904 0.9 0 0.2 0.4 0.6 0.8 1 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 Ar Spacing 0.96 0.94 0.92 0.9 0.88 0 0.2 0.4 0.6 0.8 1 Az,Ar Figure 5.6: Free energy per unit area of a ’cookie disk’ as a function of height of the cusps Az (red squares) with the other parameters at Ar = 0.7, b = 0.2, k¯ = 1.5, γ = 1.6 and the cusps are put uniformly around the edge of the disk, as a function of the size of the radial bulges Ar (green triangles) which is also shown in the left inset with the other parameters at Az = 0.6, b = 0.2, k¯ = 1.5, γ = 1.6 and the cusps are put uniformly around the edge of the disk, and as a function of the spacing between adjacent cusps (blue circles) which is also shown in the right inset with the other parameters at Az = 0.6, Ar = 0.7, b = 0.2, k¯ = 1.5 and γ = 1.6. The meaning of the numbers on the spacing axis is as follows, for example, 0.4 means that the four cusps are placed at φ = 0.0, 0.4π, 0.8π, 1.2π and vice versa. Note that for clarity, the scale on the y-axis is magnified in the insets. The lines are sketched as guidance to the eye. 84 In Fig. 5.6, the variation of the total free energy per unit area with the z-direction height of the cusps Az is shown in red squares with the other geometric parameter kept at Ar = 0.7, b = 0.2 and the cusps are distributed uniformly along the edge. The other physical parameters are k¯ = 1.5 and γ = 1.6. The free energy per unit area assumes a parabolic shape with Az in the parameter range explored. Contributions from different terms in the free energy are shown in Fig. 5.7(a). Note that the terms are re-centered around 0.0 in this figure (and also in (b) and (c) in the same figure) for a clear comparison of their relative magnitude and trend of variation. It is seen that the Frank free energy also assumes the parabolic shape as the total free energy. In this parameter range, small Az values results in a high Frank free energy because compared to flat membranes, these out-of-plane cusp structures relaxes the abrupt contact of the directors in the right-handed and left-handed region to different sides of the cusps into a more smooth structure, so the Frank free energy is reduced by the appearance of these cusps. The large value of Az is not favored from the director field perspective because it induces a too large deformation in the region near the cusps which overcomes the gain in the relaxing of the abrupt contact. It is also seen that both the absolute value of Gaussian curvature (note that Gaussian curvature is always negative in this structure) and the length of the edge are monotonically increasing with Az . In the parameter range selected, at the small Az limit, the smallness of the Gaussian curvature and the high Frank free energy dominate and at the large Az side, the longer edge and larger Frank free energy combine to overcome the gain in Gaussian curvature free energy so that both limits are unfavored. With such proper choices of k¯ and γ, a free energy valley emerges at around Az = 0.6. The behavior of the free energy when the in-plane size of the bulges is varied instead of the height of the cusps is shown in Fig. 5.6 with green triangles and with the other geometrical parameters at Az = 0.6, b = 0.2 and the cusps are distributed 85 0.15 (a) 0.015 (b) Frank Frank 0.01 Gaussian Gaussian Edge 0.1 Edge 0.005 F/Area 0.05 0 F/Area -0.005 0 -0.01 -0.05 -0.015 -0.1 -0.02 0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 Az Ar 3 (c) 2.5 Frank 2 Gaussian Edge 1.5 1 F/Area 0.5 0 -0.5 -1 -1.5 -2 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0.5 Spacing Figure 5.7: Contribution from Frank free energy (red squares), Gaussian curvature (green trian- gles) and edge energy (blue circles) in the free energy per unit area for the mixture membranes as functions of (a) the height of the cusps Az (b) the size of the radial bulges Ar and (c) the spacing between adjacent cusps. Note that for better comparison, all the terms are re-centered around 0.0. 86 uniformly. k¯ and γ are as indicated above. Compared to the out-of-plane cusps, the dependence of the free energy on the size of the radial bulges is much weaker, possibly because the cusps are sharp defected structures while the radial bulges are much smoother so that the range of the variation in the free energy density is much smaller than that as a function of Az . Thus, for clarity, it is also plotted in the left inset with a magnified scale. Again, the contributions from different terms are shown in Fig. 5.7(b). Bigger bulges introduce larger deformations in the director field far from the cusps because the directors are forced to follow a more curved edge at the boundary, but bigger bulges can also reduce the mismatch between the directors near the cusps by forcing the directors on both sides of the cusps to become more parallel. These two effects cancel each other. Thus, the Frank free energy seems to be only very weakly dependent on Ar and remains nearly unchanged aside from some numerical fluctuations. As for the edge and Gaussian curvature, the absolute value of the Gaussian curvature increases with Ar however the length of the edge first decreases with Ar and after Ar ∼ 0.15, it starts to increase with Ar . The low value of Ar is not favored due to the dominant small Gaussian curvature while large values of Ar is not favorable because of the rapid increase in the edge length. Thus, a free energy valley appears at a finite Ar which is at around 0.7 with the specific parameter choice here. It is also interesting to study how the free energy changes when the cusps are distributed equidistantly on a certain segment of the edge instead of being placed uniformly along the entire edge of the membrane. This is shown in Fig. 5.6 with blue circles and the other geometrical parameters at Az = 0.6, Ar = 0.7 and b = 0.2. k¯ and γ are as indicated above. The meaning of the number on the ”spacing” axis is as follows. Take 0.4 as an example, it means that the four cusps are placed at φ = 0.0, 0.4π, 0.8π, 1.2π respectively. Thus 0.5 corresponds to the configuration 87 where the four cusps are placed uniformly along the edge. As in the case for Ar , the dependence of the free energy on the spacing is much weaker than that on Az . So again, for clarity, it is plotted in the right inset. We show the contribution from different terms in Fig. 5.7(c). While squeezing the cusps into a narrower segment of the edge creates a more flat segment where no cusps are present which reduces the Frank free energy, it also creates a more curved geometry in the region where the cusps are squeezed, which increases the Frank free energy. From the general trend of the Frank free energy shown in the figure, we can tell that the latter effect overcomes the former so that the Frank free energy increases as the cusps are brought together. The small wiggles in the data are out of the scope of the current discussion. As expected, the length of the edge and the Gaussian curvature both increase in absolute value (note that Gaussian curvature is always negative.) upon squeezing because the membrane is more curved in the squeezed region and this region is the main contribution to these two terms compared to the flatter region. For the k¯ and γ values selected, the trend of the change in the free energy is dominated by the Frank free energy and is at minimum when the four cusps are placed uniformly as shown in Fig. 5.6. Currently, we are unable to explain the appearance of another minimum at spacing 0.25 which may be due to the numerical wiggles in the Frank free energy term. According to experimental observations, it is expected that ’cookie disks’ should be at least metastable with the size of both the out-of-plane cusps and in-plane radial bulges comparable to but smaller than the penetration depth. Recent experi- ments have indicated that these structures can become the equilibrium shape under certain physical conditions[88]. So now the most important question addressed is whether our theoretical model can predict a region in the parameter space where these ’cookie disks’ with physically reasonable geometrical parameters are stable or 88 Figure 5.8: The colored regions in the phase space parameterized by k¯ (X-axis) and γ (Y-axis) in these figures are the zones where a ’cookie disk’ is metastable (blue) and stable (yellow). (a)k¯ = 0.0, (b)k¯ = 0.5, (c)k¯ = 1.0, (d)k¯ = 2.0. Note that the complete area where ’cookie disks’ are stable or metastable should be larger but by the limitation of our computation resources we are unable to explore it in full details. at least metastable. The result is shown in Fig. 5.8. Note that we have introduced the second order Gaussian curvature modulus k¯ because the ’cookie disks’ are constructed to have negative Gaussian curvature and this second order term is needed at least in some area of the phase diagram to stabilize the theory. Shown in this figure are four slices at constant k¯ in the phase space spanned by the three parameters k, ¯ γ and ¯ The regions colored with blue are the regions where a ’cookie disk’ is metastable k. with the absolute minimum of the free energy corresponding to perfectly flat and circular disks, while the regions colored with yellow are where the ’cookie’ membranes 89 are the stable structure. In the colored zones, the size of the out-of-plane cusps is about half a penetration depth in accord with experiments. ’Cookie disks’ become metastable at the top of the yellow region because the bare line tension is too high to favor such a structure with higher edge-to-area ratio compared to flat disks. When the line tension is lowered, the stable ’cookie’ region is entered from the top in the phase diagram. But when the line tension is further lowered, structures with even larger edge-to-area ratios such as helices may become stable to replace ’cookie disks’. Also, when the Gaussian curvature modulus is too high, structures with more negative Gaussian curvature can appear. Theoretically, we have found that to the right bottom of the ’cookie’ region, one parameter (Az or Ar ) becomes very large, indicating an instability toward other geometrical structures; this corresponds to the white region in the bottom right part of the phase diagram. Comparison ¯ we see that generally, increasing its value between the diagrams with different k, shifts the ’cookie’ region to the right top. Note that the typical orders of magnitude of both k¯ and γ in previous works[68] for the TRs as shown in Chapter 3 fall in the metastable regions. In most of the metastable regions, the energy barrier between the perfect round disks and the ’cookies’ is about several Frank constant which converts to several hundreds of kB T in physical units. Given the fact that these ’cookies’ are first observed in mixture systems and possibly made by merging two disks with opposite handedness which should produce cusps at the edge due to the incompatible directions of the directors at the edge, we conclude that these ’cookie disks’ are not only possible to be seen in the stable region, but also very possible to be observed experimentally inside these metastable regions shown by blue in the phase diagrams according to our theoretical calculations. As for the stable region for the ’cookie disks’, they appear on the scale of k¯ and γ one order of magnitude higher than the value obtained for TRs. But we should point out that, in the recent experiments where the ’cookie disks’ are equilibrium structures for non-mixtures, some evidence 90 suggests a higher k¯ value than the TR systems[88]. We also recall that the critical value we calculated for the line tension γ at the phase transition between flat disks and TRs is one order of magnitude lower than the experimental value. Given the fact that we have switched to the more physical fixed BC in this computation, we argue that the results for the stable region of the ’cookie’ membranes are also in good agreement with experiments. Due to the limitation of our computation resources, we could not explore the whole phase space in every detail. So we cannot exclude the possibility that ’cookie disks’ can be stable or metastable in other regions apart from the regions obtained by us in the phase space. In other words, the ’cookie’ region obtained by our current computations is quite probably an underestimate. Nevertheless, we have proved the existence of the ’cookie disk’ structure in the phase diagram theoretically and the results can explain the experimental observations qualitatively. As a summary, our general theoretical model is applied to the ’cookie’ shape constructed by attaching out-of-plane cusps and radial bulges to the edge of a disk. A MC method is used to compute the free energy of this shape and their relative stability. The ’cookie’ membranes are found to be stable in a region where k¯ and γ are both one order of magnitude higher than the value previously obtained for TRs, and are found to be only metastable compared to flat disks in some parameter space adjacent to the stable region including the region corresponding to our results for TRs. Our theoretical results support the fact that the ’cookies’ can be observed in experiments with the possibility of either stable for unique handedness viruses and metastable for mixtures and our data are in good qualitative agreement with experimental values. 91 5.4 A Future MC Model for the Membranes In the previous part of this work, the method to compute the director field on monolayer membranes and their corresponding relative stability requires a prescribed geometry or shape of the membranes. This method is inefficient because recalculation is needed each time a geometrical parameter of the membrane is varied. We propose a MC model for the membranes where the shape of the monolayer can simultaneously be varied with the director field. The goal of proposing this model is to make possible the reproduction of various phenomena observed in experiments in one numerical MC computation without a priori knowledge of the geometry of the membrane. These phenomena can be both equilibrium as the optimal director field on a disk, or dynamical, as the emergence of ’cookie’ membranes upon coalescing two disks with opposite handedness. There are previous works producing good numerical results where shapes of nematic vesicles are varied simultaneously with the director field[92, 93]. But to our knowledge, currently there are no similar coarse grained numerical calculations on membranes with open edges. As the common method of discretizing a surface, the membrane is triangularized with 6 triangles at each vertex where directors are attached. The open edge of the membrane is characterized with the line tension related to the length of the edge. The boundary condition is that the directors are tangential at the edge. We use the chiral nearest-neighbor pairwise Lebwohl-Lasher model for the Frank free energy of the director field in the following form[94] VijF = −K(ni · nj )2 − q[(ni × nj ) · rij /rij2 ](ni · nj ) (5.4) where i, j are the indices of the neighboring vertices and rij is the vector pointing 92 from the i-th vertex to the j-th vertex. The total Frank free energy is given by the sum over all the nearest neighboring vertices. To prevent the neighboring vertices from colliding, usually a hardcore repulsive potential is assumed between vertices[95]. In our system, the area constraint is also important thus we replace the hardcore repulsion by a Lennard-Jones potential to prevent simultaneously the coalescence of neighboring vertices and the large change in the surface area. The Lennard-Jones potential on each bond connecting neighboring vertices i and j is given by[96] " 12 6 # rijm rijm  VijLJ = ǫ −2 (5.5) rij rij where rijm is the equilibrium bond length connecting the i-th vertex and the j-th vertex and ǫ is the strength of the Lennard-Jones potential. The discretized Gaussian curvature is given on each vertex as[97] 3 (2π − Σθj ) KiG = (5.6) Ai where θj is the angle at the i-th vertex in the j-th triangle sharing the vertex i and Ai is the total area of these triangles. The form of the discretized mean curvature is more complicated[98, 99] 1 X σij Hi = Ni · rji (5.7) σi rij j(i) where Ni is the surface normal at the i-th vertex given by the average over the surface normals of the triangles which share this vertex. σij is the bond length in the dual lattice given by σij = rij (cot θ1 + cot θ2 )/2 (5.8) and θ1 and θ2 are the two angles opposite the bond connecting the i-th and j-th 93 vertices. The area of the dual cell at vertex i, namely σi , is given by 1X σi = σij rij (5.9) 4 j(i) The total free energy, Gaussian curvature and mean curvature are given by the numerical integral over all the triangles. In the MC simulation, each vertex and director is allowed to move once in each MC cycle, the move is accepted or rejected by the common Metropolis algorithm. As a primitive calculation, we used equilateral triangles with a uniform equilibrium size and fixed bonding across the surface. No good results are obtained at this point. But this model may be useful when self-adaptive mesh size[100, 101, 102] and dynamical triangulation[103, 104, 105] are incorporated to model the fluidity of the composing substance and this is left for future works. Chapter Six Conclusion 95 The monolayer membranes formed by fd viruses appear in various shapes when the physical conditions are varied, which include the line tension, the concentration of the viruses in the aqueous solution and the chirality. A theoretical model for the free energy of these Sm-A membranes has been proposed in previous works[66, 71]. The major contributions to the total free energy are the de-Gennes free energy (reduced to Frank free energy) in which the deformation in the director field is considered, the Helfrich free energy involving the bending of the geometrical surfaces, and the edge energy. E-L equations have been derived from the total free energy through a standard variational method and be solved analytically for semi-infinite disks and numerically for finite disks. When the geometry of the membranes becomes complicated, it is less straight- forward to explicitly write down the exact form of the free energy as in the case of disks. We used techniques of differential geometry to rewrite the previous free energy model in a form that in principle can be applied to any prescribed geometry. We have first applied this model to the chiral twisted ribbons where the resulting E-L equations were solved by a differential equation solver in MATLAB. A first order phase transition from semi-infinite disks to TRs is found when the line tension in the theory (corresponding to the concentration of the depletant in experiments) is lowered which preempts the phase transition between semi-infinite disks and finite disks. Our calculations have shown that TRs always have lower free energy per unit area than finite disks. Thus our theory has predicted the vanishing of small disks (R < 5λt ) from the phase diagram which is verified by experiments. Also, the theory has successfully reproduced several phenomena that are observed experimentally, for instance, the size of the TRs near the phase transition, the decreasing of the tilt angle from disks to TRs and the decrease of the pitch near the phase boundary. Second, we have studied two kinds of instability of disks toward TRs. The first 96 kind is related to the phase transition happening spontaneously upon lowering the line tension where TRs grow from the rippling edge of a disk. The second kind is related to the stretching of the membrane using optical traps. For the first kind, we have proposed model shapes for the ripples attached to a flat disk. A MC method has been used to solve the model. Our MC calculations have shown that the ripples at the edge of a disk remain tiny deep inside the disks region and abruptly grow in a narrow range of the line tension near the phase transition. Our theoretical results resemble the process of the phase transition observed in experiments. For the second kind, we have proposed a simple model only involving the Helfrich free energy and treated a ribbon-like structure as a perturbation to a flat semi-infinite membrane. The free energy of this perturbation has been calculated and we have obtained the physical conditions for this quantity to be negative, i.e. this perturbation is favorable. Our results suggest that a TR will form only after certain extent of stretching, agreeing with the experiments. Finally, we have applied our model to the achiral membranes composed of viruses with opposite handedness which appear in ’cookie’ shapes with cusps and radial bulges present at the edge. The theory has correctly predicted the order of magnitude of the size of these cusps and radial bulges. The ’cookie’ membranes are predicted to be stable in a region in the phase space where both k¯ and γ are one order of magnitude higher than the values obtained for disks and TRs. And they have been calculated to be metastable in a region including the space where TRs are believed to reside in. We have argued that, although they are metastable in these regions, they are still highly possible to be observed. This is because they are possibly made from merging of disks with opposite handedness which requires cusps to reduce the director mismatch at the edge; and after these cusps are formed, the ’cookie disks’ are hard to relax to the flat ground state due to the large free energy barrier. While 97 in the stable regions for ’cookies’, it has been expected that even a membrane with unique composition can have these cusp structures and this has been verified recently in experiments. Collecting all these results for different shapes of the membranes, the success of our general theoretical model in reproducing the experimental observations of these systems is quite good. Thus we conclude that this model successfully captures the most important physical principles that govern the behavior of these fd systems and may be applied to other monolayer liquid crystalline membranes. This model may also help predict new structures that have not been observed in experiments yet. At last, we point out that due to the prescribed geometry required by our numerical method, the program is not efficient enough and assumes knowledge of the structure rather than predicting it. We have proposed a method for future works to overcome these deficiencies where the shape and the director field are subject to being varied simultaneously in MC simulation, therefore predicting structures rather than merely explaining the observed ones. Bibliography [1] P. G. deGennes and J. Prost, The Physics of Liquid Crystals (Oxford University Press, 1993) [2] M. J. Stephen and J. P. Straley, Rev. Mod. Phys. 46, 617 (1974). [3] P. M. Chaikin and T. C. Lubensky, Principles of Condensed Matter Physics (Cambridge University Press, 1995) [4] V. Frederiks and V. Zolina, Trans. Faraday Soc. 29, 919 (1933) [5] V. Frederiks and V. Zvetkoff, Sov. Phys. 6, 490 (1934) [6] H. Zocher, Trans. Faraday Soc. 29, 883 (1933) [7] V. Y. Reshetnyak, S. M. Shelestiuk and T. J. Sluckin, Mol. Cryst. Liq. Cryst. 454, 201 (2006) [8] P. Galatola, C. Oldano and M. Rajteri, Phys. Rev. E 49, 1458 (1994) [9] D. V. Makarov and A. N. Zakhlevnykh, Phys. Rev. E 81, 051710 (2010) [10] N. D. Mermin, Rev. Mod. Phys. 51, 591 (1979) [11] Jan P. F. Lagerwall and G. Scalia, Current Applied Physics 12, 1387 (2012) [12] J. Nehring and A. Saupe, J. Chem. Soc., Faraday Trans. 68, 1 (1972) [13] , T. J. Dingemans and E. T. Samulski, Liq. Cryst. 27, 131 (2000) [14] Y. Takanishi, H. Takezoe, A. Fukuda, H. Komura and J. Watanabe, J. Mater. Chem. 2, 71 (1992) [15] P. G. deGennes, Solid State Comm. 10, 753 (1972) [16] M. Kleman and O. D. Lavrentovich, Soft Matter Physics: an Introduction (Springer-Verlag New York Inc. 2003) 98 99 [17] W. Meissner and R. Ochsenfeld, Naturwissenschaften 21, 787 (1933) [18] C. Kittel, Introduction to Solid State Physics (John Wiley and Sons, 2004) [19] E. H. Brandt, Rep. Prog. Phys. 58, 1465 (1995) [20] H. Safar, P. L. Gammel and D. A. Huse, Phys. Rev. Lett. 69, 824 (1992) [21] H. F. Hess, R. B. Robinson, R. C. Dynes, J. M. Valles and J. V. Waszczak, Phys. Rev. Lett. 62, 214 (1989) [22] S. R. Renn and T. C. Lubensky, Phys. Rev. A 38, 2132 (1988) [23] T. C. Lubensky and S. R. Renn, Phys. Rev. A 41, 4392 (1990) [24] J. Toner, Phys. Rev. B 43, 8289 (1991) [25] A. J. Slaney and J. W. Goodby, Liq. Cryst. 9, 849 (1991) [26] S. R. Renn, Phys. Rev. A 45, 953 (1992) [27] J. W. Goodby, A. J. Slaney, C. J. Booth, I. Nishiyama, J. D. Vuijk, P. Styring and K. J. Toyne, Mol. Cryst. Liq. Cryst. 243, 231 (1994) [28] L. Onsager, Ann. N. Y. Acad. Sci. 51, 627 (1949) [29] S. Asakura and F. Oosawa, J. Chem. Phys. 22, 1255 (1954) [30] J. F. Joanny, L. Leibler and P. G. deGennes, J. Poly. Sci.: Poly. Phys. Ed. 17, 1073 (1979) [31] A. G. Yodh, K. H. Lin, J. C. Crocker, A. D. Dinsmore, R. Verma and P. D. Kaplan, Philos. Trans. R. Soc. London A 359 921 (2001) [32] L. Li, D. S. Zhou, J. F. Zhang and G. Xue, J. Phys. Chem. B 108, 5153 (2004) [33] T. Gibaud, E. Barry, M. Zakhary, M. Henglin, A. Ward, Y. Yang, C. Berciu, R. Oldenbourg, D. Nicastro, R. Meyer and Z. Dogic, Nature 481, 348 (2012) [34] A. J. Kim, P. L. Biancaniello and J. C. Crocker, Langmuir 22, 1991 (2006) [35] D. Baranov, A. Fiore, M. van Huis, C. Giannini, A. Falqui, U. Lafont, H. Zandbergen, M. Zanella, R. Cingolani and L. Manna, Nano Lett. 10, 743 (2010) [36] J. C. Love, L. A. Estroff, J. K. Kriebel, R. G. Nuzzo and G. M. Whitesides, Chem. Rev. 105, 1103 (2005) [37] S. G. Zhang, T. Holmes, C. Lockshin and A. Rich, Proc. Natl. Acad. Sci. USA 90, 3334 (1993) [38] J. F. Nagle and S. Tristram-Nagle, Biochim. Biophys. Acta 1469, 159 (2000) [39] D. Picot, P. J. Loll and R. M. Garavito, Nature 367, 243 (1994) [40] S. H. White and W. C. Wimley, Annu. Rev. Biophys. Biomol. Struct. 28, 319 (1999) 100 [41] B. M. Discher, Y. Y. Won, D. S. Ege, J. C. M. Lee, F. S. Bates, D. E. Discher and D. A. Hammer, Science 284, 1143 (1999) [42] D. E. Discher and A. Eisenberg, Science 297, 967 (2002) [43] H. Kellay, B. P. Binks, Y. Hendrikx, L. T. Lee and J. Meunier, Adv. Colloid and Interface Sci. 49, 85 (1994) [44] R. Atkin and G. G. War, J. Phys. Chem. B 111, 9309 (2007) [45] H. Endo, M. Mihailescu, M. Monkenbusch, J. Allgaier, G. Gompper, D. Richter, B. Jakobs, T. Stottmann, R. Strey and I. Grillo, J. Chem. Phys. 115 580 (2001) [46] U. Seifert, Adv. Phys. 46, 13 (1997) [47] S. J. Singer and G. L. Nicolson, Science 175, 720 (1972) [48] R. Goetz, G. Gompper and R. Lipowsky, Phys. Rev. Lett. 82, 221 (1999) [49] J. L. van Hemmen and C. Leibold, Phys. Rep. 44, 51 (2007) [50] C. H. Lin, M. H. Lo and Y. C. Tsai, Prog. Theo. Phys. 109, 591 (2003) [51] L. Miao, B. Fourcade, M. D. Rao, M. Wortis and R. K. P. Zia, Phys. Rev. A 43, 6843 (1991) [52] D. C. Morse and S. T. Milner, Phys. Rev. E 52, 5918 (1995) [53] F. Julicher and U. Seifert, Phys. Rev. E 49, 4728 (1994) [54] J. R. Frank and M. Kardar, Phys. Rev. E 77, 041705 (2008) [55] M. Bowick, D. R. Nelson and A. Travesset, Phys. Rev. E 69, 041102 (2004) [56] P. Biscari and E. M. Terentjev, Phys. Rev. E 73, 051706 (2006) [57] D. P. Landau and K. Binder, A Guide to Monte Carlo Simulations in Statistical Physics (Cambridge University Press, 2005) [58] A. Kaminski and S. Das Sarma, Phys. Rev. Lett. 88, 247202 (2002) [59] P. J. Reynolds, H. E. Stanley and W. Klein, Phys. Rev. B 21, 1223 (1980) [60] G. Szabo and C. Toke, Phys. Rev. E 58, 69 (1998) [61] M. Lattuada, H. Wu and M. Morbidelli, J. Colloid and Interface Sci. 268, 106 (2003) [62] M. J. Saxton, Biophysical J. 64, 1053 (1993) [63] E. Manousakis, Rev. Mod. Phys. 63, 1 (1991) [64] Z. Dogic and S. Fraden, Curr. Opin. Colloid Interface Sci. 11, 47 (2006) [65] E. Barry, D. Beller and Z. Dogic, Soft Matter 5, 2563 (2009) 101 [66] E. Barry, Z. Dogic, R. B. Meyer, R. A. Pelcovits and R. Oldenburg, J. Phys. Chem. B 113, 3910 (2009) [67] E. Barry and Z. Dogic, Proc. Natl. Acad. Sci. USA 107, 10348 (2010) [68] C. N. Kaplan, H. Tu, R. A. Pelcovits and R. B. Meyer, Phys. Rev. E 82, 021701 (2010) [69] W. Helfrich, Z. Naturforsch. C 28, 693 (1973) [70] W. Helfrich and H. J. Deuling, J. Phys. (Paris), Colloq. 36, C1-327 (1975) [71] R. A. Pelcovits and R. B. Meyer, Liq. Cryst. 36, 1157 (2009) [72] C. E. Weatherburn, Differential Geometry of Three Dimensions (the University Press, 1961) [73] Ou-Yang Zhong-can and Liu Ji-xing, Phys. Rev. Lett. 65, 1679 (1990) [74] Ou-Yang Zhong-can and Liu Ji-xing, Phys. Rev. A 43, 6826 (1991) [75] Ou-Yang Zhong-can and W. Helfrich, Phys. Rev. A 39, 5280 (1989) [76] Ou-Yang Zhong-can, Liu Ji-xing and Xie Yu-zhang, Geometric Methods in the Elastic Theory of Membrane in Liquid Crystal Phases (World Scientific Publish- ing Co. Pre. Ltd. 1999) [77] W. Helfrich and J. Prost, Phys. Rev. A 38, 3065 (1988) [78] D. P. Siegel, Biophys. J. 91, 608 (2006) [79] D. P. Siegel and B. G. Tenchov, Biophys. J. 94, 3987 (2008) [80] D. Marsh, Chem. Phys. Lipids 144, 146 (2006) [81] S. A. Safran, Advances in Phys. 48, 395 (1999) [82] private communication with C. N. Kaplan [83] I. Szleifer, D. Kramer, A. Ben-Shaul, W. M. Gelbart and S. A. Safran, J. Chem. Phys. 92, 6800 (1990) [84] G. Gompper and S. Zschocke, Phys. Rev. A 46, 4836 (1992) [85] E. M. Blokhuis and D. Bedeaux, Mol. Phys. 80, 705 (1993) [86] S. M. Oversteegen and E. M. Blokhuis, J. Chem. Phys. 112, 2980 (2000) [87] F. A. M. Leermakers, P. A. Barneveld, J. Sprakel and N. A. M. Besseling, Phys. Rev. Lett. 97, 066103 (2006) [88] private communication with Z. Dogic [89] private communication with R. B. Meyer 102 [90] T. Shemesh, A. Luin, V. Malhotra, K. N. J. Boerger and M. M. Kozlov, Biophys. J. 85, 3813 (2003) [91] M. D. Mitov, Dokl. Bulg. Akad. Nauk 31, 513 (1978) [92] H. Shin, M. J. Bowick and X. Xing, Phys. Rev. Lett. 101, 037802 (2008) [93] X. Xing, H. Shin, M. J. Bowick, Z. Yao, L. Jia and M. Li, Proc. Natl. Acad. Sci. USA 109, 5202 (2012) [94] R. Memmer and O. Fliegans, Phys. Chem. Chem. Phys. 5, 558 (2003) [95] M. P. Allen and D. J. Tildesley, Computer Simulation of Liquids (Clarendon, Oxford, 1987) [96] J. E. Lennard-Jones, Proc. R. Soc. Lond. A 106, 463 (1924) [97] G. Xu, Comput. Aided Geom. Des. 23, 193 (2006) [98] C. Itzykson, Proc. GIFT Seminar, Jaca 85, eds. J. A. Abad, M. Asorey and A. Cruz, (World Scientific, Singapore, 1986) [99] D. Nelson, T. Piran and S. Weinberg, Statistical Mechanics of Membranes and Surfaces (2nd Edition) (World Scientific, Singapore, 2004) [100] J. Fukuda, H. Stark, M. Yoneya and H. Yokoyama, Phys. Rev. E 69, 041706 (2004) [101] M. Tasinkevych, N. M. Silvestre, P. Patricio and M. M. T. da Gama, Eur. Phys. J. E 9, 341 (2002) [102] P. Patricio, M. Tasinkevych and M. M. T. da Gama, Eur. Phys. J. E 7, 117 (2002) [103] P. B. S. Kumar and M. Rao, Phys. Rev. Lett. 80, 2489 (1998) [104] H. Koibuchi, N. Kusano, A. Nidaira, K. Suzuki and M. Yamada, Phys. Lett. A, 319, 44 (2003) [105] H. Koibuchi, Phys. Rev. E 75, 051115 (2007)