Convolutional Neural Network analysis of optical texture patterns in liquid-crystal skyrmions

To effectively train any type of Machine Learning model, a sufficient amount of training data is crucial. In this study, we generate the necessary training data by minimising the Frank-Oseen elastic free energy denoted as $F$ to accurately describe the behaviour of LC under the influence of external electric fields. The system geometry and details of the numerical implementation are described in the Methods section. Under strong perpendicular boundary conditions, we obtain director profiles with so-called torons Smalyukh et al. (2010), which may be visualised as skyrmion tubes which terminate in the vicinity of the confining surfaces on point defects. This allows embedding of the skyrmion tube in a uniform background director field in three dimensions (3D). For the skyrmion tube, two-dimensional (2D) director profiles within each of the cross-sections perpendicular to the tube axis has the structure of the baby skyrmion where the director twists uniformly by 180^∘ from the tube center in all radial directions Ackerman et al. (2017).

The minimisation process renders 3D director configurations, which are then mapped to the POM images required for training of the CNNs (see Methods section for a detailed explanation). We generate single wavelength POM textures from the free energy minimising ${\bm{n}}({\bm{r}})$ by applying the standard Jones matrix method Zhao et al. (2023); Yeh and Gu (2009). The resulting optical image sums up, in an effective way, the configurational information stemming from all the skyrmions which exist within the 2D LC layers perpendicular to the direction of light propagation.

Figure 1 illustrates the impact of both LC intrinsic properties and external perturbations on the skyrmion POM textures. More specifically, Figs. 1a-d demonstrate how the texture changes by varying $\eta$ and applied voltage $U$. Figure 1e depicts the free energy of the stabilised skyrmion structure as a function of $\eta$ for two values of $U$ (see Methods for details). The observed break in slope in the blue curve at $\eta=30$ corresponds to the configurational transition where the uniform (without the torons) director profiles are stable for $\eta<30$, and non-uniform profiles with torons are observed at $\eta\geq 30$.

To predict the pitch of the chiral nematic LC we have used the POM textures around skyrmions. We have built a data set of skyrmions associated with different values of $\eta$ $\in$ {23, 24, 25, …, 40} from the numerical minimisation of the Frank-Oseen free energy. Figure 1f illustrates the network architecture used initially for this task. In this network, input images (300 × 300 pixels) pass through three blocks of two 4 × 4 convolutions and one 3 × 3 max-pooling layers, followed by three fully connected layers (with 32, 32, and 16 nodes, respectively) and an output layer. We have used rectified linear unit (ReLU) activation functions in all convolutional and fully connected layers, while the output layer uses a softmax activation function with 11 nodes. We have separated 15% of the data for a final evaluation (test set) of the model and have used the remaining 85% as a training set. The network parameters have been optimised using the Adam algorithmKingma and Ba (2014) (learning rate of 0.001), and the loss function is categorical cross-entropy (commonly used in multinomial classification). To avoid overfitting, we have applied an early stopping regularisation procedure (with patience set to 5 epochs) over all convolutional and fully connected layers. In Fig. 1g, we show how the loss function (given by the mean square error) evolves with the number of epochs. From comparing the training and the validation loss functions, we conclude that the CNN is not overfitted for this data. Figure 1h shows the confusion matrix obtained by applying the trained network to the 15% of data never exposed to the algorithm. These results demonstrate the excellent accuracy that our network can achieve in identifying the pitch of a skyrmion in the test set, where the average accuracy is $\approx 0.96$.

We now demonstrate the feasibility of using POM textures of isolated skyrmions to predict the intensity of an applied electric field. Differently from the pitch classification which is an integer number, it is now a regression problem where the network output is a continuous variable representing the applied voltage $U$. Figure 2a presents a series of POM images of a single skyrmion at a fixed $\eta$, and for increasing $U$. A transformation in the texture surrounding the skyrmion is evident, with even the size of the skyrmion exhibiting a dependence on $U$. This observation highlights the sensitivity of the skyrmion textures to external stimuli, suggesting their potential as probes for characterising such influences.

We have employed the CNN architecture depicted in Fig. 2b to predict $U$ from the POM textures similar to those in Fig. 2a. The network consists of convolutional layers to extract features, followed by fully connected layers for regression. We have increased the number of blocks of two 4 × 4 convolutions and one 3 × 3 max-pooling layers by one, reduced the number of fully connected layers by one and replaced the softmax activation function of the output layer by a linear activation function. A ReLU activation functions have been used after all convolution operations. We have trained this network by optimising the mean square error (loss function), following the same procedures used for the pitch classification. Figure 2c shows the evolution of the mean square error as the loss function during the training and testing phases, indicating that the CNN effectively learns the relationship between the POM textures and the applied voltage without overfitting. The performance of the trained CNN is illustrated in Fig. 2d, where the predicted $U$ is plotted against the true values. The close agreement with the 1:1 line and a coefficient of determination of approximately 0.94 demonstrate the high accuracy of our method in recovering the true voltage. This result evidences the fact that the CNNs can utilise effectively the information encoded within LC solitonic defects such as skyrmions to evaluate external stimuli applied to the system.

Predicting the free energy of a liquid crystal system is of importance for understanding the stability and dynamics of topological defects such as skyrmions. By accurately estimating the free energy, we can anticipate how these defects will respond to changes in experimental parameters such as temperature, applied voltage, or confinement. This capability is relevant not only for fundamental research but also for the development of future skyrmion-based devices. The ability to extract the free energy directly from POM textures using machine learning can aid the optimisation of device parameters and accelerate the design process.

We have employed the CNN architecture shown in Fig. 2e, which is similar to the one (Fig.  1f) used for classifying the skyrmion pitch. This architecture enables the effective capture of intricate features within the POM images and their relation to the system’s free energy landscape. The network comprises four blocks of two 4x4 convolutional layers with a 3x3 max-pooling layer, ReLU activation functions after all convolution operations, and a linear activation function in the output layer. We have used 80% of the data for training, 15% for validation, and 5% for final testing. The dynamics of the training process is presented in Fig. 2f, which displays the loss functions dependence on the training epochs. Both the training and testing losses decrease steadily, indicating successful learning without overfitting. In Fig. 2g we highlight the accuracy of this approach, with the predicted free energy values closely matching the true values with a coefficient of determination of $\approx 0.96$.

The ability to predict liquid crystal properties in systems with multiple skyrmions holds significant promise for advancing both fundamental research and technological applications. ”Schools” of skyrmions, exhibiting collective behaviour and intricate interactions, have garnered considerable attention due to their potential in areas such as novel displays with dynamic light control and efficient drug delivery systems Sohn et al. (2019a, 2018, b). It is essential to extend our machine learning approach to handle scenarios with multiple skyrmions to accurately characterise and exploit these complex systems.

In addition to its relevance for applications, predicting properties in multi-skyrmion systems is crucial for understanding the fundamental nature of emergent collective skyrmion dynamics governed by out-of-equilibrium interactions. The later appear due to the elastic director distortions and topological constraints, and play a key role in determining the self-assembling behaviour and stability of skyrmion ensembles. To investigate these interactions, we have generated a dataset of POM images featuring pairs of skyrmions at various relative positions and pitch values, and under different applied voltages. This dataset enables us to explore the interplay between skyrmion separation, material properties, and external stimuli.

Figures 3a-d display the POM textures of two skyrmions at fixed separation, while varying both $\eta$ and $U$. These images reveal a complex dependence of the textures on both intrinsic material properties and external influences. Figure 3e plots the free energy as a function of the relative distance between two skyrmions. Notably, the interaction behaviour undergoes a shift from repulsive at zero voltage to attractive at $U=3.5V$. This shows significant role of external fields in modulating the skyrmion interactions.

To predict the free energy of a system containing two interacting skyrmions, we have use the same CNN architecture as for the case of an isolated skyrmion free energy learning (see Fig. 2e). The efficacy of our approach is evident in the training and testing performance. Figure 3f shows that the loss function converges rapidly to a low value for both the training and test data, suggesting efficient learning and generalisation capabilities of the CNN. Furthermore, Fig.3g demonstrates the remarkable accuracy of our method, where the predicted free energy values exhibit a near-perfect 1:1 relationship with the true values. This exceptional agreement, reflected in a coefficient of determination of approximately 0.999, highlights the power of CNNs in learning to accurately predict the free energy of complex multi-skyrmion systems from their POM textures.

A possibility to combine data from both single and multiple skyrmion configurations is not only a demonstration of the generalisability and versatility of our method but also a necessity for its practical application. In realistic scenarios, LC devices may contain a varying number of skyrmions, and our model must be capable of predicting accurately properties across different configurations. Figure 4a illustrates this concept, showcasing the use of data from both one- and two-skyrmion configurations for the CNN training. This diverse dataset encompasses a range of LC properties, including the pitch, and external effects like the applied voltage, ensuring the model’s versatility and robustness.

We demonstrate the effectiveness of this approach by predicting the free energy of different skyrmion configurations, regardless of whether they involve single or multiple skyrmions. Figure 4b depicts the loss functions depending on the training epochs, revealing a fast and clear convergence to a low value. This convergence indicates successful learning and generalisation across distinct skyrmion configurations. In Fig. 4c, we present the predicted free energy values against the true values, demonstrating a very good agreement. This finding emphasises the capability of our method to accurately estimate the free energy in diverse skyrmion systems.

The versatility of our approach extends beyond the free energy prediction. Figures 4d and e highlight the model’s ability to predict both external perturbations and liquid crystal properties using data from multiple skyrmion configurations. Figure 4d shows the predicted applied voltage as a function of the true one, exhibiting a strong 1:1 agreement. This accuracy demonstrates the model’s effectiveness in capturing the influence of external stimuli on skyrmion textures. Furthermore, the confusion matrix for the classification of the LC pitch in Fig. 4e reveals a very good agreement as well, further confirming the model’s ability to extract intrinsic material properties from complex POM images. These results demonstrate collectively the potential of our approach for characterising and optimising LC systems with varying skyrmion configurations.

We use the Frank-Oseen theory of liquid crystals where we write $F$ as the free energy for chiral nematic LC in the following form:

	$\displaystyle F=\frac{1}{2}\int_{V}\biggl{[}K_{1}(\nabla\cdot{\bm{n}})^{2}+K_{2}\left({\bm{n}}\cdot\nabla\times{\bm{n}}-\frac{2\pi}{P}\right)^{2}+$
	$\displaystyle K_{3}({\bm{n}}\times\nabla\times{\bm{n}})^{2}+\varepsilon_{0}\Delta\varepsilon({\bm{E}}\cdot{\bm{n}})^{2}\biggr{]}d^{3}r,$		(1)

where ${\bm{n}}({\bm{r}})$ is the nematic director field, $K_{1},K_{2}$ and $K_{3}$ are positive elastic constants describing splay, twist and bend director distortions, respectively, and $P$ is the cholesteric pitch, and the integral is taken over the three-dimensional (3D) domain $V$ occupied by the liquid crystal. We assume $V=L\times L\times L_{z}$, where $L_{z}$ is the separation between the confining rigid surfaces of the area $L\times L$. External electric field ${\bm{E}}$ couples to the nematic director according to the last term under the integral, where $\varepsilon_{0}$ is the vacuum permittivity and $\Delta\varepsilon$ is the dielectric anisotropy. The last term in the integrand of Eq. (1) approximates the local electric field by the constant external field, which is valid for weak fields, as used in the experiments Ackerman et al. (2017); Sohn and Smalyukh (2020). In this study we assume that the electric field is applied in the direction $\hat{\bm{z}}$ perpendicular to the confining cell walls and has a form ${\bm{E}}=-U\hat{\bm{z}}$, where $U$ is constant voltage drop across the cell.

The free energy (1) is minimised by using the following relaxation dynamics of the director field:

$$\partial_{t}n_{\mu}=-\frac{1}{\gamma}\frac{\delta F}{\delta n_{\mu}}.$$

(2)

where $\gamma$ is the rotational viscosity, and $\mu=x,y,z$. Minimisers of (1) corresponds to steady state solutions of (2), which we compute numerically on a square lattice. The spatial derivatives on the r.h.s of Eq. (2) are approximated by using finite-differences and the integration over time is performed using the fourth-order Runge-Kutta method. The director field at the confining surfaces $z=0$ and $z=L_{z}$ is kept fixed perpendicular to the surfaces, and periodic boundary conditions are applied in the $x$ and $y$ directions. Equation (2) is integrated subject to the constraint $({\bm{n}}\cdot{\bm{n}})=1$. The values of the model parameters used in this study are provided in table 1.

To facilitate the formation of skyrmionic configurations we adopt as an initial condition the following axially symmetric Ansatz for the director field

	$\displaystyle n_{x}({\bm{r}})=\sin\left(a({\bm{r}})\right)\sin\left(b({\bm{r}})+\frac{\pi}{2}\right)$
	$\displaystyle n_{y}({\bm{r}})=\sin\left(a({\bm{r}})\right)\cos\left(b({\bm{r}})+\frac{\pi}{2}\right)$
	$\displaystyle n_{z}({\bm{r}})=-\cos\left(a({\bm{r}})\right),$		(3)

where (for $z^{\prime}\equiv z-C_{z}$)

$$a({\bm{r}})=\begin{cases}\frac{\pi}{2}\left[1-\tanh\left(\frac{B}{2}(r-\sqrt{R^{2}-(z^{\prime})^{2}})\right)\right],R>|z^{\prime}|\\ 0,R\leq|z^{\prime}|\end{cases}$$

(4)

	$\displaystyle b({\bm{r}})=\tan^{-1}\left(\frac{x-C_{x}}{y-C_{y}}\right)$		(5)
	$\displaystyle r=\sqrt{(x-C_{x})^{2}+(y-C_{y})^{2}}.$		(6)

$(C_{x},C_{y},C_{z})^{\mathrm{T}}$ above is the skyrmion center, $R$ set the skyrmion extention, $B$ controls the width of the twist wall of the skyrmionic tube, $r$ is the distance from a given point $(x,y,z)^{\mathrm{T}}$ to the skyrmion symmetry axis, and $b({\bm{r}})$ is the $2D$ polar angle. The values of the parameters used in the simulation are (in the simulation units defined in Table 1): $R=0.45L_{z}$, $B=0.5$, $C_{x}=C_{y}=L/2$, $C_{z}=Lz/2$.

We have implemented a CNN approach, which is a special type of neural network designed for processing data with relevant spatial patterns Goodfellow et al. (2016). Unlike traditional methods whihc require manual feature extraction, such as textures and shapes, a CNN directly processes raw pixel data from the image to extract pertinent features Simonyan and Zisserman (2015).

At the core of a CNN is the convolution operation. The network takes an input feature map, represented as a three-dimensional matrix, and extracts tiles from it. Filters are applied to these tiles to compute new features, resulting in an output feature map (convolved feature) that may differ in size and depth from the input. These filters move across the input feature map grid both horizontally and vertically, one pixel at a time, extracting the corresponding tiles. Each convolutional layer performs this operation and passes the result to the next layer. The filter, or kernel, can be adjusted by modifying hyperparameters such as its size or stride. This process is analogous to the response of a neuron in the visual cortex, where each convolutional neuron processes information only from its receptive field. While fully connected feedforward neural networks could theoretically be used to extract features and classify data of this nature, they become impractical for large inputs like high-resolution images due to the massive number of neurons required, even in a shallow architecture.

Following each convolution operation, the CNN applies an activation function to the convolved feature, introducing non-linearity into the model. This non-linearity allows the network to adapt to various data types and to differentiate between different outcomes effectively.

A pooling layer is typically applied after the convolution operation to downsample the convolved feature map, reducing its dimensionality while preserving the most relevant feature information. This layer operates similarly to the convolutional layer by sliding over the feature map, extracting tiles, and performing a simple calculation on each, such as the maximum (max pooling) or average value (average pooling). Pooling helps summarize the presence of features, making representations invariant to small translations in the input data and improving computational efficiency. Pooling layers are often inserted periodically between successive convolutional layers in a CNN architecture. The rationale behind pooling is that the precise location of a feature is less important than its general location relative to other features.

The fully connected (FC) layer functions similarly to those found in conventional neural network models, operating on a flattened input where each input is connected to all neurons. The output from the convolution and pooling layers is flattened and fed into the fully connected layer, where numerous dot product operations are performed to produce the final output. FC layers are typically found toward the end of a CNN architecture and are used to optimize objectives such as class scores.

There are numerous ways to design a convolutional neural network, involving choices related to depth, width, number and size of filters, stride, and more. These decisions are often empirical, based on the input data type and inspired by successful architectures for similar tasks. Nonetheless, some design patterns are commonly found across CNN architectures Sigaki et al. (2020). These include simplicity, symmetry, incremental feature construction, normalisation, and downsampling as the network depth increases Szegedy et al. (2015). Our design choices, though guided by these principles, were largely determined through trial-and-error and cross-validation across several network parameters.