Abstract
We studied the impact of the cell thickness on configurations of line disclinations within a plane-parallel nematic cell. The Lebwohl-Lasher semimicroscopic approach was used and (meta)stable nematic configurations were calculated using Brownian molecular dynamics. Defect patterns were enforced topologically via boundary conditions. We imposed periodic circular nematic surface fields at each confining surface. The resulting structures exhibit line defects which either connect the facing plates or remain confined within the layers near confining plates. The first structure is stable in relatively thin cells and the latter one in thick cells. We focused on structures at the threshold regime where both structures compete. We demonstrated that “history” of samples could have strong impact on resulting nematic configurations.
1. Introduction
Line topological defects are ubiquitous in nematic liquid crystals (NLC) which is fingerprinted even in their name [1]. They could be stabilised topologically by appropriate boundary conditions [2, 3] or due to energy reasons [4]. They have a strong impact on optical NLC properties and are therefore of potential interest for various electro-optic applications.
Nematic uniaxial liquid crystals exhibit simultaneously liquid properties and local orientational order [5]. The latter is at macroscopic level commonly presented by the nematic director field . It points along the local uniaxial order, where states with are physically equivalent. In bulk equilibrium is spatially homogeneous and aligned along a single symmetry breaking direction. NLC can exhibit line dislocations, which are characterised by the winding number m = ±1/2 [6, 7]. It reveals the total reorientation of on encircling the defect in counter-clockwise direction. Furthermore, one can assign the total topological charge q to a line defect [6] by enclosing it by a surface. This charge is integer and reveals how many realizations of all orientations are realised in the nematic order parameter space [6]. Note that the core structure of line m = ±1/2 line defects is biaxial [8], and the center of the core exhibits negative uniaxiality.
Line defects could be stabilised in different ways. For example, they could be enforced by AFM imprinted patterns to plates enclosing NLC in plane-parallel geometry as illustrated in [2]. It has been shown [9] that, in such geometries the line defects could either span the opposite plates or are confined to the vicinity of the bounding plates. In this paper we focus on these competing structures. We henceforth refer to the defect configurations that (i) connect the bounding plates and (ii) remain confined close to the planes, as the (i) connected and (ii) confined defect configurations, respectively. In our study we consider networks of line defects in a plane-parallel cell of thickness h. We use the Lebwohl-Lasher semimicroscopic lattice model [10, 11] where the local orientational order is presented by nematic pseudospins. We assume that the bounding plates are patterned by a lattice of concentric circles, enforcing the circular planar nematic alignment, which give rise to line defects. We focus on the impact of h and history of samples on defect patterns. The plan of the paper is as follows. In Section 2 we present model and in Section 3 results. In the final section we summarize our results.
2. Materials and Methods
We consider nematic structures within a plane-parallel cell of thickness h. At the bounding plates we enforce spatially varying nematic patterns and calculate the corresponding nematic structures within the cell. In our modelling we use the semimicroscopic Lebwohl-Lasher lattice approach [10, 11]. In this modelling the nematic orientational ordering is described in terms of nematic spins (states with are physically equivalent) residing at each lattice site. We henceforth refer to unit vectors as spins. The simulation lattice is cubic, characterised by the lattice constant a0.
The simulation lattice consists of M × M × L sites in the Cartesian coordinate system (x, y, z). Here L is proportional to the cell thickness, i.e., h = La0, while M L stands for its lateral dimensions. In the following we set a0 = 1, so that we identify L with the cell thickness. Individual sites are denoted by a set of indices (i, j, k): 1 i ≤ M, 1 j ≤ N, and 1 k L. Each site is occupied by a spin , which tends to orient in parallel direction with its nearest neighbors. The bounding plates have imprinted a two-dimensional (2D) circular nematic “surface field” (see Figures 1 and 2, top left), which enforces line defects to the LC body. Therefore, we impose a lattice of two-dimensional m 1 topological point defects. This surface imprinted structure is positioned symmetrically on both plates. In the model the surface field is determined by frozen-in spins at the bounding plates.
The coupling between neighboring spins inside the cell is given by the coupling constant J 0, while the coupling between the plates and the nearby spins is parametrized with the anchoring strength W 0. The bulk interaction energy F is the sum of the terms over all lattice sites:where The subscript n.n. in the inner sum denotes the 6 nearest neighbors of the site with indices (i, j, k). The bottom and upper boundary plates are determined artificially by k = 0 and k = L + 1, respectively, and impose the nematic structure determined by frozen-in spins. The corresponding nearest inner spin planes (with k = 1 and k = L, respectively) interact with the boundary planes with the LL interaction written in (1), where . At the lateral boundaries we impose the periodic boundary condition. Besides the torques due to coupling between neighboring spins the random thermal fluctuations according to the Gaussian distribution of deviations from equilibrium are included. For numerical purpose, we introduce the dimensionless time step = D∆t and the dimensionless temperature = kT/J = kT. Here D is the effective rotational diffusion coefficient, ∆t is the physical time step, and k is the Boltzmann constant. Simulation details are described in [11]. Furthermore, we scale W in units of J.
A surface imposes spin configuration, consisting of a lattice of 2D circular defects with the winding number m = 1 on the boundary plate. Each circular pattern, centered at (x0, y0), is determined byIn (2a) and (2b) φ is the azimuthal angle in the (x, y) plane and θ is the angle between the spin and the x-axis; therefore θ = φ + π/2 imposes a circular pattern. We place four equal and symmetrically positioned circular defects on each plate as shown in Figure 1. In between the circles we impose a homogeneously aligned spin field along the x-axis. Both confining plates impose the same spin field, i.e., .
In simulations we vary the cell thickness L. The remaining parameters are fixed: W = 5 and = 0.5. These conditions correspond to a relatively strong anchoring deeply in the nematic phase.
For given conditions we calculate spin configurations using Brownian molecular dynamics, where details are described in [11]. We present key structural characteristics using the following quantities. We define the traceless tensor order parameter:Here Sijk,m denotes the m-th component of the spin , while the triangular bracket denotes the average over the cell and δmn is the Kronecker symbol. The scalar order parameter S is conventionally defined as the largest eigenvalue of the tensor.
In order to measure impact of the bounding plates on the average LC configuration within the cell we define the quantity S(0, k) which measures the average matching of spins in the k-th plane inside the cell with the frozen spins in the confining boundaries, i.e.,If the matching is perfect it holds S(0, k) = 1. Another useful macroscopic parameter is the mean value of the square of the spin component Sz, again with averaging over horizontal planes:Note that the bounding plates enforce zero z-component of the spins.
We study also the impact of sample “history” on spins. For this purpose we impose four different initial spin configurations: spin configuration imposed by the plates (which we label by “D”), homogeneous configuration with all the spins along the x-axis (“XH”), homogeneous configuration with all the spins along the z-axis (“ZH”), and random configurations, where the directions of the spins are randomly distributed according to isotropic distribution (“R”).
3. Results and Discussion
Of our interest is to analyse the thickness h=La0 driven transformation between connected and confined defect configurations. In thin enough cells the surface boundary conditions stabilise a connected defect configuration. When the cell thickness is increased, the defect structure, as enforced by the plates, tends to unwind into preferable homogeneous structure in the bulk of the cell. This is indicated by 2D cross-sections parallel to the plate in Figures 1 and 2, where the initial spin configuration was “D” and the equilibrium spin configurations for two different thicknesses of the cell are compared: L = 10 (very thin cell) and L = 30 (thicker cell), respectively. Because of the symmetry it is sufficient to illustrate the spin configurations for the lower half of the cell, i.e., k ≤ L/2. Thus the four configurations presented in either of Figures 1 and 2 represent the spin configuration for the planes with increasing distance from the bottom plate, indicated by increasing parameter k. In order to have a definite (albeit superficial) direction of the spin also in the center points of the defects, the spin was set to be aligned in z-direction, i.e., perpendicular to the plates there (dots at the defect centers in the figure).
It is evident from Figure 1 (L = 10) that the defect structure is essentially preserved through the cell for small thickness, but the centers of the defects are smoothed. Thus the line disclinations essentially connect the confining plates. On the contrary, in Figure 2 (L = 30) the spin configuration in the middle of the cell is roughly homogeneous along the x-axis.
Different structures for cells with different thicknesses can also be forecast theoretically and revealed experimentally with the use of polarizing optical microscopy [12, 13]. The essence of the optical polarizing microscopy is the following. The NLC is positioned between a pair of crossed polarizers (perpendicular polarization directions). In the case of optically inactive medium in the cell no light passes through the second polarizer. However, the optically active medium, such as NLC, rotates the polarization direction because of the difference in ordinary and extra-ordinary refraction indices, and so the intensity of the transmitted light is nonzero. The degree of polarization rotation depends on nematic director field. In numerical simulation of the optical pattern, which mimics the experimental situation, the path of the beam through the NLC cell is divided into short segments, and at each segment the small rotation angle of the polarization is calculated by using the Jones matrix. If the direction of the light beam is z (perpendicular to the cell plates), the intensity of the transmitted light generally depends on (x, y) coordinates, resulting in characteristic optical pattern. Each pattern in our calculation corresponds to a single configuration (snapshot). Theoretically modelled patterns for L = 10 and L = 30 are compared in Figure 3. The left image (L = 10) is characteristic for the winding number m = 1. However, some traces of the optical patterns of the m = 1 line defect are still present for thicker cell (L = 30, right image).
The structures (equilibrium spin configurations) presented in Figures 1 and 2 above indicate that there exists a threshold value of the cell thickness above which the defect structure unwinds in the middle part of the cell (thicker cells). In our study is between the values 20 and 30, where the thickness of the cell is approximately equal to the size of the imposed circles. The direction of the spins in the predominantly homogeneous part of the cell is parallel to the x-axis because this direction is favored by the alignment of frozen spins between the circles.
Furthermore, we have studied the effect of different initial spin configurations. In the case of thin cells (L = 10 or slightly larger) the final equilibrium configurations for all different initial configurations (“D”, “XH”, “ZH”, and “R”) are similar. This was shown by comparison of 2D horizontal cross-sections after relaxation, similarly as for Figures 1 and 2, and it has been also confirmed by similar values of parameters S(0, k) and in (4) and (5). The value of is of the order of 0.02 in the middle of the cell (k = L/2), indicating that the spins predominantly lie in the (x, y) plane. The small Sz component is present because of thermal fluctuations.
Only the “ZH” initial spin configurations result in the equilibrium configuration which tends to deviate from the (x, y) plane configuration to larger extent, at least for larger values of L. For instance, for L = 30 we have found that the “ZH” initial configurations lead to equilibrium configuration with ≈ 0.6 at the middle of the cell, while the other three initial configurations result in ≈ 0.02, as in the case L = 10. However, other differences for initial spin configurations “D”, “XH”, and “R” have been noticed. As expected, the final configuration in the cell has the largest similarity to the defect configuration on the plate, when the initial configuration is “D”. This is also confirmed quantitatively by the largest value of the parameter S(0, k) in (4) for all 2D planes (all values of k). Figures 4 and 5 show the plane averages and S(0, k) as functions of the index k for all initial spin configurations and for thicknesses L = 10 and 30, respectively.
The degree of alignment of spins along the x-axis could be quantified by in a similar way as was defined in (5) for the z-component. The “R” initial configuration results in the equilibrium spin configuration between those corresponding to “D” and “XH”. We have also noticed that for thin cells (L = 10) the quantitative parameters, such as S(0, k), for “R” initial spin configuration are closer to the corresponding values for “D” as compared to “XH”. The opposite is true for thicker cells. Our simulations reveal that for initial conditions “D”, “XH”, and “R” the structural differences decrease on increasing the cell thickness. However, the “HZ” initial spin configuration could lead to significantly different structures as Figure 5 reveals. In this case we also observe a pronounced statistical behavior of the system with respect to the “escape” of spin configurations along the z-direction. By repeating calculations several times for the same set of parameters we observe that the system either “escapes” or does not and the extent of this bifurcation behavior increases on approaching the threshold thickness . The reason behind this is thermal fluctuations included in the simulations. As a quantitative measure of the bifurcation extent we set that the “escape” is realized when (k = L/2) exceeds = 0.1. Our statistical analysis in the cell 62 × 62 × 30 reveals that the probability for the escape along the z-direction is roughly equal to 40 %.
We also calculated the interaction energy per site of the equilibrium states = /(M2L) (see (1)). We have found that it increases slightly with the cell thickness: , , , and for L = 10, 15, 20, and 30, respectively. The uncertainty of these values due to thermal fluctuations in our simulations is δ ~ 0.01. The quoted energy values must be compared to the energy per site for completely homogeneous state: since J = 1. The smaller energy per site for very thin cells results from relatively strong coupling between the surface and neighbouring spins: W = 5. The surface anchoring energy has been included into the total energy and its contribution to the sites at both plates compensates the rise of energy due to inhomogeneity. On the other hand differences in the equilibrium energy for the four different initial configurations for the same L are roughly comparable to fluctuation δ.
4. Conclusions
We investigated nematic structures within plane-parallel cells of thickness h=La0 whose bounding plates enforced topologically line-defect patterns. The circular nematic fields at each bounding surface favour surface topological defects bearing m = 1 winding number at the centers of the circular patterns. Due to the finite anchoring strength these nucleation sites split and as a result two line defects with m = 1/2 emerge from each nucleation site. In order to achieve winding number neutral structure two m = are formed to compensate the total winding number within each (x,y) plane. In thin cells these defect lines span the bounding plates and their winding number remains the same across the cells. On the other hand in thick cells the defect lines remain confined to the vicinity of the bounding plates. Consequently, they switch the sign of the winding number while linking m = 1/2 and daughter nucleation sites. We analysed structures in the critical region, L ~ , where these two different realizations of defect patterns compete. The resulting patterns depend relatively strongly on initial conditions. We imposed four different histories to the nematic cells. We either imposed the surface imposed pattern through whole the cell, homogeneously aligned pattern along the x-axis or z-axis, and isotropic profiles. We label these histories as “D”, “XH”, “ZH”, and “R”, respectively. These patterns could be experimentally imposed. For example, the first pattern is realised in NLC cells below the critical thickness. The second and the third patterns could be realised by using a strong enough spatially uniform electrical field aligned along the x-axis and z-axis, respectively. The fourth profile could be realised by quenching the sample from the isotropic phase. In the case of “ZH” we find behaviour which is reminiscent to the percolation one. In our future study we intend to analyse external electric field driven rewiring of defect patterns in thin (L ) and thick (L ) cells.
In the studied case we considered substrates imposing four circular “easy axis surface fields”. Note that larger number of circles and different distance between their centers or different arrangement (e.g., hexagonal) would quantitatively affect the results. However, we claim that qualitative structure would be similar. Namely, the key parameter influencing the switching mechanism is the length of disclinations. The crucial parameters are therefore the cell thickness h and the shortest distance between termination ends of “confined” disclinations.
Furthermore, we considered conditions for which each m = 1 surface imposed singularity nucleated a pair of m = 1/2 line defects. However, in general also other structural solutions might appear depending on h and the strength of the surface anchoring interactions. Namely, in our modelling we assume that the confined surfaces locally favour unsplit m = 1 defect structure at the defect nucleating sites. On the other hand, bulk elasticity favours formation of defect lines bearing “elementary” (smallest possible) winding numbers, i.e., . A simple modelling [5] in the approximation of equal nematic elastic constants reveals that the director field dependent elastic free energy cost of a single disclination of strength m within a finite region (characterised by a linear length R) scales as , where Consequently, in bulk a pair of m = 1/2 dislocations is energetically more favourable than a single m = 1 line defects. Our preliminary simulations (using Landau-de Gennes continuum approach in terms of the tensor nematic order parameter), where we focus on structures emerging from a single m = 1 surface imposed singularity, reveal three qualitatively different nematic structures on varying h. For a single m = 1 singularity is most stable. Within the window a pair of two m = 1/2 disclinations appears instead. Furthermore, for an “escaped” structure corresponds to the global minimum, where the singularity in the director field is avoided by a nonsingular escape of the nematic director field along the direction [14]. Note that the described behaviour on crossing was already observed experimentally [15, 16]. A more detailed analysis will be published elsewhere.
The observed behavior is reminiscent to the percolation phenomenon. We illustrate this analogy for the simplest case of cubic cell (network) in percolation theory. The cell is divided into small cubes (sites) with two significantly distinct values of some physical quantity. Let us take the usual case of either electrically conducting or insulating sites. The conducting and insulating sites are randomly spatially distributed, so that there is a fraction p of conducting sites. In the vicinity of the percolation threshold and for a finite network there is a probability P that a randomly distributed conducting sites form a conducting (percolation) path between two opposite faces of the entire cube. The percolation threshold is usually defined as the fraction of conducting sites (i.e., the volume fraction of conducting phase) where the percolation probability has the value P = 50 %. The percolation threshold is roughly independent on the network size, but the steepness of the monotonically increasing function P(p) around increases with increasing network size [17]. It holds P = 0 for p , while P = 100 % for p >> . The fraction p in percolation theory plays similar role as the thickness of the NLC cell in our study, and corresponds to the threshold thickness .
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.