Curvature-driven stability of defects in nematic textures over spherical disks

TL;DR

This paper analytically investigates how curvature influences the stability and arrangement of defects in nematic liquid crystal textures on spherical disks, revealing phase-transition-like behavior and defect interactions driven by curvature effects.

Contribution

It introduces a curvature-driven stability mechanism for nematic defects on spherical disks, showing how curvature affects defect positions and interactions without topological constraints.

Findings

Curvature prevents +1 and +1/2 defects from forming boojum textures.

A narrow area window causes +1 defect to move from boundary to center, showing phase transition behavior.

Curvature induces alternating repulsive and attractive interactions between +1/2 defects.

Abstract

Stabilizing defects in liquid-crystal systems is crucial for many physical processes and applications ranging from functionalizing liquid-crystal textures to recently reported command of chaotic behaviors of active matters. In this work, we perform analytical calculations to study the curvature driven stability mechanism of defects based on the isotropic nematic disk model that is free of any topological constraint. We show that in a growing spherical disk covering a sphere the accumulation of curvature effect can prevent typical +1 and +1/2 defects from forming boojum textures where the defects are repelled to the boundary of the disk. Our calculations reveal that the movement of the equilibrium position of the +1 defect from the boundary to the center of the spherical disk occurs in a very narrow window of the disk area, exhibiting the first-order phase-transition-like behavior. For the pair of +1/2 defects by splitting a +1 defect, we find the curvature driven alternating repulsive and attractive interactions between the two defects. With the growth of the spherical disk these two defects tend to approach and finally recombine towards a +1 defect texture. The sensitive response of defects to curvature and the curvature driven stability mechanism demonstrated in this work in nematic disk systems may have implications towards versatile control and engineering of liquid crystal textures in various applications.

Equations

$$F={\int }_{D}fdA+\lambda ({\mathbf{n}}^2-1),$$

$$f=\frac{1}{2}{K}_1(\text{div}\mathbf{n}{)}^2+\frac{1}{2}{K}_3(\mathbf{n}\times \text{curl}\mathbf{n}{)}^2,$$

$${\partial }_j\frac{\partial f}{\partial (\frac{\partial n_i}{\partial u^j})}+\frac{{\partial }_j\sqrt{g}}{\sqrt{g}}\frac{\partial f}{\partial (\frac{\partial n_i}{\partial u^j})}-\frac{\partial f}{\partial n_i}=-\lambda n_i,$$

$${\oint }_{\Gamma }d\theta =k\pi ,$$

$$\mathbf{n}=\cos (\phi +{\theta }_0){\mathbf{e}}_1+\sin (\phi +{\theta }_0){\mathbf{e}}_2,$$

$$F_{+1,p}=\frac{K_1}{2}{\iint }_{x^2+y^2\le r_p^2}\frac{1}{x^2+y^2}dxdy,$$

$$F_{+1,p}(c)=\frac{K_1}{2}{\iint }_{x^2+y^2\le r_p^2}\frac{1}{(x-c{)}^2+y^2}dxdy.$$

$$F_{+1,p}^{\prime }(c)=K_1{\iint }_{x^2+y^2\le r_p^2}\frac{x-c}{[(x-c{)}^2+y^2{]}^2}dxdy$$

$$=K_1{\iint }_{(x^{\prime }+c{)}^2+y^2\le r_p^2}\frac{x^{\prime }}{(x^{\prime 2}+y^2{)}^2}d{x}^{\prime }dy,$$

$$F_{+1/2,p}(c)=2K_1{\int }_0^{r_p-c}dx{\int }_0^{\sqrt{r_p^2-(x+c{)}^2}}dy\frac{1}{x^2+y^2},$$

$$(n_r,n_{\phi })=\left(\text{sin}\theta (\alpha ),\text{cos}\theta (\alpha )\right),$$

$$F_{\mathrm{spiral}}(\alpha )=\left(1+{\left(\frac{\alpha }{\text{ln}\frac{r_p}{r_0}}\right)}^2\right)F_{+1,p},$$

$$\text{div}\mathbf{n}=\frac{\text{cos}\theta }{\text{sin}\theta }\frac{n_1}{R}+\frac{1}{R}\frac{\partial n_1}{\partial \theta }+\frac{1}{R\text{sin}\theta }\frac{\partial n_2}{\partial \phi },$$

$$\text{curl}\text{n}=\frac{1}{R\text{sin}\theta }(-\frac{\partial n_1}{\partial \phi }+n_2\text{cos}\theta +\text{sin}\theta \frac{\partial n_2}{\partial \theta }){\mathbf{e}}_{\mathbf{r}}$$

$$-\frac{n_2}{R}{\mathbf{e}}_{\theta }+\frac{n_1}{R}{\mathbf{e}}_{\phi },$$

$$f$$

$$F_{+1,s}$$

$$(u,v)=(\frac{Rx}{z+R},\frac{Ry}{z+R}),$$

$$(u,v)=(\frac{R\text{sin}\theta \text{cos}\phi }{\text{cos}\theta +1},\frac{R\text{sin}\theta \text{sin}\phi }{\text{cos}\theta +1}).$$

$$(u-u_0{)}^2+v^2\le r_{eq}^2,$$

$$u_0=\frac{2c{R}^2}{-r_s^{\prime 2}+2R(R+\sqrt{R^2-c^2})},$$

$$r_{eq}^2=\frac{r_s^{\prime 2}{R}^2(4R^2-r_s^{\prime 2})}{[r_s^{\prime 2}-2R(R+\sqrt{R^2-c^2}){]}^2}.$$

$$\displaystyle\frac{\partial(u,v)}{\partial(\theta,\varphi)}=\left(\begin{array}[]{cc}\frac{R\textrm{cos}\varphi}{1+\textrm{cos}\theta}&-\frac{R\textrm{sin}\theta\textrm{sin}\varphi}{1+\textrm{cos}\theta}\\ \frac{R\textrm{sin}\varphi}{1+\textrm{cos}\theta}&\frac{R\textrm{sin}\theta\textrm{cos}\varphi}{1+\textrm{cos}\theta}\\ \end{array}\right),$$

$$dudv=\left|\frac{\partial (u,v)}{\partial (\theta ,\phi )}\right|d\theta d\phi =\frac{u^2+v^2}{\text{sin}\theta }d\theta d\phi ,$$

$$\frac{dudv}{u^2+v^2}=\frac{d\theta d\phi }{\text{sin}\theta }.$$

$$F_{+1,s}=\frac{K_1}{2}{\iint }_D\frac{1}{u^2+v^2}dudv,$$

$$\Delta F$$

$$F_{+1,s}(c)=\frac{K_1}{2}{\iint }_{x^2+y^2\le r_{eq}^2}\frac{1}{(x+u_0{)}^2+y^2}dxdy$$

$$=K_1{\int }_{-r_{eq}}^{r_{eq}}dx\frac{1}{x+u_0}\text{arctan}\frac{\sqrt{r_{eq}^2-x^2}}{x+u_0},$$

$$F_{+1,s}^{\prime }(c)=K_1{\int }_{-r_{eq}}^{r_{eq}}(G_1+G_2+G_3)dx,$$

$$\frac{F_{+\frac{1}{2},s}}{2K_1}={\iint }_{D_1}\frac{1}{x^2+y^2}dA+{\iint }_{D_2}\frac{1}{R^2-y^2}dA,$$

$$\frac{F_{+\frac{1}{2},s}^{\prime }(c)}{2K_1}={\int }_0^{\frac{b-c}{R}}{F}_1(x)dx+{\int }_0^{\sqrt{b^2-c^2}}{F}_2(y)dy,$$

Full text

I Introduction

Functionalizing the rich variety of self-assembled liquid-crystal (LC) structures represents a trend in LC research Bisoyi and Kumar (2011); Alexander et al. (2012); Umadevi et al. (2013); Urbanski et al. (2017). Confining LCs in various geometries in the form of droplets Lopez-Leon et al. (2011); Pairam et al. (2013); Yamamoto and Sano (2016), shells Fernández-Nieves et al. (2007); Lopez-Leon and Fernandez-Nieves (2009); Liang et al. (2012) and fibers Yoshino et al. (2010); Fleischmann et al. (2014) using modern microfluidic technology and characterization methods opens the prospect of many application opportunities, and brings new scientific problems related to the creation and engineering of complex director arrangements de Gennes and Prost (1993); Lopez-Leon and Fernandez-Nieves (2011); Sengupta et al. (2012); Urbanski et al. (2017). LC textures can be strongly affected by the distribution and type of topological defects, which are singularities in the otherwise continuous LC director field Kleman (1983); Chaikin and Lubensky (2000); Nelson (2002a). The extraordinary responsiveness of LC makes the manipulation of defects a challenge in applications. Stabilizing defects in two-dimensional LC systems is directly related to arrangement of LC textures Kleman (1983); de Gennes and Prost (1993); Mbanga et al. (2014); Darmon et al. (2016), fabrication of controllable valency in colloid-LC-based artificial atoms Nelson (2002b); Gharbi et al. (2013); Koning et al. (2016), modulation of coupled geometries where LC lives Frank and Kardar (2008); Xing et al. (2012); Chen et al. (2013); Ramakrishnan et al. (2013); Pismen (2014); Mostajeran (2015); Leoni et al. (2017), and relevant applications in active matter systems Sanchez et al. (2012); Zhou et al. (2014); Keber et al. (2014); Lavrentovich (2016); Peng et al. (2016). A prototype model to study the stability mechanism of defects in LC is the isotropic two-dimensional LC disk model with a single elastic constant Langer and Sethna (1986); Rudnick and Bruinsma (1995); Pettey and Lubensky (1999). In a flat freestanding LC disk, defects tend to move swiftly to the boundary to form a boojum texture, which is a two-dimensional version of its namesake in superfluid helium-3 Bhattacharyya et al. (1977); Stein et al. (1978); Langer and Sethna (1986). A “virtual boojum” texture with a topological defect outside the sample has been predicted in planar circular LC domains by Langer and Sethna  Langer and Sethna (1986), and it has been found to be a local energy extremal Rudnick and Bruinsma (1995); Pettey and Lubensky (1999). Sufficiently strong pinning boundary conditions can stabilize a defect within a circular LC domain Langer and Sethna (1986); Rudnick and Bruinsma (1995); Pettey and Lubensky (1999).

Exploring other stability mechanisms of defects in LC samples in addition to imposing boundary conditions constitutes an underlying scientific problem towards versatile control and engineering of LC director arrangement. Confining LC over spherical surfaces can generate various regularly arranged stable defect patterns Vitelli and Nelson (2006); Fernández-Nieves et al. (2007); Shin et al. (2008); Bates (2008); Bowick and Giomi (2009); Xing et al. (2012); Seč et al. (2012); Koning et al. (2016). Vitelli and Nelson have studied two-dimensional nematic order coating frozen surfaces of spatially varying Gaussian curvature, and found the instability of a smooth ground-state texture to the generation of a single defect using free boundary conditions Vitelli and Nelson (2004). These results of LC order on closed spheres and topographies with varying curvature show that curvature suffices to provide a stability mechanism for defects even without imposing any pinning boundary condition. However, for LC order on a closed sphere, it is unknown to what extent the appearance of defects is energetically driven, while they must appear as a consequence of the spherical topology. To remove the topological constraint, we study nematic order, the simplest LC order, on a spherical disk. Here we emphasize that, due to the fundamentally distinct topologies of sphere and disk, the appearance of defects on spherical disks is not topologically required; the emergence of defects therein is purely geometrically driven. According to the continuum elasticity theory of topological defects in either LC or crystalline order, the stress caused by defects can be partially screened by Gaussian curvature Chaikin and Lubensky (2000); Nelson (2002a); Bowick and Giomi (2009); Koning and Vitelli (2016). Therefore, one expects the appearance of defects on a sufficiently curved spherical disk. It is of interest to identify the transition point for a defect to depart from the boundary of the disk, and illustrate the nature of the transition by clarifying questions such as: Will the defect move rapidly or gradually with the accumulation of curvature effect? Will the defect split as the nematic texture becomes more and more frustrated by the curvature? Once split, will the resulting defects become stable on the spherical disk?

We perform analytical calculations based on the isotropic nematic disk model to address these fundamental problems. This theoretical model may be realized experimentally in Langmuir monolayers Knobler and Desai (1992); Fang et al. (1997); Pettey and Lubensky (1999); Gupta and Manjuladevi (2012) and liquid-crystal films Langer and Sethna (1986); Xue et al. (1992) deposited at the surface of water droplets whose curvature is controllable by tuning the droplet size Urbanski et al. (2017). Flat space experiments in these two-dimensional monolayer systems at air-water interface have revealed stable liquid-crystal phases Knobler and Desai (1992); Xue et al. (1992); Gupta and Manjuladevi (2012).

In this work, we first discuss the two instability modes of a +1 defect over a flat disk, either sliding to the boundary or splitting to a pair of +1/2 defects. By depositing the nematic order over a spherical surface, we analytically show that bending deformation of a director field is inevitable everywhere, which implies the appearance of defects to release the curvature-driven stress. By comparing a flat and a spherical nematic disk of the same area, both containing a +1 defect at the center, we derive for the analytical expression for the difference of the Frank free energy, and show that the spherical disk always has higher energy. However, when the +1 defect deviates from the center of the disk, the free energy curves become qualitatively different for flat and spherical disks when the disk area exceeds some critical value. Specifically, the equilibrium position of the +1 defect rapidly moves from the boundary to the center of the spherical disk in a narrow window of the disk area, exhibiting the first-order phase-transition-like behavior. For the pair of +1/2 defects by splitting a +1 defect, we further show the curvature-driven alternating repulsive and attractive interactions between the two defects. When the spherical cap occupies more area over the sphere, the pair of +1/2 defects tend to approach until merging to a +1 defect texture. The recombination of the pair of +1/2 defects into a +1 defect is consistent with the result of the +1 defect case. These results demonstrate the fundamentally distinct scenario of defects in a spherical disk from that on a planar disk. We also briefly discuss the cases of nematic order on hyperbolic disks. In this work, the demonstrated distinct energy landscape of LC defects created by curvature is responsible for the stability of defects, and may have implications in the design of LC textures with the dimension of curvature.

II Model and Method

In the continuum limit, the orientations of liquid-crystal molecules lying over a disk are characterized by a director field $\mathbf{n}(\mathbf{x})$ that is defined at the associated tangent plane at $\mathbf{x}$. The equilibrium nematic texture is governed by minimizing the Frank free energy Chaikin and Lubensky (2000)

[TABLE]

where the integration is over the disk $D$. The Frank free energy density

[TABLE]

where $K_{1}$ and $K_{3}$ are the splay and bending rigidities, respectively. The Lagrange multiplier $\lambda$ is introduced to implement the constraint of $\mathbf{n}\cdot\mathbf{n}=1$. In general, $\lambda$ is a function of coordinates. The twist term $(\mathbf{n}\cdot\textrm{curl}\,\mathbf{n})^{2}$ vanishes in nematics confined on a sphere (see Appendix B). Equation (2) has been widely used to analyze the deformation in nematic phases. For nematics on curved surfaces, the operators of divergence and curl in the Frank free energy are promoted to be defined on the curved manifold and carry the information of curvature. Note that the curl operator relies on the extrinsic geometry of the surface  Napoli and Vergori (2012). Note that the Frank free energy model in Eq. (2) describes the distortion free energy of uniaxial nematics. A formalism based on the tensorial nematic order parameter has been proposed to characterize the distortion of both uniaxial and biaxial nematics and defects therein Kralj et al. (2011); Rosso et al. (2012).

We work in the approximation of isotropic elasticity with $K_{1}=K_{3}$. Under such an approximation, one can show that the free energy is invariant under the local rotation of the director field by any angle, whether the disk is planar or curved (see Appendix A). In other words, the energy degeneracy of the system becomes infinite when $K_{1}=K_{3}$. Such configurational symmetry is broken when the ratio $K_{1}/K_{3}$ is deviated from unity. While the states selected by the differential in the values for $K_{1}$ and $K_{3}$ are of interest in other contexts such as in the ground states of spherical nematics Bowick and Giomi (2009), here we work in the isotropic regime to highlight the curvature effect of substrates on the configuration of nematics.

The general Euler-Lagrange equation of the Frank free energy on a curved surface $\mathbf{x}(u^{1},u^{2})$ is

[TABLE]

where $i,j=1,2$, and $g$ is the determinant of the metric tensor. The second term in Eq.(3) is due to the spatially varying $g$. The nematic textures studied in this work are solutions to Eq.(3).

To characterize defects that are named disclinations in a two-dimensional director field, we perform integration of the orientation $\theta$ of the director $\mathbf{n}$ with respect to any local reference frame along any closed loop $\Gamma$:

[TABLE]

where $k$ is nonzero if $\Gamma$ contains a defect. Unlike in a vector field where $k$ can only be integers, two-dimensional nematics supports both integer and half-integer disclinations due to the apolarity of liquid-crystal molecules, i.e., $\mathbf{n}\equiv-\mathbf{n}$.

III Results and discussion

We first discuss the case of nematics on a planar disk. It is straightforward to identify the following solution to the Euler-Lagrange equation:

[TABLE]

where $\varphi=\textrm{arctan}(y/x)$ is the polar angle, $\theta_{0}$ is a constant, and $\mathbf{e}_{i}$ is the unit basis vector in Cartesian coordinates. The strength of the defect located at the origin of the coordinates is +1. The associated Lagrange multiplier is $\lambda=K_{1}/(x^{2}+y^{2})$. The contributions to the splay and bending terms in the free energy density are $K_{1}\textrm{cos}^{2}\theta_{0}/[2(x^{2}+y^{2})]$ and $K_{3}\textrm{sin}^{2}\theta_{0}/[2(x^{2}+y^{2})]$, respectively. When $\theta_{0}$ increases from 0 to $\pi/2$, the +1 defect transforms from the radial (pure splay) to the azimuthal (pure bending) configurations. In this process, the sum of the splay and bending energies is an invariant under the isotropic elasticity approximation. The total free energy of the configuration in Eq.(5) is

[TABLE]

where $r_{p}$ is the radius of the planar disk.

We show that the +1 defect at the center of the planar disk in Eq.(5) is unstable and tends to slide to the boundary of the disk. For simplicity, we employ free boundary condition. Consider a +1 defect like in Eq.(5) at $(c,0)$, where $c\leq r_{p}$. Its free energy is

[TABLE]

To avoid the singularity point at $(c,0)$ in the evaluation for $F_{+1,p}(c)$, we take the derivative of $F_{+1,p}(c)$ with respect to $c$. Physically, this procedure returns the force on the defect. While the free energy may diverge, a physical force must be finite. After some calculation, we have

[TABLE]

where variable substitution is applied in the last equality. The integral domain is shown in Fig. 1(a). The defect is located at $x^{\prime}=0$ (i.e., $x=c$) and $y=0$. We see that the integration in the red region returns zero, since the integrand $x^{\prime}/(x^{\prime 2}+y^{2})^{2}$ is an odd function of $x^{\prime}$. In the rest region where $x^{\prime}<0$, the integrand is negative. Therefore, $F^{\prime}_{+1,p}(c)$ is negative when the defect is deviated from the center of the disk. $F^{\prime}_{+1,p}(c=0)=0$. In other words, once deviated from the center of the disk, the defect will slide to the boundary to reduce the free energy of the system. Figure 1(b) shows the numerical result on the dependence of $F^{\prime}_{+1,p}(c)$ on $c/r_{p}$.

An alternative instability mode of the central +1 defect in the planar disk is to split into two +1/2 defects. Such a process may occur when the interaction energy of the two repulsive +1/2 defects dominates over the core energy of the defects. To analyze the energetics of the +1/2 defects, we construct the director field containing two +1/2 defects by cutting and moving apart an azimuthal configuration as shown in Fig.1(c), where the +1/2 defects are represented by red dots. The region between the two half azimuthal configurations is filled with a uniform director field. The Frank free energy of such a configuration is

[TABLE]

where the separation between the two defects is $2c$. $F^{\prime}(0)=-2K_{1}/r_{p}<0$. In Fig. 1(d), we plot $F^{\prime}_{+1/2,p}(c)$ versus $c$. The negative sign indicates the repulsive nature of the two +1/2 defects. The resulting +1/2 defects are ultimately pushed to the boundary of the disk under the repulsive interaction.

In the preceding discussions, we employ the free boundary condition where directors at the boundary do not have preferred orientations. Another important class of boundary condition is to fix the orientation of the molecules at the boundary. Homeotropic and planar liquid-crystal samples are two typical cases, where the directors are perpendicular and parallel to the boundary, respectively. Imposing these pinning boundary conditions over the aster configuration can lead to spiral deformations de Gennes and Prost (1993). Note that a recent study has demonstrated a dynamic consequence of the radial-to-spiral transition of a +1 defect pattern in the system of swimming bacteria in a liquid-crystal environment Peng et al. (2016). It is observed that the swimming mode of bacteria changes from bipolar to unipolar when the +1 defect pattern becomes spiral. For the general pinning boundary condition that the angle between $\mathbf{n}$ and the tangent vector at the boundary is $\alpha$ ($\alpha\in[0,\pi/2]$), we obtain the solution to Eq.(3):

[TABLE]

where $n_{r}$ and $n_{\varphi}$ are the components of $\mathbf{n}$ in polar coordinates $(r,\varphi)$, $\theta(\alpha)=\alpha\ \textrm{ln}(r/r_{0})/\textrm{ln}(r_{p}/r_{0})$, $r_{p}$ and $r_{0}$ are the outer and inner radius of the planar disk as shown in Fig. 2. The magic spiral solution in Ref. de Gennes and Prost (1993) is a special case of $\alpha=\pi/2$. The associated Lagrange multiplier is $\lambda=(K_{1}/r^{2})\{1+\alpha^{2}/[\textrm{ln}(r_{p}/r_{0})]^{2}\}$. The configuration of the solution in Eq.(10) is plotted in Fig. 2. The originally straight radial lines deform to spiral curves to satisfy the boundary condition. The Frank free energy of the spiral configuration is

[TABLE]

where $F_{+1,p}$ is the free energy of an aster configuration in a planar disk given in Eq.(6). Eq.(11) shows that the boundary effect does not enter the integral of $F_{+1,p}$. The energy cost associated with the spiral deformation conforms to a quadratic law with respect to the angle $\alpha$. And its dependence on the size of disk is relatively weak in a logarithm relation.

Now we discuss two-dimensional nematic texture confined on spherical disks. Consider a director field $\mathbf{n}$ on a sphere $\mathbf{n}=n_{1}(\theta,\varphi)\mathbf{e}_{\theta}+n_{2}(\theta,\varphi)\mathbf{e}_{\varphi}$, where $\mathbf{e}_{\theta}$ and $\mathbf{e}_{\varphi}$ are the unit tangent vectors in spherical coordinates. $\theta$ and $\varphi$ are the polar and azimuthal angles, respectively. We first show that on spherical geometry a director field without any splay and bending deformations is impossible. Topology of the two-dimensional sphere dictates that a harmonic vector field on a sphere is impossible Nakahara (2003). A vector field is called harmonic if it is divergence-free, irrotational, and tangent to the spherical surface. A director field is a vector field with the extra constraints of $|\mathbf{n}|=1$ and $\mathbf{n}\equiv-\mathbf{n}$. Therefore, it is a topological requirement that one cannot completely eliminate both bending and splay deformations in a director field living on a sphere.

In addition to the above global analysis, we will further show that an irrotational director field is impossible at any point on a sphere. In other words, bending of a director field is inevitable everywhere on a spherical surface. We first present the general expressions for the divergence and curl of a director field over a smooth surface: $\textrm{div}\ \mathbf{n}=\frac{1}{\sqrt{g}}\partial_{i}(\sqrt{g}n_{i})$ and $\textrm{curl}\ \textbf{n}=(\star d\textbf{n}^{\flat})^{\sharp}$, where $\star$ is the Hodge dual, $\scriptstyle\flat$ and $\scriptstyle\sharp$ are the musical isomorphisms, $d$ is exterior derivative (see Appendix B). Applying these expressions on a sphere, we have

[TABLE]

and

[TABLE]

where $\mathbf{e_{r}}$ is the unit normal vector. According to Eq.(13), we clearly see that at least one of the last two terms must be nonzero. In contrast, Eq.(12) shows that a divergence free director field with vanishing splay deformation without any bend deformation is possible. The simplest example is the direction field with only the azimuthal component: $\mathbf{n}=\mathbf{e}_{\varphi}$. Such a director field is divergence free but with bending deformation. Note that in the calculation for the curl of the director field, we use the condition that the sphere is embedded in three-dimensional Euclidean space. The divergence of the director field does not depend on how the sphere is embedded in the Euclidean space. One can check that the twist term $(\mathbf{n}\cdot\textrm{curl}\,\mathbf{n})^{2}=0$.

The stability analysis of defects in nematic textures over spherical disks is based on the following expression for the Frank free energy density in spherical coordinates:

[TABLE]

Note that the first term $K_{3}/(2R^{2})$ in the bending part represents the irremovable bending deformation of a director field over spherical substrates. This term vanishes in the limit of $R\rightarrow\infty$. One can check that for a divergence-free director field $\mathbf{n}=\mathbf{e}_{\varphi}$, $f=K_{1}/(2R^{2}\sin^{2}\theta)$. The singularities at $\theta=0$ and $\theta=\pi$ correspond to the two +1 defects at the north and south poles.

We first discuss if the +1 defect can be supported by spherical geometry. All the degenerate nematic configurations containing a +1 defect at the center of the spherical cap are characterized by the director field $\mathbf{n}=(c_{1},c_{2})$, where the constants $c_{1}$ and $c_{2}$ satisfy $c_{1}^{2}+c_{2}^{2}=1$. These degenerate states have the same Frank free energy:

[TABLE]

where the integration is over a spherical cap $D$ with spherical radius $R$ and geodesic radius $r_{s}$. And these states are solutions to the Euler-Lagrange equation (see Appendix C).

In order to derive for $F_{+1,s}-F_{+1,p}$, the free energy difference of a +1 defect configuration on spherical and planar disks, we introduce the following coordinates transformation. For generality, the Cartesian coordinates of the center of the spherical cap are $(c,0,\sqrt{R^{2}-c^{2}})$ as shown in Fig. 3. The center of the spherical cap is located at the north pole for $c=0$. The region of the spherical cap is $D=\{(x,y,z)|x^{2}+y^{2}+z^{2}=R^{2},(x-c)^{2}+y^{2}+(z-\sqrt{R^{2}-c^{2}})^{2}\leq r_{s}^{\prime 2}\}$. $r_{s}^{\prime}$ is the Euclidean distance from the center to the boundary of the spherical cap. The area of such a spherical cap is $S=\pi r_{s}^{\prime 2}$. Now we construct the stereographic projection from the spherical cap to the plane of equator. Specifically, we draw a line connecting the south pole of the sphere and any point at $(x,y,z)$ or $(\theta,\varphi)$ on the spherical cap. The point on the spherical cap is thus projected to the intersection point $(u,v)$ of this line and the equator plane. The projection is described by the formula

[TABLE]

or, in terms of spherical coordinates,

[TABLE]

The stereographic projection has a convenient geometric property that any spherical cap not containing the point of projection (south pole) is projected to a circular disk on the equator plane:

[TABLE]

where

[TABLE]

and

[TABLE]

To guarantee that the spherical cap contains the north pole, it is required that $r_{s}^{\prime 2}\geq 2R(R-\sqrt{R^{2}-c^{2}})$. Alternatively, $c\leq r_{s}^{\prime}\sqrt{1-[r_{s}^{\prime}/(2R)]^{2}}$ for given $r_{s}^{\prime}$. On the other hand, the spherical cap occupies no more than half of a sphere, so $r_{s}^{\prime}\leq\sqrt{2}R$.

From the Jacobian of the coordinates transformation in Eq.(17)

[TABLE]

and

[TABLE]

we finally have

[TABLE]

We therefore obtain the desired expression for Eq.(15) in the $(u,v)$ coordinates:

[TABLE]

where the integral domain $D=\{(u,v)|u^{2}+v^{2}\leq r_{s}^{\prime 2}R^{2}/(4R^{2}-r_{s}^{\prime 2})\}$. Note that now the integrands in Eq.(20) and Eq.(6) have the same functional form and can be conveniently compared. A subtle point worth mentioning is that the direct subtraction of Eq.(6) from Eq.(20) will lead to a wrong expression of $\Delta F=F_{+1,s}-F_{+1,p}=-(\pi/2)K_{1}\ln[4-(r_{s}^{\prime}/R)^{2}]$. One can check that $\Delta F$ fails to converge to the expected zero in the limit of $R\rightarrow\infty$. Here, the subtlety is from the fact that the integrands in Eq.(20) and Eq.(6) have singularity at the origin point. To eliminate this singularity, one has to cut off the small defect core. The integral domain of Eq.(20) should be $D=\{(u,v)|(a/2)^{2}\leq u^{2}+v^{2}\leq r_{s}^{\prime 2}R^{2}/(4R^{2}-r_{s}^{\prime 2})\}$, where $a$ is the radius of the defect core. The prefactor of $1/2$ is due to the shrink of the defect size in the previously introduced stereographic projection. The integral domain in Eq.(6) also becomes $a^{2}\leq x^{2}+y^{2}\leq r_{p}^{2}$. To conclude, the change of the total free energy in the deformation of the planar to the spherical nematic disk in the constraint of fixed disk area $A_{d}$ is

[TABLE]

We check that $\Delta F$ approaches zero in the limit of $R\rightarrow\infty$, as expected. Equation (21) shows that $F_{+1,s}$ is always larger than $F_{+1,p}$.

However, it will be shown that a +1 defect can be stabilized within a sufficiently curved spherical disk despite the higher energy in comparison with the planar disk case. We analyze the stability of the +1 defect from the derivative of the free energy with respect to its position in the disk. The expression for the free energy is rewritten in the new coordinates $\{x,y\}$, where $x=u-u_{0}$ and $y=v$:

[TABLE]

where $u_{0}$ and $r_{eq}$ are given in Eq.(18). From Eq.(22), we have

[TABLE]

where $G_{1}=-[u_{0}^{\prime}(c)/(x+u_{0})^{2}]\textrm{arctan}[\sqrt{r_{eq}^{2}-x^{2}}/(x+u_{0})]$, $G_{2}=-u_{0}^{\prime}(c)\sqrt{r_{eq}^{2}-x^{2}}/[(x+u_{0})(u_{0}^{2}-2u_{0}x+r_{eq}^{2})]$, and $G_{3}=r_{eq}r_{eq}^{\prime}(c)/[(u_{0}^{2}-2u_{0}x+r_{eq}^{2})\sqrt{r_{eq}^{2}-x^{2}}]$. The $G_{3}$ term can be integrated out: $K_{1}\int_{-r_{eq}}^{r_{eq}}G_{3}dx=K_{1}\pi r_{eq}r_{eq}^{\prime}(c)/|u_{0}^{2}-r_{eq}^{2}|$. Local analysis around the defect at $x=-u_{0}$ shows that both the $G_{1}$ and the $G_{2}$ terms are odd functions of $x$, and can be canceled in the integration of $x$ near the defect. The singularity associated with the defect is therefore removed. Note that $F^{\prime}_{+1,s}(c)$ is negative in the large $R$ limit, which is consistent with the planar disk case.

Now we analyze zero points of $F^{\prime}_{+1,s}(c)$. The defect is stable at a zero point where the slope of the $F^{\prime}_{+1,s}(c)$ curve is positive. With the increase of $c$, numerical analysis shows that the $G_{1}$ term decreases and the $G_{3}$ term increases, both starting from zero at $c=0$. While the $G_{1}$ and the $G_{3}$ terms are comparable, the $G_{3}$ term is much smaller than either of them. The competition of the $G_{1}$ and the $G_{3}$ terms may lead to another zero point at the $F^{\prime}_{+1,s}(c)$ curve in addition to the unstable zero point at $c=0$.

In Figs. 4(a)–4(c), we plot $F^{\prime}_{+1,s}(c)$ versus $c$ at typical values for $r_{s}^{\prime}$. We see that the $F^{\prime}_{+1,s}(c)$ is negative and monotonously decreasing when the spherical cap is smaller than a critical value. With the increase of $r_{s}^{\prime}$, a second zero point appears at $c=c^{*}$, where a perturbed defect will be restored to the original equilibrium position. It indicates that the equilibrium position of the defect starts to depart from the boundary of the disk. We introduce the quantity $r_{d}^{\prime}/r_{s}^{\prime}$ to characterize the equilibrium position of the defect over the spherical cap, where $r_{d}^{\prime}$ is the Euclidean distance between the center of the disk and the defect. The variation of the optimal position of the +1 defect with the size of the spherical cap is summarized in Fig. 4(d). A pronouncing feature of the $r_{d}^{\prime}/r_{s}^{\prime}$ vs $r_{s}^{\prime}/R$ curve is the rapid decrease from unity to zero when $r_{s}^{\prime}/R$ varies by only about $0.1\%$. It corresponds to the movement of the defect from the boundary to the center of the disk. Such a transition occurs in the narrow window of $r_{s}^{\prime}$ when the spherical cap occupies about half of the sphere. Note that the spherical cap becomes a half sphere when $r_{s}^{\prime}=\sqrt{2}R$.

Here, it is of interest to compare a +1 defect in nematics and a five-fold disclination in a two-dimensional hexagonal crystal on a sphere. Both nematic and crystalline order are frustrated on a sphere, leading to the proliferation of defects. The resulting defects in condensed matter orders are to screen the geometric charge of the substrate surface, which is defined to be the integral of Gaussian curvature. Over a spherical crystal, the topological charge of a five-fold disclination can be screened by a spherical cap of area $A_{0}/12$ ($A_{0}$ is the area of sphere), since 12 five-fold disclinations are required over a spherical crystal by topological constraint Nelson (2002a). Topological analysis of a spherical nematics shows that a sphere can support two +1 defects, so the topological charge of a +1 defect can be screened by a spherical cap of area $A_{0}/2$. Our energetics calculation is consistent with such topological analysis; it is when the spherical cap becomes as large as a half sphere that a +1 defect will be energetically driven to move to the center of the disk.

We proceed to discuss the split of a +1 defect into two +1/2 defects over a spherical cap. Like the case of the planar disk, we first construct the director field containing two +1/2 defects by cutting an azimuthal +1 defect configuration. As shown in Fig. 5(a), the resulting director field on the spherical cap is composed of three parts: the middle uniform region where $\mathbf{n}=(-z/\sqrt{x^{2}+z^{2}},0,x/\sqrt{x^{2}+z^{2}})$, and the symmetric azimuthal configurations at the two sides. The origin of the Cartesian coordinates is at the center of the sphere, and the z-axis passes through the north pole. The two +1/2 defects are indicated by red dots in Fig.5. Their x-coordinates are $x=\pm c$. The center of the spherical cap is at the north pole. The Frank free energy density of the middle uniform configuration is $f=K_{1}/[2(x^{2}+z^{2})]$. By putting them together and working in the Cartesian coordinates over the equator plane, we have

[TABLE]

where the surface element of the spherical cap $dA=(R/\sqrt{R^{2}-x^{2}-y^{2}})dxdy$, $D_{1}=\{(x,y)|y\in[0,\sqrt{b^{2}-(x+c)^{2}}],x\in[0,b-c]\}$, and $D_{2}=\{(x,y)|y\in[0,\sqrt{b^{2}-x^{2}}],x\in[0,c]\}$. $b$ is the radius of the circular boundary of the spherical cap. $b=r_{s}^{\prime}\sqrt{1-(r_{s}^{\prime}/2R)^{2}}$.

From Eq.(24), we have

[TABLE]

where $F_{1}(x)=R^{2}(c+Rx)/\{\sqrt{-b^{2}+c^{2}+R^{2}+2Rcx}\\ \sqrt{b^{2}-(c+xR)^{2}}[-b^{2}+c(c+2Rx)]\}$, and $F_{2}(y)=R/[(R^{2}-y^{2})\sqrt{R^{2}-c^{2}-y^{2}}]$. It is straightforward to show that $F^{\prime}_{+1/2,s}(0)=-\sqrt{R^{2}-b^{2}}/(Rb)<0$. It indicates the repulsive interaction between two infinitely close +1/2 defects. Numerical evaluation of Eq.(25) shows that when the spherical disk is sufficiently large, the departing +1/2 defects can be stabilized within the disk. The plots of $F^{\prime}_{+1/2,s}(c)$ at typical values for $b/R$ are shown in Figs. 5(b)–5(d). We see that when $b/R>0.93$, the $F^{\prime}_{+1/2,s}(c)$ curve starts to hit the horizontal zero line, leading to the two zero points indicated by the blue and the green dots in Fig.5(d). When the separation between the two defects is smaller than the value at the blue dot or larger than the value at the green dot, they repel with each other. In the regime between the two zero points, the defects attract with each other. The curvature-driven alternating repulsive and attractive regimes in the $F^{\prime}_{+1/2,s}(c)$ curve are indicated by the arrows in Fig.5(d). The left zero point (blue dot) represents the equilibrium configuration of the +1/2 defects. In Fig. 5(e), we show the variation of the equilibrium position of the +1/2 defects with the size of the spherical disk. When the spherical cap occupies more area over the sphere, the distance between the two +1/2 defects in the equilibrium configuration shrinks. In the limit of a half sphere, the two +1/2 defects merge together, becoming a +1 defect. This result is consistent with our previous analysis of the +1 defect case, where the optimal position of the +1 defect over a half sphere is at the center of the disk.

We proceed to discuss nematic order on $Poincar\acute{e}$ disk with constant negative Gaussian curvature Modes and Kamien (2007). The associated metric over a hyperbolic disk with Gaussian curvature $K_{G}$ is characterized by $ds^{2}=4(dx^{2}+dy^{2})/(1+K_{G}r^{2})^{2}$, where $r^{2}=x^{2}+y^{2}$. The area element $dA=4dxdy/(1+K_{G}r^{2})^{2}$. For the director field $\textbf{n}=n_{1}(x,y)\mathbf{e}_{1}+n_{2}(x,y)\mathbf{e}_{2}$, where $\mathbf{e}_{1}$ and $\mathbf{e}_{2}$ are the orthogonal unit basis vectors, its divergence and curl are div $\textbf{n}=(1/2)(1+K_{G}r^{2})(\partial n_{1}/\partial x+\partial n_{2}/\partial y)-K_{G}(n_{1}x+n_{2}y)$, and $\textrm{curl}\ \textbf{n}=(1/2)(1+K_{G}r^{2})(\partial n_{2}/\partial x-\partial n_{1}/\partial y)-K_{G}(n_{2}x-n_{1}y)$, respectively (see Appendix B for the derivation of curl n). We first consider a defect-free uniform director field $(n_{1},n_{2})=(\textrm{cos}\theta_{0},\textrm{sin}\theta_{0})$ whose associated Lagrange multiplier is $\lambda=K_{1}K_{G}(1+K_{G}r^{2})$. $\theta_{0}\in[0,\pi/2]$. The associated Frank free energy density is independent of $\theta_{0}$: $f=K_{1}K_{G}^{2}r^{2}/2$. We see that the uniform state in $Poincar\acute{e}$ disk has a non-zero energy density that increases with $r$ in a power law. It is due to the special metric structure of the $Poincar\acute{e}$ disk. Now we consider a +1 defect configuration in the nematic texture on $Poincar\acute{e}$ disk. It can be characterized by $\textbf{n}=[(c_{2}x-c_{1}y)\mathbf{e}_{1}+(c_{1}x+c_{2}y)\mathbf{e}_{2}]/\sqrt{x^{2}+y^{2}}$, where $c_{1}$ and $c_{2}$ are both constants satisfying $c_{1}^{2}+c_{2}^{2}=1$, such that the magnitude of $\mathbf{n}$ is unity. Varying the value of $c_{1}$ from zero to unity, we obtain director fields from radial to azimuthal configurations. The associated Frank free energy density is

[TABLE]

Since $f^{\prime}_{+1,h}(r)=(K_{1}/2)(-1+K^{2}r^{4})/(2r^{3})<0$, the Frank free energy density decreases with $r$. On the other hand, due to the homogeneity of the $Poincar\acute{e}$ disk, the optimal position of a +1 defect is always at the boundary of the disk.

Finally, we discuss some effects that are not taken into consideration in our calculations. First, by introducing anisotropy in the elastic constants, the free energy varies with the local rotation of the director field. Despite the reduced energy degeneracy arising from the elasticity anisotropy, both radial and azimuthal configurations based on which our calculations are performed are still solutions to the Euler-Lagrange equation (see Appendix C). Therefore, introducing elasticity anisotropy does not change the major conclusions about the optimal positions of both $+1$ and $+1/2$ defects. Second, in addition to curvature, the thickness of liquid-crystal shells is an important parameter to control the number and orientation of defects Lopez-Leon et al. (2011); Koning et al. (2016). It has been experimentally observed that thickness variation can produce a number of novel defect configurations over a spherical liquid-crystal shell Lopez-Leon et al. (2011). It is of great interest to include the effect of thickness in a generalized Frank free energy model to account for these new experimental observations Koning et al. (2016). This is beyond the scope of this study. Third, spatial variations in nematic order parameter within defect cores contribute to the condensation free energy of topological defects Kleman (1983); Chaikin and Lubensky (2000). Notably, nematic textures in defect core regions can exhibit featured patterns and energy profiles, such as highly biaxial nematic order in the cores of $+1/2$ defects Schopohl and Sluckin (1987) and local melting of the nematic ordering Kralj et al. (2011). A recent study has demonstrated that the condensation energy associated with the defect core plays an important role in the formation of defects triggered by strong enough curvature Mesarec et al. (2016). In our study, we focus on the optimal locations of pre-existent defects. They are determined by the variation of the free energy with the positions of the defects, where the contribution from the defect core structures is canceled without considering the boundary effect of defects.

IV CONCLUSION

In summary, we investigate the curvature-driven stability mechanism of LC defects based on the isotropic nematic disk model where the appearance of defects is not topologically required, and present analytical results on the distinct energy landscape of LC defects created by curvature. We show that with the accumulation of curvature effect both +1 and +1/2 defects can be stabilized within spherical disks. Specifically, the equilibrium position of the +1 defect will move abruptly from the boundary to the center of the spherical disk, exhibiting the first-order phase-transition-like behavior. We also find the alternating repulsive and attractive regimes in the energy curve of a pair of +1/2 defects, which leads to an equilibrium defect pair separation. The sensitive response of defects to curvature and the curvature-driven stability mechanism demonstrated in this work may have implications in the control of LC textures with the dimension of curvature.

Appendix A: Infinite degree of degeneracy in the one elastic constant

approximation

Let us consider a planar nematic disk. $\mathbf{n}=(\textrm{cos}\theta(x,y),\textrm{sin}\theta(x,y))$ in Cartesian coordinates. In the one elastic constant approximation, the Frank free energy is $F=\frac{1}{2}K_{1}\iint_{x^{2}+y^{2}\leq r_{p}^{2}}|\triangledown\theta(x,y)|^{2}dxdy$. It is easily seen that rotating a director by a constant angle does not change the Frank free energy.

For a nematic field on a sphere, by inserting $n=(\textrm{cos}\Psi(\theta,\varphi),\textrm{sin}\Psi(\theta,\varphi))$ in spherical coordinates into Eq.(14), we obtain the expression for the Frank free energy density

[TABLE]

where $\nabla\Psi=\frac{1}{R}\frac{\partial\Psi}{\partial\theta}\textbf{e}_{\theta}+\frac{1}{R\textrm{sin}\theta}\frac{\partial\Psi}{\partial\varphi}\textbf{e}_{\varphi}$, and the notation $\Psi_{\varphi}$ is an abbreviation for $\partial\Psi/\partial\varphi$. Obviously, the Frank free energy density is invariant under the transformation $\Psi\rightarrow\Psi+c$.

The conclusion that the nematic texture has infinite degree of degeneracy in the one elastic constant approximation can be generalized to any generally curved surface by writing the Frank free energy under the one constant approximation in the form of

[TABLE]

where the integration is over an area element $dS$ on the surface $\mathbf{x}(u^{1},u^{2})$, $\alpha(u^{1},u^{2})$ is the angle between $\mathbf{n}(u^{1},u^{2})$ and any local reference frame, and $A_{i}$ is the spin connection Vitelli and Nelson (2004). The free energy is invariant under the rotation $\alpha(u^{1},u^{2})\rightarrow\alpha(u^{1},u^{2})+c$.

Appendix B: Calculating curl $\mathbf{n}$ on spherical geometry

In a coordinates independent expression, curl $\mathbf{n}$=$(\star dn^{\flat})^{\sharp}$ Nakahara (2003); Berger (2012); Mikusinski and Taylor (2012). The operators $\star$, $\scriptstyle\flat$ and $\scriptstyle\sharp$ are to be explained below. $\star$ is an operator called Hodge dual. When applied on an antisymmetric tensor $\alpha=\frac{1}{k!}\alpha_{i_{1},\cdot\cdot\cdot,i_{k}}e^{i_{1}}\wedge\cdot\cdot\cdot\wedge e^{i_{k}}$, where $e^{i_{1}},\cdot\cdot\cdot,e^{i_{n}}$ are dual bases,

[TABLE]

$\scriptstyle\flat$ and $\scriptstyle\sharp$ are the musical isomorphisms. $X^{\flat}=g_{ij}X^{i}dx^{j}$, where $X=X^{i}\partial_{i}$. $\omega^{\sharp}=g^{ij}\omega_{i}\partial_{j}$, where $\omega=\omega_{i}dx^{i}$.

Consider a vector field $\mathbf{n}$ defined on a two-dimensional sphere. $\mathbf{n}=n_{1}(\theta,\varphi)\mathbf{e}_{\theta}+n_{2}(\theta,\varphi)\mathbf{e}_{\varphi}$, where $\mathbf{e}_{\theta}$ and $\mathbf{e}_{\varphi}$ are the unit tangent vectors in spherical coordinates. Applying the above formulas on such a vector field, we have

[TABLE]

[TABLE]

and

[TABLE]

We finally obtain Eq.(13). It is of interest to note that the curl of a director field $\mathbf{n}$ on a generally curved surface is $\textrm{curl}\ \mathbf{n}=-\tau_{n}\mathbf{n}-c_{n}\mathbf{t}+\kappa_{n}\boldsymbol{\nu}$, where $\{\mathbf{n},\mathbf{t},\boldsymbol{\nu}\}$ constitute the Darboux basis Napoli and Vergori (2012). $\tau_{n}$ and $c_{n}$ are the components of the extrinsic curvature tensor $\mathbf{L}$. $\mathbf{L}_{nn}=c_{n}$, and $\mathbf{L}_{nt}=\mathbf{L}_{tn}=-\tau_{n}$. In general, the extrinsic curvature influences the Frank free energy of nematics on a curved surface. It is only on a flat or spherical surface $\tau_{n}=0$ and $c_{n}$ is a constant. So the extrinsic curvature effect only contributes a constant term in the Frank free energy Napoli and Vergori (2012).

Appendix C: Euler-Lagrange equations in Cartesian and spherical

coordinates

In this appendix, we present the Euler-Lagrange equations in Cartesian and spherical coordinates derived from Eq.(3), and show that both radial and azimuthal configurations are solutions to the Euler-Lagrange equations. We also show that the anisotropic elastic constants will not change the main result of curvature-driven alternating repulsive and attractive interactions between the two +1/2 defects due to the fact that the elastic modulus $K_{1}$ plays no role in the energy expression.

In two-dimensional Cartesian coordinates, $\mathbf{n}=[n_{1}(x,y),n_{2}(x,y)]$. The components of the director field in equilibrium nematic textures satisfy the following Euler-Lagrange equations:

[TABLE]

and

[TABLE]

It is found that both radial $(n_{1},n_{2})=(x/\sqrt{x^{2}+y^{2}},y/\sqrt{x^{2}+y^{2}})$ and azimuthal $(n_{1},n_{2})=(-y/\sqrt{x^{2}+y^{2}},x/\sqrt{x^{2}+y^{2}})$ configurations satisfy the above Euler-Lagrange equations with $\lambda=K_{1}/(x^{2}+y^{2})$ and $\lambda=K_{3}/(x^{2}+y^{2})$, respectively. The spiral configuration $\mathbf{n}=c_{1}(x/\sqrt{x^{2}+y^{2}},y/\sqrt{x^{2}+y^{2}})+c_{2}(-y/\sqrt{x^{2}+y^{2}},x/\sqrt{x^{2}+y^{2}})$ ($c_{1}^{2}+c_{2}^{2}=1$ and neither $c_{1}$ nor $c_{2}$ is 0) is the solution to the Euler-Lagrange equations only in the one elastic constant approximation.

In spherical coordinates, $\mathbf{n}=n_{1}\mathbf{e}_{\theta}+n_{2}\mathbf{e}_{\varphi}$. In equilibrium nematic textures, $n_{1}$ and $n_{2}$ satisfy the following Euler-Lagrange equations:

[TABLE]

We remark that the equilibrium equations in spherical coordinates are invariant under uniform local rotation of the director field. Similarly, one can show that both radial and azimuthal configurations are solutions to the Euler-Lagrange equations with $\lambda=K_{1}/(R^{2}\textrm{sin}^{2}\theta)$ and $\lambda=K_{3}/(R^{2}\textrm{sin}^{2}\theta)$, respectively. The spiral configuration of $n_{1}=c_{1},n_{2}=c_{2}$ ($c_{1}$ and $c_{2}$ are non-zero constants satisfying $c_{1}^{2}+c_{2}^{2}=1$) satisfies the equilibrium equation only when $K_{1}=K_{3}$. For the two $+1/2$ defects configurations discussed in the main text, we show that introducing elasticity anisotropy does not change the curvature-driven alternating repulsive and attractive interactions between the defects. For the two $+1/2$ defects configuration on a spherical disk where an azimuthal configuration is separated by a uniform configuration, the associated Frank free energy is

[TABLE]

where $D_{1}$ and $D_{2}$ are given below Eq.(24). We see that since the entire defect configuration is divergence free, the parameter $K_{3}$ does not appear in the expression for the Frank free energy. Therefore, anisotropy in elastic constants does not change the featured interaction between the $+1/2$ defects.

Acknowledgement

This work was supported by NSFC Grant No. 16Z103010253, the SJTU startup fund under Grant No. WF220441904, and the award of the Chinese Thousand Talents Program for Distinguished Young Scholars under Grant No. 16Z127060004.