Flowing Liquid Crystal Torons Around Obstacles

Abstract

Liquid crystal torons, localized topological structures, are known for their stability and dynamic behaviour in response to external stimuli, making them attractive for advanced material applications. In this study, we investigate the flow of torons in chiral nematic liquid crystals around obstacles. We simulate the fluid flow and director field interactions using a hybrid numerical method combining lattice Boltzmann and finite difference techniques. Our results reveal that the toron dynamical behaviour depends strongly on the impact parameter from the obstacle. At impact parameters smaller than half cholesteric pitch, the flowing toron is destabilized by the interaction with the obstacle; otherwise, the flowing toron follows a trajectory with a deflection which decays exponentially with the impact parameter. Additionally, we explore the scattering of torons by multiple obstacles, providing insights into how the dynamics of these structures respond to complex environments.

1. Introduction

Liquid crystals (LCs) represent a fascinating class of soft materials characterized by their fluidity combined with anisotropic properties, making them highly responsive to external fields and boundary conditions. This unique combination has led to widespread applications, particularly in display technologies, and has also spurred interest in fundamental research on the stability of their phases and topological structures [1]. Among these structures, solitons in chiral nematic LCs have emerged as a significant area of study due to their rich topological properties and potential for novel technological applications [2,3]. Besides skyrmions, which are two-dimensional configurations, many other topological-protected structures were realized in experiments, including torons, hopfions, and solitonic macromolecules. We will refer to torons as the three-dimensional elementary skyrmion terminating at two-point defects to satisfy the uniform surface boundary conditions and match the topologically nontrivial skyrmion tube with the uniform far-field background and reserve skyrmion to describe the two-dimensional mid-plane structure.

Solitons in liquid crystals are localized, non-singular distortions of the director field that exhibit robust particle-like behaviour and complex interactions. These structures are of particular interest due to their ability to be manipulated by external stimuli, such as electric fields, which can drive their motion and induce collective behaviour [4,5]. The study of soliton dynamics in LCs has provided insights into fundamental aspects of topology in soft matter and has highlighted their potential in advanced material applications, including reconfigurable optical devices and microfluidic systems [4,6,7,8,9,10,11,12,13].

While much of the research on LC torons has focused on their behaviour under uniform external fields, little is known under non-uniform conditions such as static colloids or other obstacles trapped in the LC matrix. Furthermore, only a handful of studies have addressed the interaction between torons and mass flows [14,15,16,17], which has thus remained an overlooked problem of relevance in various fields. Understanding the dynamical response of LC torons in non-uniform or heterogeneous media is relevant for practical applications, where flowing torons may encounter various types of obstacles, such as geometrically patterned substrates, dispersed colloidal particles or defects in the LC matrix. This response will significantly alter the torons’s initial trajectory, in addition to its shape and stability, potentially driving novel dynamical phenomena.

In the theoretical study of torons, three approaches are commonly used in the literature. The first approach involves the numerical solution of continuum equations resulting from the minimization of the elastic free energy, which is typically implemented using finite difference or finite volume methods [4,18,19,20]. This approach has focused either on the dynamic evolution of torons or on the stable configuration of static torons. The second approach is based on analytical studies of the continuum equations, under simplified conditions required to make the problems tractable [21,22]. This method often provides valuable insights into the toron behaviour in idealized situations. The third approach, a more recent development, uses particle-based modelling, where torons are treated as individual particles [23]. Langevin dynamics is employed to simulate their motion, and potential functions for the toron interactions are inferred from experimental data or from simulations using the first method. Moreover, there are attempts to use neural networks to measure system parameters using optical images obtained from experiments on skyrmions [24]. In this study, we use the first approach, which is combined with the lattice Boltzmann method, to simulate material flow, providing a comprehensive and versatile framework to study the interaction between topology and external mass flows that may be perturbed by an arbitrary array of obstacles.

In what follows, we study the flow of liquid crystal torons, driven by an external mass flow, around solid obstacles. The toron’s stability and trajectory are found to depend strongly on the impact parameter of the obstacle. When the impact parameter is less than $h<5\mathsf{\mu }\mathrm{m}$, the toron becomes unstable and disintegrates, and flows where torons would otherwise remain stable. However, when the impact parameter is larger than $h\ge 5\mathsf{\mu }\mathrm{m}$, the stability of the toron is preserved while the toron’s trajectory is deflected by an amount that decays exponentially with the impact parameter. Notably, the toron’s velocity is found to increase as it approaches the obstacles due to flow constriction. We have also examined the dynamics and stability of torons flowing around one and between two obstacles, providing additional insights into how multiple obstacles may influence the dynamics of flowing torons.

2. Method

2.1. Hydrodinamic Equations

The dynamics of the liquid crystal (LC) director field are described by the Ericksen–Leslie model [25,26,27,28]. This model couples two dynamical equations: one for the material flow and another for the director field. These equations are particularly adequate to describe the behaviour of LCs in the nematic or cholesteric phases.

For the velocity field, we use the Navier–Stokes equation along with the continuity equation [29]:

(1)$\begin{array}{ccc}& & \rho {\partial }_t{u}_{\alpha }+\rho u_{\beta }{\partial }_{\beta }{u}_{\alpha }={\partial }_{\beta }\left[-p{\delta }_{\alpha \beta }+{\sigma }_{\alpha \beta }^v+{\sigma }_{\alpha \beta }^e\right]+\rho g\hfill \end{array}$

(2)$\begin{array}{ccc}& & {\partial }_{\alpha }{u}_{\alpha }=0,\hfill \end{array}$

where the viscous stress tensor is defined as

(3)$\begin{array}{cc}\hfill {\sigma }_{\alpha \beta }^v=& {\alpha }_1{n}_{\alpha }{n}_{\beta }{n}_{\mu }{n}_{\rho }{D}_{\mu \rho }+{\alpha }_2{n}_{\beta }{N}_{\alpha }+{\alpha }_3{n}_{\alpha }{N}_{\beta }\\ & +{\alpha }_4{D}_{\alpha \beta }+{\alpha }_5{n}_{\beta }{n}_{\mu }{D}_{\mu \alpha }+{\alpha }_6{n}_{\alpha }{n}_{\mu }{D}_{\mu \beta }.\end{array}$

In these equations, $\rho$ denotes the fluid density, g is the external acceleration that drives the fluid flow, p is the pressure, $\mathit{u}$ is the fluid velocity, $\mathit{n}$ is the director field describing the direction of alignment of the LC molecules, and ${\alpha }_n$ is the Leslie viscosity. The kinematic transport, which describes the influence of the macroscopic flow field on the microscopic structure, is given by

(4)$\begin{array}{c}\hfill N_{\beta }={\partial }_t{n}_{\beta }+u_{\gamma }{\partial }_{\gamma }{n}_{\beta }-W_{\beta \gamma }{n}_{\gamma },\end{array}$

while the shear rate and vorticity tensors are defined as

(5)$\begin{array}{c}\hfill D_{\alpha \mu }={\displaystyle \frac{1}{2}}\left({\partial }_{\alpha }{u}_{\mu }+{\partial }_{\mu }{u}_{\alpha }\right),W_{\alpha \mu }={\displaystyle \frac{1}{2}}\left({\partial }_{\mu }{u}_{\alpha }-{\partial }_{\alpha }{u}_{\mu }\right).\end{array}$

The elastic stress tensor is

(6)$\begin{array}{c}\hfill {\sigma }_{\alpha \beta }^e=-{\partial }_{\alpha }{n}_{\gamma }{\displaystyle \frac{\delta E}{\delta \left({\partial }_{\beta }{n}_{\gamma }\right)}},\end{array}$

where E represents the Frank–Oseen elastic free energy:

(7)$\begin{array}{c}\hfill E=\int dV\left({\displaystyle \frac{K_{11}}{2}}{(\nabla ·\mathit{n})}^2+{\displaystyle \frac{K_{22}}{2}}{\left(\mathit{n}·[\nabla \times \mathit{n}]+q_0\right)}^2+{\displaystyle \frac{K_{33}}{2}}{[\mathit{n}\times [\nabla \times \mathit{n}]]}^2\right).\end{array}$

where $K_{11}$, $K_{22}$ and $K_{33}$ are the Frank elastic constants, and $q_0=2\pi /P$, with as P the cholesteric pitch.

The second set of equations describes the evolution of the director field:

(8)$\begin{array}{cc}\hfill & {\partial }_t{n}_{\mu }={\displaystyle \frac{1}{{\gamma }_1}}{h}_{\mu }-{\displaystyle \frac{{\gamma }_2}{{\gamma }_1}}{n}_{\alpha }{D}_{\alpha \mu }-u_{\gamma }{\partial }_{\gamma }{n}_{\mu }+W_{\mu \gamma }{n}_{\gamma },\hfill \end{array}$

where ${\gamma }_1={\alpha }_3-{\alpha }_2$ is the rotational viscosity determining the relaxation rate of the director, and ${\gamma }_2={\alpha }_3+{\alpha }_2$ is the torsion coefficient characterizing the contribution to the viscous torque from the velocity field gradients. The ratio ${\gamma }_2/{\gamma }_1$ is known as the alignment parameter, with $|{\gamma }_2/{\gamma }_1| > 1$ indicating flow-aligning and $|{\gamma }_2/{\gamma }_1| < 1$ flow-tumbling systems. The molecular field is given by

(9)$\begin{array}{c}\hfill h_{\mu }=-{\displaystyle \frac{\delta E}{\delta n_{\mu }}}.\end{array}$

The simulations employed a hybrid numerical method. The velocity field is resolved using the lattice Boltzmann method [30,31], summarized in the next subsection, with the elastic and viscous stress tensors (except the term proportional to ${\alpha }_4$) introduced as force terms. The director field equation, Equation (8), was solved using a predictor–corrector finite difference algorithm [32]. The derivatives were calculated using central moment differences:

(10)$\begin{array}{c}\hfill {\displaystyle \frac{d\mathcal{F}}{dx}}={\displaystyle \frac{\mathcal{F}(x+\mathsf{\Delta }x)-\mathcal{F}(x-\mathsf{\Delta }x)}{2\mathsf{\Delta }x}}+\mathcal{O}\left(\mathsf{\Delta }{x}^2\right),\end{array}$

where $\mathcal{F}$ is a generic function and $\mathsf{\Delta }x$ is the step. Note that the gradients are calculated only in the fluid nodes while keeping the quantities in the solid nodes constant. On solid boundaries, including the obstacles, infinite homeotropic anchoring and no-slip conditions were applied using the bounce-back condition [30].

The simulations started with the liquid at rest, and the directors mostly aligned perpendicular to the plates except near the toron, whose configuration was obtained by minimizing its free energy starting from the Ansatz of Ref. [14]. The material parameters were chosen to be close to those of MBBA at 22 °C [33], except for the absolute viscosity (or equivalently, ${\alpha }_4$), which was doubled to ensure reasonable simulation times while keeping its value of the same order of magnitude of the material one. Note that the physical behaviour is mainly controlled by the Ericksen number [14]: $Er\equiv {\displaystyle \frac{\mu v_{t}P}{K}}$, where $v_t$ is the toron velocity, P is the cholesteric pitch, $\mu$ is the viscosity, and K is the average elastic constant. If we multiply $\mu$ by two, $v_t$ must be divided by two in order to have the same Ericksen number and, thus, the same behaviour. The code was parallelized in CUDA-C, and the simulations were executed on GPUs [34]. A typical performance of the simulation is ∼328 MLUPS (Mega Lattice Updates Per Seconds) in a Nvidia Tesla V100 using double precision calculations.

The use of GPUs allowed us to perform 3D simulations in a feasible time frame. By comparison, our previous C++ implementation, parallelized with OpenMP and running on 8 CPU cores, was 2 orders of magnitude slower. Additionally, we found that the current method, based on the Ericksen–Leslie theory, is more suitable for simulating flowing liquid crystal torons than the Beris–Edwards model (as in Ref. [16]). It uses less memory (storing only the 3 components of the director field instead of the 6 components of the Q-tensor) and assumes a constant scalar order parameter, reducing spurious effects such as non-uniform far fields.

2.2. Lattice Boltzmann Method

We obtain the velocity field by solving the discretized Lattice Boltzmann Equation (LBE) using the single relaxation time (SRT) collision approximation and Guo’s forcing scheme [35], as follows:

(11)$\begin{array}{cc}\hfill & f_i(\mathit{x}+{\mathit{e}}_i\mathsf{\Delta }t,t+\mathsf{\Delta }t)-f_i(\mathit{x},t)=\hfill \\ & -{\displaystyle \frac{\mathsf{\Delta }t}{\tau }}(f_i(\mathit{x},t)-f_i^{\mathrm{eq}}(\mathit{x},t))+\mathsf{\Delta }t\left(1-{\displaystyle \frac{1}{2\tau }}\right)S_i^F(\mathit{x},t),\hfill \end{array}$

where $f_i$ is the particle distribution function, and $f_i^{\mathrm{eq}}$ is the equilibrium distribution function, which is determined by

(12)$\begin{array}{c}\hfill f_i^{\mathrm{eq}}\left(\mathit{x}\right)=w_i\rho \left(\mathit{x}\right)\left[1+{\displaystyle \frac{{\mathit{e}}_i·\mathit{u}\left(\mathit{x}\right)}{c_s^2}}+{\displaystyle \frac{{({\mathit{e}}_i·\mathit{u}\left(\mathit{x}\right))}^2}{2c_s^4}}-{\displaystyle \frac{\mathit{u}\left(\mathit{x}\right)·\mathit{u}\left(\mathit{x}\right)}{2c_s^2}}\right],\end{array}$

where $c_s$ is the lattice speed of sound, defined as $c_s={\displaystyle \frac{\mathsf{\Delta }x}{\sqrt{3}\mathsf{\Delta }t}}$. $\mathsf{\Delta }x$, $\mathsf{\Delta }t$, and the reference density ${\rho }_0$ are set to unity (lattice units). $w_i$ and ${\mathit{e}}_i$ are the lattice weights and the discrete velocity vector in the $i^{\mathrm{th}}$ direction, respectively. In this work, we use the $D_3{Q}_{19}$ lattice arrangement; hence, $w_i$ and ${\mathit{e}}_i$ are given by

(13)${\mathit{e}}_i=\left(\begin{array}{ccccccccccccccccccc}0& 1& -1& 0& 0& 0& 0& 1& -1& 1& -1& 1& -1& 1& -1& 0& 0& 0& 0\\ 0& 0& 0& 1& -1& 0& 0& 1& 1& -1& -1& 0& 0& 0& 0& 1& -1& 1& -1\\ 0& 0& 0& 0& 0& 1& -1& 0& 0& 0& 0& 1& 1& -1& -1& 1& 1& -1& -1\end{array}\right),$

(14)$w_i=\left(\begin{array}{ccccccccccccccccccc}{\displaystyle \frac{1}{3}}& {\displaystyle \frac{1}{18}}& {\displaystyle \frac{1}{18}}& {\displaystyle \frac{1}{18}}& {\displaystyle \frac{1}{18}}& {\displaystyle \frac{1}{18}}& {\displaystyle \frac{1}{18}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}& {\displaystyle \frac{1}{36}}\end{array}\right).$

The macroscopic fluid density is defined as

(15)$\rho \left(\mathit{x}\right)=\sum _i{f}_i\left(\mathit{x}\right).$

The relaxation time ($\tau$) is related to the kinematic viscosity of the fluid ($\nu$) through

(16)$\nu =c_s^2\left(\tau +{\displaystyle \frac{1}{2}}\right)\mathsf{\Delta }t.$

The forcing term in the LBE equation $S_i^F$ is defined as

(17)$\begin{array}{c}\hfill S_i^F\left(\mathit{x}\right)=w_i\left[{\displaystyle \frac{{\mathit{e}}_i·\mathit{F}\left(\mathit{x}\right)}{c_s^2}}+{\displaystyle \frac{({\mathit{e}}_i{\mathit{e}}_i-c_s^2\mathit{I}):(\mathit{u}\left(\mathit{x}\right)\mathit{F}\left(\mathit{x}\right)+\mathit{F}\left(\mathit{x}\right)\mathit{u}\left(\mathit{x}\right))}{2c_s^4}}\right],\end{array}$

where $\mathit{F}$ represents the resultant force acting on the flow field (external force and stresses arising from the liquid crystal distortions), and $\mathit{u}$ is the macroscopic fluid velocity, which is defined as

(18)$\mathit{u}\left(\mathit{x}\right)={\displaystyle \frac{1}{\rho \left(\mathit{x}\right)}}\sum _i{\mathit{e}}_i{f}_i\left(\mathit{x}\right)+{\displaystyle \frac{\mathsf{\Delta }t}{2\rho \left(\mathit{x}\right)}}\mathit{F}\left(\mathit{x}\right).$

The bounce-back boundary condition is used to model the no-slip boundary condition over the solid nodes [30]. The idea behind the bounce-back boundary condition is to reflect the distribution functions hitting a solid node with the same magnitude in the reverse direction. The exact location of the solid boundary is assumed to be between two lattice points [30].

3. Results

We report the results of numerical simulations of a flowing toron driven by an external force, around a cylindrical pillar near its initial trajectory. The liquid crystal was initialized at rest, with alignment along the vertical direction, except at the toron, where we used the relaxed configuration obtained from the Ansatz [14]. The system dimensions were set to $L_x\times L_y\times L_z=96\times 56\times 16$, and the material parameters are listed in Table 1. The pillar, with a radius of $R=3.13\mathsf{\mu }\mathrm{m}$ (or 5 $\mathsf{\Delta }x$), homeotropic anchoring and no-slip boundary, was placed at $x=L_x/2$ and $y=L_y/2+h$. The fluid was driven by an external body force of magnitude $g=312.5$ m/s² (or $2\times {10}^{-9}$ l.u.) in the x direction.

As expected, the pillar distorts the velocity field and perturbs the toron initial trajectory, as shown in Figure 1. In the mid-plane, Figure 1a, the toron shrinks slightly as it flows around the pillar but later recovers its original configuration. To confirm the toron’s stability post-pillar interaction, we relaxed the final configuration in Figure 1 without the external force, confirming that it remains stable (not shown). The toron cross-section in Figure 1b reveals that the toron becomes asymmetric in the direction of the flow, consistent with previous observations in Poiseuille-like flows without the pillar [16]. The velocity field is perturbed by the toron, while the pillar introduces large velocity gradients in the flow field (Figure 1c). In the cross-sectional view, Figure 1d, the flow exhibits a Poiseuille-like profile, except near the toron and the pillar.

The deflection of the initial trajectory of the toron depends on the impact parameter h, as shown in Figure 2. Trajectories are deflected significantly when the impact parameter is less than the cholesteric pitch. However, below a threshold, $h=3.75 \mathsf{\mu }\mathrm{m}$ or $0.43$ P, in our simulations, the toron becomes unstable and disintegrates. The toron in Poiseuille flows already becomes unstable without the obstacle at higher speeds [17], and here, it occurs due to the local increase in the velocity field around the obstacle. When the toron collides directly with the pillar ($h=0$), it disintegrates before flowing around the obstacle. For small but non-zero values of h, the toron flows around the pillar but later becomes unstable, as depicted in Figure 3.

The inset on the left of Figure 2 depicts the deflection $\mathsf{\Delta }y$ of the toron initial trajectory as a function of the impact parameter, h. As h increases, $\mathsf{\Delta }y$ decays exponentially due to elastic distortions of the director field. Exponential decays are usually observed around defects or distortions in liquid crystal textures [36]. However, it is clear that the decay reported here does not tend to zero as h tends to infinity, as expected for elastic distortions only. This is likely due to the finite size of the flow domain in the y direction. The constriction in the velocity field, caused by the obstacle at $x=L_x/2$, leads to a velocity increase in this region, which in turn leads to a deflection. The inset on the right shows the toron velocity as a function of position, where we note that a stable toron accelerates near the pillar, with its velocity increasing by about 26% before returning to its initial value, due to the flow constriction by the obstacle. Notice that this increase in the velocity field occurs even at the largest value of h.

The free energy integrated over the entire volume, calculated using Equation (7), is plotted in Figure 4 as a function of time for various impact parameters, h. At small h, the free energy rises sharply, reaching a peak when the toron loses stability and disintegrates. At larger impact parameters, h, the toron remains stable, and the free energy stabilizes, with slight perturbations occurring as the toron flows around the pillar, particularly evident in the orange curve ($h=5\mathsf{\mu }\mathrm{m}$) around $t=1$ s.

Until now, we considered the flow of a toron around a single pillar. To explore the effects of multiple obstacles, we simulated the toron flow between two pillars, which act as a flow constriction (see Figure 5). In this configuration, the toron accelerates by 44% as it flows between the pillars, as shown in Figure 5c. Similar to the single-pillar case, the toron becomes unstable at small values of h. In our simulations, the toron remains stable for $h>6.25 \mathsf{\mu }\mathrm{m}$, which is almost twice the impact parameter required for toron stability in the single-pillar case. This is likely due to the toron’s size reduction when passing around the pillars, an effect that is amplified by the presence of two obstacles. Additionally, in this configuration, the toron does not deviate from its initial trajectory.

In addition to simulating the interaction between the toron and a cylindrical pillar with homeotropic anchoring, we investigated other obstacle geometries (not shown). For example, we replaced the pillar with a solid sphere, mimicking a colloid, and applied different anchoring conditions (both homeotropic and toron-like). However, no qualitative differences were observed. Thus, we conclude that the primary factor determining the flow of torons around obstacles, as well as their stability, appears to be the effect of the large velocity gradients introduced by the obstacles rather than their detailed form.

4. Conclusions

In this study, we investigated the interaction between flowing liquid crystal torons and solid obstacles. Through numerical simulations using a combination of the lattice Boltzmann method and finite differences, we demonstrated how the dynamics of flowing torons respond to the presence of obstacles in the flow domain. Our findings reveal that the toron’s stability and trajectory depend strongly on the impact parameter of the obstacle. At short impact parameters (smaller than half cholesteric pitch), the toron becomes unstable and disintegrates, while at larger impact parameters, the toron’s trajectory is deflected around the obstacle. The deflection of the trajectory and the increase in velocity as the toron approaches the obstacle depend strongly on the impact parameter. In addition, we explored the dynamics of torons flowing between two obstacles, providing insights into how multiple obstructions influence the toron´s stability and motion. An extension of this work could address the study of flowing torons in a lattice or in disordered arrays of obstacles. Moreover, it is also interesting to investigate the dynamics of flowing torons in gradient velocity fields. Symmetry breaking resulting in deflections of the toron trajectories in directions perpendicular to the drive, as reported in Couette-like flows [15], may be a general phenomenon that is also observed here. An extension of this work could involve the study of flowing torons in a lattice of obstacles and in a disordered distribution of them. Moreover, it would be interesting to investigate the dynamics of torons in gradients of velocity field in the mid-plane. Possibly, they may break the symmetry and move to transverse directions like in Couette-like flows [15].

This study contributes to our understanding of toron interactions in complex environments which could be used in applications in reconfigurable optical devices and microfluidic systems.

Author Contributions

Conceptualization, R.C.V.C. and M.M.T.d.G.; Methodology, R.C.V.C.; Software, M.S. and R.C.V.C.; Investigation, J.P.A.S. and R.C.V.C.; Writing—original draft, R.C.V.C.; Writing—review & editing, M.S. and M.M.T.d.G.; Visualization, J.P.A.S.; Supervision, R.C.V.C. and M.M.T.d.G.; Project administration, R.C.V.C. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data that support the findings of this study are available from the corresponding author, upon reasonable request.

Conflicts of Interest

The authors declare no conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Table

Figure 1. Snapshots of a toron flowing close to a pillar: (a) z-component of the director field in the mid-xy plane; (b) z-component of the director field in the xz plane that intercepts the geometric centre of the toron; (c) magnitude of the velocity field in the mid-xy plane; and (d) the magnitude of the velocity field in the xz plane that intercepts the geometric centre of the toron. The white arrow indicates the direction of the flow. The impact parameter of the obstacle is [Forumla omitted. See PDF.].

Figure 2. Trajectories of torons deflected by pillars at different impact parameters, i.e., the lateral distance between the toron’s initial trajectory and the pillar’s centre, h (legend in [Forumla omitted. See PDF.]). The colour of the pillars is the same as that of the corresponding trajectories. Torons that are destabilized and destroyed by the impact are indicated in the figure. The inset on the left depicts the dependence of the vertical deflection of the toron trajectory, [Forumla omitted. See PDF.], on the impact parameter, h. The solid line is an exponential fit to the data. The inset on the right shows the toron velocity as a function of the position.

Figure 3. Unstable toron flowing around a pillar: (a) z-component of the director field in the xy mid-plane; and (b) z-component of the director field in the xz plane that follows the geometric centre of the toron. The white arrow indicates the direction of the flow. The impact parameter of the obstacle is [Forumla omitted. See PDF.].

Figure 4. Time evolution of the free energy of flowing torons around pillars, for different impact parameters h in [Forumla omitted. See PDF.].

Figure 5. Snapshots in the mid-plane of a toron flowing through two obstacles. (a) Z-component of the director field. (b) Magnitude of the velocity field. The white arrow indicates the direction of the applied force. Each pillar has an impact parameter [Forumla omitted. See PDF.] [Forumla omitted. See PDF.]. (c) Toron velocity as a function of the position for different values of h in [Forumla omitted. See PDF.]. The vertical lines on the top are the times when the toron flows between the pillars.

Table 1

Parameters used in the simulation and physical units.

Symbol	Simulation Units	Physical Units
$\rho$	1	1088 Kg/m3
$\mathsf{\Delta }x$	1	0.625 $\mathsf{\mu }\mathrm{m}$
$\mathsf{\Delta }t$	1	2$\times {10}^{-9}$ s
$K_{11}$	$1.67 \times {10}^{-7}$	$6.4 \times {10}^{-12}$ N
$K_{22}$	$7.88 \times {10}^{-8}$	$3.0 \times {10}^{-12}$ N
$K_{33}$	$2.62 \times {10}^{-7}$	$9.98 \times {10}^{-12}$ N
${\alpha }_1$	0.0373	0.0036 Pa·s
${\alpha }_2$	−0.4496	−0.044 Pa·s
${\alpha }_3$	−0.0203	−0.0020 Pa·s
${\alpha }_4$	0.9318	0.091 Pa·s
${\alpha }_5$	0.3084	0.030 Pa·s
${\alpha }_6$	−0.1617	−0.016 Pa·s
P	14	8.75 $\mathsf{\mu }\mathrm{m}$
$L_x$, $L_y$, $L_z$	56, 56, 16	35, 35, 10 $\mathsf{\mu }\mathrm{m}$

References

1. Chen, H.W.; Lee, J.H.; Lin, B.Y.; Chen, S.; Wu, S.T. Liquid crystal display and organic light-emitting diode display: Present status and future perspectives. Light Sci. Appl.; 2018; 7, 17168. [DOI: https://dx.doi.org/10.1038/lsa.2017.168] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/30839536]

2. Smalyukh, I.I.; Lansac, Y.; Clark, N.A.; Trivedi, R.P. Three-dimensional structure and multistable optical switching of triple-twisted particle-like excitations in anisotropic fluids. Nat. Mater.; 2010; 9, pp. 139-145. [DOI: https://dx.doi.org/10.1038/nmat2592] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/19966792]

3. Ackerman, P.J.; Trivedi, R.P.; Senyuk, B.; van de Lagemaat, J.; Smalyukh, I.I. Two-dimensional skyrmions and other solitonic structures in confinement-frustrated chiral nematics. Phys. Rev. E; 2014; 90, 012505. [DOI: https://dx.doi.org/10.1103/PhysRevE.90.012505] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/25122322]

4. Ackerman, P.J.; Boyle, T.; Smalyukh, I.I. Squirming motion of baby skyrmions in nematic fluids. Nat. Commun.; 2017; 8, 673. [DOI: https://dx.doi.org/10.1038/s41467-017-00659-5]

5. Sohn, H.R.O.; Smalyukh, I.I. Electrically powered motions of toron crystallites in chiral liquid crystals. Proc. Nat. Acad. Sci. USA; 2020; 117, pp. 6437-6445. [DOI: https://dx.doi.org/10.1073/pnas.1922198117]

6. Fukuda, J.i.; Žumer, S. Quasi-two-dimensional Skyrmion lattices in a chiral nematic liquid crystal. Nat. Commun.; 2011; 2, 246. [DOI: https://dx.doi.org/10.1038/ncomms1250]

7. Posnjak, G.; Čopar, S.; Muševič, I. Points, skyrmions and torons in chiral nematic droplets. Sci. Rep.; 2016; 6, 26361. [DOI: https://dx.doi.org/10.1038/srep26361]

8. Ackerman, P.J.; Smalyukh, I.I. Static three-dimensional topological solitons in fluid chiral ferromagnets and colloids. Nat. Mater.; 2017; 16, pp. 426-432. [DOI: https://dx.doi.org/10.1038/nmat4826]

9. Guo, Y.; Afghah, S.; Xiang, J.; Lavrentovich, O.D.; Selinger, R.L.B.; Wei, Q.H. Cholesteric liquid crystals in rectangular microchannels: Skyrmions and stripes. Soft Matter; 2016; 12, pp. 6312-6320. [DOI: https://dx.doi.org/10.1039/C6SM01190J]

10. Tai, J.S.B.; Smalyukh, I.I. Three-dimensional crystals of adaptive knots. Science; 2019; 365, pp. 1449-1453. [DOI: https://dx.doi.org/10.1126/science.aay1638]

11. Zhao, H.; Tai, J.S.B.; Wu, J.S.; Smalyukh, I.I. Liquid crystal defect structures with Möbius strip topology. Nat. Phys.; 2023; 19, pp. 451-459. [DOI: https://dx.doi.org/10.1038/s41567-022-01851-1]

12. Sohn, H.R.O.; Ackerman, P.J.; Boyle, T.J.; Sheetah, G.H.; Fornberg, B.; Smalyukh, I.I. Dynamics of topological solitons, knotted streamlines, and transport of cargo in liquid crystals. Phys. Rev. E; 2018; 97, 052701. [DOI: https://dx.doi.org/10.1103/PhysRevE.97.052701] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/29906829]

13. Dkhar, D.S.; Kumari, R.; Malode, S.J.; Shetti, N.P.; Chandra, P. Integrated lab-on-a-chip devices: Fabrication methodologies, transduction system for sensing purposes. J. Pharm. Biomed. Anal.; 2023; 223, 115120. [DOI: https://dx.doi.org/10.1016/j.jpba.2022.115120] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/36343538]

14. Coelho, R.C.V.; Tasinkevych, M.; da Gama, M.M.T. Dynamics of flowing 2D skyrmions. J. Phys. Cond. Matter; 2021; 34, 034001. [DOI: https://dx.doi.org/10.1088/1361-648X/ac2ca9] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34607323]

15. Coelho, R.C.V.; Zhao, H.; Amaral, G.N.C.; Smalyukh, I.I.; da Gama, M.M.T.; Tasinkevych, M. “Hall” transport of liquid crystal solitons in Couette flow. arXiv; 2024; arXiv: 2405.10850

16. Coelho, R.C.V.; Zhao, H.; Tasinkevych, M.; Smalyukh, I.I.; Telo da Gama, M.M. Sculpting liquid crystal skyrmions with external flows. Phys. Rev. Res.; 2023; 5, 033210. [DOI: https://dx.doi.org/10.1103/PhysRevResearch.5.033210]

17. Amaral, G.N.C.; Zhao, H.; Sedahmed, M.; Campante, T.; Smalyukh, I.I.; Tasinkevych, M.; da Gama, M.M.T.; Coelho, R.C.V. Liquid crystal torons in Poiseuille-like flows. arXiv; 2024; arXiv: 2409.09486

18. Duzgun, A.; Saxena, A.; Selinger, J.V. Alignment-induced reconfigurable walls for patterning and assembly of liquid crystal skyrmions. Phys. Rev. Res.; 2021; 3, L012005. [DOI: https://dx.doi.org/10.1103/PhysRevResearch.3.L012005]

19. Duzgun, A.; Nisoli, C.; Reichhardt, C.J.O.; Reichhardt, C. Directed motion of liquid crystal skyrmions with oscillating fields. New J. Phys.; 2022; 24, 033033. [DOI: https://dx.doi.org/10.1088/1367-2630/ac58b8]

20. Duzgun, A.; Nisoli, C.; Reichhardt, C.J.O.; Reichhardt, C. Commensurate states and pattern switching via liquid crystal skyrmions trapped in a square lattice. Soft Matter; 2020; 16, pp. 3338-3343. [DOI: https://dx.doi.org/10.1039/C9SM02312G]

21. Alvim, T.; Gama, M.M.T.d.; Tasinkevych, M. Collective variable model for the dynamics of liquid crystal skyrmions. Commun. Phys.; 2024; 7, 2. [DOI: https://dx.doi.org/10.1038/s42005-023-01486-5]

22. Long, C.; Selinger, J.V. Coarse-grained theory for motion of solitons and skyrmions in liquid crystals. Soft Matter; 2021; 17, pp. 10437-10446. [DOI: https://dx.doi.org/10.1039/D1SM01335A] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/34761790]

23. Teixeira, A.W.; Tasinkevych, M.; Dias, C.S. Particle-based model of liquid crystal skyrmion dynamics. Soft Matter; 2024; 20, pp. 2088-2099. [DOI: https://dx.doi.org/10.1039/D3SM01422C] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/38348527]

24. Terro, J.; Tasinkevych, M.; Dias, C.S. Convolutional Neural Network analysis of optical texture patterns in liquid-crystal skyrmions. arXiv; 2024; arXiv: 2410.15246

25. Ericksen, J.L. Hydrostatic theory of liquid crystals. Arch. Ration. Mech. Anal.; 1962; 9, pp. 371-378. [DOI: https://dx.doi.org/10.1007/BF00253358]

26. Leslie, F.M.; Green, A.E. Some thermal effects in cholesteric liquid crystals. Proc. R. Soc. Lond. Ser. A Math. Phys. Sci.; 1968; 307, pp. 359-372. [DOI: https://dx.doi.org/10.1098/rspa.1968.0195]

27. Stewart, I.W. The Static and Dynamic Continuum Theory of Liquid Crystals: A Mathematical Introduction; CRC Press: Boca Raton, FL, USA, 2019.

28. Federico Toschi, M.S. Flowing Matter; Springer International Publishing: Cham, Switzerland, 2019; [DOI: https://dx.doi.org/10.1007/978-3-030-23370-9]

29. Faber, T.E. Fluid Dynamics for Physicists; Cambridge University Press: Cambridge, UK, 1995.

30. Krüger, T.; Kusumaatmaja, H.; Kuzmin, A.; Shardt, O.; Silva, G.; Viggen, E.M. The Lattice Boltzmann Method—Principles and Practice; Springer International Publishing: Cham, Switzerland, 2016.

31. Succi, S. The Lattice Boltzmann Equation: For Complex States of Flowing Matter; Numerical Mathematics and Scientific Computation Series; Oxford University Press: Oxford, UK, 2018.

32. Grossmann, C.; Roos, H.G.; Stynes, M. Numerical Treatment of Partial Differential Equations; Springer: Berlin/Heidelberg, Germany, 2007; [DOI: https://dx.doi.org/10.1007/978-3-540-71584-9]

33. Turk, J.; Svenšek, D. Backflow-mediated domain switching in nematic liquid crystals. Phys. Rev. E; 2014; 89, 032508. [DOI: https://dx.doi.org/10.1103/PhysRevE.89.032508]

34. Kirk, D.B.; Hwu, W.M.W. Programming Massively Parallel Processors: A Hands-On Approach; 1st ed. Morgan Kaufmann Publishers Inc.: San Francisco, CA, USA, 2010.

35. Guo, Z.; Zheng, C.; Shi, B. Discrete lattice effects on the forcing term in the lattice Boltzmann method. Phys. Rev. E; 2002; 65, 046308. [DOI: https://dx.doi.org/10.1103/PhysRevE.65.046308]

36. Doi, M. Soft Matter Physics; Oxford University Press: Oxford, UK, 2013; [DOI: https://dx.doi.org/10.1093/acprof:oso/9780199652952.001.0001]