Controlling motile disclinations in a thick nematogenic material with an electric field

Manipulating topological disclination networks that arise in a symmetry-breaking phase transformation in widely varied systems including anisotropic materials can potentially lead to the design of novel materials like conductive microwires, self-assembled resonators, and active anisotropic matter. However, progress in this direction is hindered by a lack of control of the kinetics and microstructure due to inherent complexity arising from competing energy and topology. We have studied thermal and electrokinetic effects on disclinations in a three-dimensional nonabsorbing nematic material with a positive and negative sign of the dielectric anisotropy. The electric flux lines are highly nonuniform in uniaxial media after an electric field below the Fréedericksz threshold is switched on, and the kinetics of the disclination lines is slowed down. In biaxial media, depending on the sign of the dielectric anisotropy, apart from the slowing down of the disclination kinetics, a nonuniform electric field filters out disclinations of different topology by inducing a kinetic asymmetry. These results enhance the current understanding of forced disclination networks and establish the presented method, which we call fluctuating electronematics, as a potentially useful tool for designing materials with novel properties in silico.

Topological singularities such as points, lines and walls are ubiquitous in phases with broken symmetry. Canonical examples include dislocations in solids 1 , vortex lines and rings 2 in superfluid 3 He and 4 He, Abrikosov vortex lines in superconductors 3 , vortex lines in Bose-Einstein condensates 4 , umbilic lines 5 and π solitons 6 (disclinations) in nematic liquid crystals (NLC) that provide a testing ground for theories of cosmology 2,7 , λ lines in cholesteric fluids 8 , Bloch and Néel lines in ferromagnets 9 , walls in lipid membranes 10 and string networks in ecology 11 . NLC phases display rich birefringence under a polarizing microscope during phase ordering from a disordered state after a rapid quench in pressure or temperature, resulting in the formation of disclinations with integer and fractional topological charge. These singularities proliferate after nucleation and form contractile loops after intercommutation 12 . Unlike dislocations, disclinations possess intricate kinetics, microstructure, and equivalence with an electric charge. Strings are charge neutral with either topological charge ±1 or ±1/2 residing at the two segments or end points of the string to form topological dipoles. Higher multipoles and integer charged dipoles also nucleate within the charge neutral strings at the early stage of kinetics. Subsequently, these structures rupture into fractionally charged dipolar strings. Similar to electrodynamics, like topological charges repel and unlike charges attract and annihilate in pairs while monopoles are nonexistent to retain charge neutrality unless created by symmetry-breaking boundaries, an inclusion of impurity or external drive with a laser beam 13 .
Existence, classification and recombination rules of disclinations in equilibrium, which play an important role in the material design, is governed by the energy landscape as well as the geometry (topology) of the order parameter space 12 . Strings in uniaxial NLC displayed in Fig. 1 (frames a-d) are topologically defined by RP  π = ( ) 1 2 2 with homotopy group π 1 in the projective plane 2 RP resulting in the abelian group  2 with topological charge ±1/2 12 . After a theoretical proposal 14 , π solitons have been seen in fluorescence confocal polarized-light microscopy of pentylcyanobiphenyl (5CB) NLC 6,7 , molecular simulations 15 and field theoretic computations 13,16 . Likewise, biaxial disclinations displayed in Fig. 1 (frames e-h) are defined by (g) Uniaxial order and director distribution on a portion of the xy-slice plane of (e). Note the similarity in the microstructure of ±1/2 C z defects with uniaxial defects in (c) and no variation in n for ±1/2 C y defects. (h) The corresponding variation of uniaxial and biaxial order displaying core structure of disclination segment of both class in-plane. Although both scalar orders decrease in the core of C y defect, biaxiality increases for a decreasing uniaxiality in C z defect. Parameters are defined in Methods (section 3) and material (computation) parameters are tabulated in Table 1.  20,26 . By investigating beyond the uniform field assumption 27 , in this article, we have developed a fluctuating electronematics method based on the thermal description of the GLdG theory, with the physical control over the role of each forcing to accurately describe thermal and electrokinetic phenomena in three dimensional NLC media. Using simple analytical argument, we provide an understanding of the underlying mechanism responsible for the outcome in both uniaxial and biaxial media in the free draining limit at moderate to small electric field intensity, where the advective flow of the anisotropic media can be neglected. The interdependency between the topology in the orientational order of the media and morphology of the external field is elucidated through the measurement of disclination kinetics and morphology. It turns out that small magnitude of nonuniform electric field can significantly dilate the coarsening kinetics of disclination network and can eradicate certain topological class of disclinations in biaxial media.

Results
Thermoelectrokinetic effects of disclination network in a coarsening uniaxial NLC. Here we systematically probe on the role played by the elasticity of the medium and various external forcing (k B T, E) on the kinetics and microstructure of the string disclination assembly.
Role of thermal fluctuations and anisotropic elasticity. The disclination kinetics is immensely influenced by external agents like thermal and electric forces, anisotropic elastic effects or shear. To illustrate the role of thermal fluctuations, we estimate the decay of disclination density per unit area in the degenerate elastic constant approximation without an electric field (κ = Θ = |E| = 0 in equation (7)). Using surface triangulation method 28 , we calculate the total surface area of the disclinations in the media. Figure 2: frame a plots the same for three different values of k B T including the athermal scenario and in frame (b) we show a portion of the three-dimensional volume that distinguishes the disclination kinetics between the athermal and thermal scenario. In the athermal scenario, three different regime with marked exponents emerges out in the evolution process, that had been previously quantified as diffusive regime → Porod's law regime → diffusive regime 16 . The early diffusive regime corresponds to domain coarsening before nucleation of disclinations (t = 0.84 ms), while Porod's law scale designates the defect annihilation kinetics (t = 1.4 ms) and finally, the late stage diffusion is attributed to contraction of isolated loops (t = 11 ms). Clearly, thermal fluctuations tend to increase the disclination surface density without affecting the scaling laws, that is important for materials (e.g. PAA) having transition temperature way above the room temperature. Fluid viscosity η can be obtained from the Stokes-Einstein relation k B T/Kη = constant. However, frame c shows that disclination density in the athermal media not only increases for an increase in κ, but also stretches the Porod's law regime. The slope changes from 0.8 to 1 as elastic anisotropy is increased. While the slope of unity is also obtained in experiments with 5CB 7,29 , we re-establish our earlier claim 20 about the crucial contribution of the anisotropic elasticity of the medium. Intuitively, anisotropic elastic constant results into asymmetric diffusion constants in Cartesian directions and thus brings asymmetry in the speed of ±1/2 integer point defects 23,30 . The loss of area in forming a contractile loop is greatly reduced for higher anisotropy.
Before we embark on discussing coarsening in the presence of an electric field, we revisit the coarsening kinetics when no electric field is present. As seen in Supplementary Movie S1, disclinations with higher line tension and curvature are energetically disfavoured, culminating in stretched strings that form contractile loops after intercommuting with neighbouring strings. By taking into account the effect of viscous drag and length decrement of a string in forming a loop at time t, the disclination surface density ρ is found 7 to scale as ρ ∝ t −1 . An identical scaling law is obtained for a planar disclination when equating the rate of change of the line tension per unit volume with the energy density loss rate 1 . Our method accurately reproduces the slope of 1 ± 10 −3 within the   Table 1. Parameters values excercised to mimic uniaxial (upper row) and biaxial (lower row) thermotropic NLC. A rectangular simulation box of size 80 2 × 160 μm 3 with grid spacing Δx = Δy = Δz = 1 μm and time step Δt = 1μs is considered. We use material parameters for 5CB at T = 33.65 °C, ε 0 = 1,ε a = ±1, ε s = 0.74 ε a 51 and use earlier excercised material parameters for biaxial media 16,32 . Using equation (4), we estimate E F = 2 × 10 −3 V/μm for 5CB and with |E| = E F × (10 −2 , 10 −1 ,1), nondimensional parameters are ε 1 = 1.5 × (10 −4 , 10 −3 , 10 −2 ), ε 2 = 1.8 × (10 −3 , 10 −2 , 10 −1 ). In biaxial media, we use |E| = 1.5 V/μm < E F , ε 1 = 0.14, ε 2 = 3.35.  Fig. 3 ('no field' curve in frame g) with the material parameters of 5CB. Thermal fluctuations and elastic anisotropy lead to an increment in the disclination surface density at a given time and the period of disclination annihilation kinetics is extended without affecting the physical laws.
Role of an electric field. Next, we elaborate on the thermal and electrokinetic effects on coarsening in uniaxial NLC after the onset and cessation of an unidirectional voltage pulse. The electric field is applied at an instant t on when only disclinations with dipolar charge ±1/2 are present, and, the medium is free of ±1 dipolar strings or point charges. The switch-off time of the field is set at an instant t off when the electric flux lines within the medium do not change substantially. For a shallow quench below the supercooling line, for both signs of the dielectric anisotropy constant and below the Fréedericksz threshold of the electric field, the orientational order is comparable in magnitude with the non-dimensionalized electric potential. Thus within the medium, the flux lines are very much distorted as seen in Fig. 3. Instead of being aligned along the field direction, the isotropic cores of the strings displaying reduced uniaxial order are little deformed by the electric force 6 and have a little contribution in distorting the flux lines. This is attributed to the weak coupling of scalar order to the electric field, unlike the director that strongly couples to the electric field. In no switch-off scenario (t off → ∞), the nonuniformity of the flux lines is retained in the electrically forced nematic phase devoid of disclinations around t∼35 ms (not shown). The flux lines in Fig. 3 (frames a-f), however, tend to gain uniformity as the electric field is increased towards the Fréedericksz threshold for both signs of the dielectric anisotropy constant. As observed in frame g, electrically forced strings are little thinner due to a reduction of the surface density, but they are long lived due to a dilated kinetics. After the field is switched off, disclination kinetics and surface density do not immediately return to the zero-field behavior, but lag for an interval of ∼ ms 10 . Thus, the material retains a memory of the field onset in the process of exhibiting an elastic response. The electric field agitates the isotropic background towards an uniaxial medium, thus increasing the uniaxial order (frame h) and non-monotonically decreasing the total free energy (frame i). Though we obtain Brochard-Legér lines connecting π solitons under an intense field above the Fréedericksz threshold, the lack of backflow in our method cannot reproduce a ceasing motility of disclinations 6 . Rather, they annihilate at a faster rate due to the uncompensated electric drag force. Electrokinetic effects under intense forcing, such as electroconvection or rheochaos, can be quantitatively captured only if backflow is systematically included. Unlike in colloidal suspensions 24 , the question of the existence of correlations between fluctuations in orientation and velocity has to be answered from experiments 31 before attempting a numerical study in which backflow effects are included.
We qualitatively argue about the slowing down of disclination kinetics for a planar isolated disclination loop. An in-plane estimate is valid for 5CB while twist constant is much smaller than splay or bend constants 20 . Also deep within the uniaxial phase where the director is aligned uniformly, Frank constant K can sufficiently define the medium elasticity 26 . In such situations, planar disclination energy per unit length is  = being the topological charge,  the bounding plane and c a constant, defines a planar defect configuration n = [cosf, sinf]. Performing the surface integral, the reduced disclination energy per unit length is A portion of the three-dimensional volume is shown for clarity, that displays domain decomposition before nucleation of disclinations → intercommuting disclinations → contractile disclination loops. (c) Increased elastic anisotropy leads to increased disclination density with a prolonged Porod's law regime. Note that both thermal fluctuations and elastic anisotropy aids in early nucleation of the isotropic domain, but prolongs the disclination kinetics. The rendered colours in field values are indicated in Fig. 1  where  discl c is the disclination core energy. Equating the elastic energy ξ / discl  per unit area with the drag force −η∂ t ξ per unit length yields contributions from (i) the equilibrium kinetics ξ ∼ − t 1/2 and (ii) the electric field effect ξ ∼ ν e t with ν = πε 0 ε a E 2 /8kη. While ν is independent of the sign of k or ε a , ν > 0 implies of a temporal reduction of the loop extinction kinetics. Physically this can be interpreted as a reduction in speed of approach between ±1/2-charged topological dipole within a charge-neutral loop due to the external forcing. As observed in Fig. 3 (frame g) and in supplementary Figure S1, nonuniform electric field substantially prolong the kinetics when compared to the uniform field scenario (see Supplementary Movie S4).
To shed light on the effect of electric forces on the disclination core structure in thermal uniaxial media with ε a > 0, in Fig. 4 we sketch the spatial variation of the surface of S around a planar defect for different field intensity and compare with the equilibrium scenario. The fluctuation amplitude at S ueq is reduced due to the application of an electric field, resulting in a reduction of the disclination surface density (Fig. 3, frame g). However, we do not find any significant distortion of the core for an increasing field strength which indicates that the field, unlike the director orientation, cannot influence the sufficiently isotropic core other than a complete melting of the disclina- Thermoelectrokinetic effects in coarsening biaxial NLC. Next, we examine the role played by the isotropic elasticity of the medium and various external forcing (k B T, E) on the kinetics and microstructure of the string disclination assembly of different homotopy class. We do not find excitingly different outcome when investigating the role of anisotropic elasticity and, thus, here we restrict ourselves in reporting results in κ = 0 limit, in par with other investigations 16,32 .  Table 1.
Role of thermal fluctuations. We estimate the consequence of thermal fluctuations on the biaxial disclinations of class {C y , C z }. In Fig. 3, we plot the evolution of disclination density for different values of k B T and compare with the athermal scenario. To remind, in the inset to the left panel, we portray the early, intermediate and late  stage of the disclination kinetics. Similar to the uniaxial disclinations, here we also find that thermal fluctuation increases the disclination density per unit area with a comparable slope. The increase in slope for an increase in k B T during early stage of the kinetics for C y class hints for a delayed emergence of Porod's regime, that is absent in C z class. More interestingly, we observe an equivalence of the C z class of biaxial disclinations with that of the ±1/2-integer disclinations in uniaxial nematics, both (i) in morphology ( Fig. 1: frame c and frame g) and (ii) kinetics, as observed in the identical slope of 0.8 in the Porod's law scaling regime ( Fig. 5: right panel and Fig. 2: frame a). Also, the mismatch of the slope between the C y and C z class of disclinations (Fig. 5) suggests of minor influence within each other in the course of annihilation.
Role of an electric field. To conclude, we examine the thermal and electrokinetic effects on the coarsening of biaxial NLC when the sample is rapidly cooled from a disordered phase in the presence of a steady voltage pulse. Figure 6 plots the instantaneous snapshots as well as the kinetic evolution of the medium. Although the topological structure and kinetic pathway of disclinations in biaxial NLC have been predicted for long 1,18,33 , experimental advance to stabilize disclinations by avoiding crystallization continues to be the holy grail of research on thermotropic biaxial mesogens 34 . Similar to the dilated kinetics of electrically forced uniaxial disclinations (see Fig. 3), we find in Fig. 6 (frame g) that the onset of an electric field increases the lifetime of biaxial disclinations of homotopy class {C y , C z } at the initial stage for both sign of dielectric anisotropy. For either class of disclinations, ν(= πε 0 ε a E 2 /8kη) >0 qualitatively explains the slowing down. For both signs of the dielectric anisotropy constant of the material, flux lines are massively distorted in the presence of disclinations (frames a,d). After an interval of ~5 ms, a clear asymmetry between the disclination kinetics of different topology becomes evident. As shown in frames b,c,e,f,g, this results in long-lived disclinations of either class with uniform electric field lines. This is attributed to an increment (decrement) of the total free energy with a positive (negative) value of the dielectric anisotropy constant (frame h). Thus, the dielectric energy has a strong influence in selecting disclinations of the desired class as the bulk and elastic energies increase negligibly from the no field scenario. Physically, the acceleration (or retardation) in the loop extinction kinetics at the late stage can be interpreted as effective acceleration (or retardation) in speed of approach between ±1/2-charged topological dipole within a charge-neutral loop of

Discussions
Electrorheology of line defects, with [35][36][37] or without 23 particulate inclusion, under an intense electric field have established the coupling of orientation tensor with hydrodynamics and uniformity in the electric field, though the effect of nonuniformity in the electric field 19,38-40 and the effect of thermal fluctuations are less explored. The effect of hydrodynamics is assumed to be negligible under moderate to weak electric field intensity 41 . We have examined the role of thermal fluctuations and nonuniformity of electric field in this limit and have shown that fluctuating electronematics is a robust tool to mimic laboratory experiments 6,7,13 in silico for anisotropic NLC in three dimensions. From the structure of the orientation tensor, we present a simple way to identify and classify the line defects and to compute physical quantities from the geometry of disclinations. We have shown how the spatial uniformity in electric field is gained in approaching the Fréedericksz limit. Apart from modifying the kinetic pathway of the coarsening of athermal disclination network, the external stimuli in terms of temperature fluctuations and local electric field essentially probe two emergent length scales: (i) interfacial correlation length between isotropic and nematic phase and (ii) radius of curvature of disclination loop. Neglecting any local heating effects due to the variation of temperature, any change in the correlation length is attributed to the disclination core size as well its geometric position within the three-dimensional volume. Although the evolution is temporally dilated compared to the no-field scenario, the external field cannot sufficiently modify the radius of curvature of isotropic disclinations -thus the lines are not stretched along the direction of the electric field, rather they retain their shape even when the director gets aligned along or perpendicular to the field direction depending on the sign of material's dielectric anisotropy. The inhomogeneity of the nematic orientation is manifest in the inherent nonuniformity of the local electric field -resulting in the highly nonuniform electric flux lines within the sample. The electric field induces a memory to the material that exhibit an elastic response and also induces a kinetic asymmetry within disclinations of the different class. On the other hand, increase in thermal fluctuations tends to increase the disclination surface density.
This complex interaction can be intelligently engineered to yield a fascinating outcome in a more complex scenario, for example, fractal nematic colloids 42 , metadevices 43 and photonic applications 44 . The electric field induced kinetic asymmetry leading to the class selection of biaxial disclinations can develop into novel materials in topologically similar systems. Other than NLC, the presented work has resemblance with line defects in passive 2,3,9 and active 10,45 soft matter including conducting microwires 46 and self-assembled resonators 47 , and thus has the potential to bring exciting applications in diverse systems.

Methods
Fluctuating electronematics: model energy and thermal kinetics. Instantaneous orientational order, that distinguishes between the disordered liquid state and partially ordered nematic state, is characterized by a symmetric traceless second rank tensor 1 Figure 7. A portion of the xy-slice plane of an 80 2 × 160 μm 3 thermotropic biaxial NLC (portrayed in frame 1 e,g) displaying ±1/2 integer defects corresponding to two different class of disclinations. Note that singular points of homotopy class C y are markedly different than C z by ∂n = 0, otherwise for both classes, (∂S,∂B 2 ,∂l) ≠ 0. Material (computation) parameters are tabulated in Table 1 2 . Y l m are the spherical harmonics of degree l with order m, 〈⋅〉 t denotes the ensemble average at that instant t and symmetric traceless is symbolized with′ ⋅ ′. While B 2 cease to null as S attains its maximum, other moments R 2 generated from the planar projection of {n, l} is also exercised to define biaxiality 33 . B 2 is numerically less expensive in considering one less degree of freedom, nevertheless, R 2 plays an important role in understanding the field-induced switching kinetics 48 .
The ground state free energy including the excitations due to the spatial distortions and dielectric coupling are represented in the phenomenological Ginzburg-Landau-de Gennes (GLdG) free energy functional beq b eq beq b eq beq 2 3 1/2 2 1/2 and clumsy algebraic expression for S beq is omitted for brevity.
Inhomogeneities due to the excitations above the ground state is concealed within the symmetry allowed lowest order terms in  , where ratio κ = L 2 /L 1 and Θ = L 3 /L 1 of the elastic constants can be mapped to the Frank constants splay, bend and twist 49 . For MBBA, (κ,Θ) ≈ 1 while for 5CB, κ ≈ 40,Θ ≈ 1 which designate nearly equal splay and bend constants for both materials (Θ ≠ 0) but the twist constant is an order smaller than splay or bend. This results into nucleation of integer topological charged nematic droplets in the metastable isotropic medium of 5CB, while topologically uncharged nematic droplets nucleate in the metastable medium of MBBA 20,50 .
The optical dielectric permittivity tensor is related to the orientation tensor when separated into symmetric and antisymmetric part ε = ε s δ + ε a Q, where δ is the Kronecker delta, ε s = Trε/3 and ε ε ε = − ⊥ 2( )/3 a with ε ε ⊥ ( ) being the permittivity along (orthogonal to) the director. The material parameters are ε s = 0.74ε a , where for 5CB ε a = 5.8 at T = 33.65 °C and ε a = −0.7 for MBBA at T = 25 °C 51,52 . Application of a spatially varying electric field E =−∂Ψ leads to an electric displacement D = ε ⋅ E and therefore to a dielectric (free) energy term with ε 0 being the vacuum permittivity and Ψ the electric potential. We characterize the intensity of the electric field with respect to the order 19,53 and thermal energy by the nondimensional ratios ε 1 = (ε 0 ε a E 2 /8πAS ueq ) 1/2 and ε 2 = (ε 0 ε a E 2 /8πk B T) 1/2 . We estimate the Fréedericksz threshold to orient a director in a uniformly oriented nematic state along (orthogonal to) an electric field in a twist geometry, calculated by minimizing the free energy for director distortion and field coupling, to be 54 where L x is the spacing between the electrodes. In the present study, confinement effects e.g. centrosymmetry breaking geometric restriction at the boundary by coverslips, patterned or chemically active walls, curvature induced polarization and the presence of free ions are not considered. For bent-core molecules, curvature induced polarization can be incorporated by adding where c 1,…,4 are coefficients 55 . Free ions can also be neglected by retaining x 0 2 , where n is the free ion density and e the electric charge 19 .
When the medium is sufficiently dry so that the long ranged hydrodynamic interaction produced by the motile disclinations do not interfere the kinetics and fluid inertia plays no role -thus restricting to an overdamped relaxational kinetics without convection of momentum, the Langevin equation displaying the time evolution of the electric potential together with the orientation tensor can be written as where Γ = Γ[δ ik δ jl + δ il δ jk − 2δ ij δ kl /3] is a 4 th rank tensor that maintains the symmetric-traceless property on the right-hand side of second equation (5). The rotational diffusion constant Γ is approximated to be independent of Q and the stochastic term ξ satisfies the property of the orientation tensor 〈ξ(x, t)〉 = 0, 〈ξ(x, t)ξ(x′, t′)〉 = 2k B TΓδ(x − x′)δ(t − t′) and is constructed as a summation of Wiener process to keep discrete fluctuation dissipation (FDT) spectrum intact over all Fourier modes and thus to sample Gibbs distribution in thermal equilibrium 26 . We stress at this point that the functional derivative of with respect to Ψ is the divergence of the electric displacement D 38 . This legitimate the imposition of Maxwell's equation along with the Q-tensor equation. Unlike ref. 40 , we neglect cross-coupling between terms proportional to δ δψ / total  to the right-hand side of ∂ t Q equation and  δ δ Q / total to the right-hand side of ∂ t ψ equation for simplicity. By substituting equation (3) in equation (5), the coupled Maxwell-GLdG equation in expanded form reads and the electric field can remain spatially uniform by decoupling from the Q-tensor equation (6) to yield, t ij i j i j ij mn m n ij i qr j qr i l jl a i j ij 2 32 An extensive numerical route of investigation is presented next. The deterministic part of the equation (7) has been widely exercised in two-dimensional monolayered thin films in one elastic constant approximation (κ = Θ = 0) to examine electrokinetic effects of point defects in switching experiments 27,56 . However, an impeccable role of backflow and elastic anisotropy is found to decipher asymmetric speed of ±1/2 integer defects under intense electric field 23,30 , that also leads to ceasing motility of disclinations due to the backflow 6 . A quantitative measure in three-dimensional media under intense electric field, albeit experimentally posed for single disclination in deep uniaxial state 6 , is yet to be sought by including Beris-Edwards model 41 for fluid flow to the presented fluctuating electronematics model. Here instead, we focus on moderate to the small magnitude of electric forces including thermal fluctuations at temperature close and below T * where advective effects can be neglected for simplicity.
Stochastic method of lines for fluctuating electronematics. We consider a thick rectangular slab of insulating thermotropic NLC material in thermal equilibrium, with an 80 μm × 160 μm base and 80 μm height. At equilibrium, periodic boundaries in three Cartesian directions are retained, that can be realized as a free standing anisotropic thick film from a groove. The electric potential at x = 0 is fixed at zero and at x = L x is held according to the desired field magnitude with Dirichlet boundary condition and periodicity is maintained in the yz-directions. This is mimic by suspending the thick slab within two planar laser beams kept at a different potential or placing within an electrode without the notion of an easy axis.
An elegant way to numerically integrate equations (6,7) while retaining the symmetric-traceless property of {Q, ξ} tensor is by projection on a basis of five 3 × 3 matrices as 26,57 Q = ∑ l a l T l ;ξ = ∑ l ζ l T l (l = 1, …, 5), so that the kinetics is projected into the basis coefficients a l . The fluctuating force in the projected basis has the property, 〈ζ l (x, t)〉 = 0,〈ζ l (x, t)ζ m (x′, t′)〉 = 2k B TΓδ lm δ(x − x′)δ(t − t′). After projection, equation (6)  while equation (7)   To numerically integrate the above equations (8 or 9), we adopt a central finite differencing to spatially discretize the Laplacian and mixed derivates. An explicit stochastic method of lines (SMOL) integrator is exercised for seamless temporal integration 58 . SMOL is a stochastic generalization of the deterministic method of lines approach 59 that relies on discretizing the spatial derivates without a temporal discretization, thus yielding to a set of ordinary differential equations in time, that can be easily integrated using the standard numerical libraries 60 . Though spatial accuracy can be increased by using spectral collocation method 61 , obtained solution is usually limited by the temporal accuracy of the integrator. SMOL is numerically stable without computational hindrance with convergence, is less computationally overloaded, spatiotemporally second order accurate, satisfies discrete FDT and can faithfully reproduce lab-based experiments in silico 7,13,20,29,30 .
In coarsening kinetics, a kinetic length scale is extracted from the scattered light intensity inscribed within the structure function S(q,t) as ξ = [∫d 3 qS(q,t)q 2 /∫d 3 qS(q,t)] −1/2 , where S(q, t) = Q ij (q, t)Q ij (−q, t)/∫d 3 qQ ij (q, t)Q ij (−q, t) with the orientation tensor defined as  62 . At equilibrium, the Ornstein-Zernicke form of the spatial correlation 〈Q(0)Q(x)〉 ∝ Γe −|x|/ζ /|x| defines the coherence length ζ = [{1 + 2(κ + Θ)/3}L 1 /A] 1/2 . By definition, ζ determines the core size of a disclination which is distinct from the length scale ξ that denotes the mutual separation between two strings and also the separation between ±1/2 integer dipoles within a shrinking disclination loop (or the radius of curvature). Thus, the length scale obtained from the disclination surface density ρ ξ − 2 shown in frame 3 g is identical with the length scale obtained from correlation functions 16 . In our numerics, time, length, and energy scales are resolved by non-dimensionalizing equations (8,9) 20 and strictly retaining (i) Courant-Friedrichs-Lewy condition for timestep to avoid stiffness 63 , (ii) ζ Δ  x to resolve the grid spacing and distortion length embedded in the disclination core, (iii) ∂S < S/ζ for the validity of the GLdG method and, (iv)  k T B barrier height between isotropic-nematic state at clearing point to avoid any spurious oscillations between the two phases, thus, is the non-dimensionalized distortion free energy 16 . With the chosen values of {ε 0 , ε s , ε a }, electric field relaxation is rapid compared to the orientational kinetics, so as a steady electric field is obtained for each step of Q-evolution. Grid size independence, numerical accuracy, and validity of physical tests were confirmed for each presented graphics and for both equations (8,9).
We finally discuss difference of our method with existing methodology in connection with the electrorheology of disclinations. Within Leslie-Ericksen (LE) theory in one dimension using free energy minimization technique, both the effect of rheology 64 and nonuniform electric field 40 has been exercised. Traditionally, LE theory is suitable in deep within the nematic phase where it is assumed that {S, B 2 , l} are constant in the orientation tensor and elastic distortion is concealed in n only. Earlier attempts to include hydrodynamics 65 to LE theory were through the finite element method (FEM) 66 . Also, the finite volume method (FVM) applied to Smoluchowski equation for the orientational probability density function is yet at a preliminary stage 67 . On the other hand, after the advent of large scale computation, numerical solution of the athermal orientation tensor equation retaining all independent degrees of freedom is exercised through cell-dynamic scheme (CDS) 18 , finite difference methods (FDM) 35,68 as well as the method of lines (MOL) 59 approach to the GLdG theory. Attempts to include hydrodynamics and thermal fluctuations to the FDM via fluid particle dynamics (FPD) 35 were replaced with the Lattice-Boltzmann method (LBM) 37,69 , while inclusion of thermal fluctuations retaining second-order numerical accuracy is extended as a stochastic generalization of the method of lines (SMOL) 26 . In this manuscript, SMOL complements the inclusion of Maxwell's equation to the existing formulation and motivates the association of hydrodynamics 24 for intense electric field studies and Fourier's law of latent heat conduction 70 for confined NLC systems.
Identification and topological classification of disclinations. Disclinations of different homotopy class are obtained after extracting {S,B 2 ,n,l} from the basis coefficients a i on each space point by a similarity transformation. The scalar values are colour rendered according to the indicated bars in Fig. 3. Disclinations are identified by sketching the isosurfaces with specific isovalue of the scalar fields. At a late stage of the kinetics when few isolated disclinations are existent, a plane from the three-dimensional volume is sliced in which the lateral section of disclinations appear as points. Once the physical location of disclinations are identified, the homotopy class, charge and sign of defect is calculated from the spatial distribution of {n,l}. Similar to the Volterra construction in crystal dislocation, a Burgers circuit γ displayed in Fig. 7 is constructed using neighboring lattice points that encircle the defect 18 . The angular shift of {n,l} is measured while traversing one complete loop and the charge and sign of defects are estimated using the hodograph method 1 . If u = (sinθ cosφ, sinθ sinφ, cosθ) is directionally equivalent with −u = (sinθ′cosφ′, sinθ′sinφ′, cosθ′) where (θ, φ) are polar and azimuthal angle in an arbitrary frame, then following transformation retains the centrosymmetry, θ π θ φ π φ ′ → − ′→ + . ; (10) For ±1/2 integer disclinations in a slice plane of uniaxial NLC shown in frame 1c, while traversing γ, n rotates by ± π. In biaxial NLC, homotopy classes are identified using the following recipe: 18 (i) for C x class of disclinations, n rotates by ±π but l does not rotate, (ii) for C y class of disclinations, n does not rotate but l rotates by ± π and, (iii) for C z class of disclinations, both {n, l} rotates by ±π. We do not find any C x class of disclinations, which is consistent with the analytic prediction on two-dimensional nonabelian vortices 17 and numerical computations 16,18 . This algorithm not only supersedes the traditional defect classifying approaches using vector field, tensor glyph or hyperstreamline seeding through Mueller and Westin matrices 71,72 but can uniquely determine disclinations from the structure of orientation tensor rather than rely on the vectorial information.