# A Multi-Branched Radial Basis Network Approach to Predicting Complex Chaotic
Behaviours

Aarush Sinha
<EMAIL_ADDRESS>
###### Abstract

In this study, we propose a multi branched network approach to predict the
dynamics of a physics attractor characterized by intricate and chaotic
behavior. We introduce a unique neural network architecture comprised of
Radial Basis Function (RBF) layers combined with an attention mechanism
designed to effectively capture nonlinear inter-dependencies inherent in the
attractor’s temporal evolution. Our results demonstrate successful prediction
of the attractor’s trajectory across 100 predictions made using a real-world
dataset of 36,700 time-series observations encompassing approximately 28
minutes of activity. To further illustrate the performance of our proposed
technique, we provide comprehensive visualizations depicting the attractor’s
original and predicted behaviors alongside quantitative measures comparing
observed versus estimated outcomes. Overall, this work showcases the potential
of advanced machine learning algorithms in elucidating hidden structures in
complex physical systems while offering practical applications in various
domains requiring accurate short-term forecasting capabilities.

##  1 Introduction

In traditional mathematics, a radial basis function is a function that is
based on the distance between the input and a specified point, such as the
origin or a center point. A radial function is any function that meets this
property [1].
A radial function is a function $\varphi:[0,\infty)\to\mathbb{R}$. When paired
with a metric on a vector space $\|\cdot\|:V\to[0,\infty)$ a function
$\varphi_{\mathbf{c}}=\varphi(\|\mathbf{x}-\mathbf{c}\|)$ is said to be a
radial kernel centered at $\mathbf{c}$. A Radial function and the associated
radial kernels are said to be radial basis functions if, for any set of nodes
$\\{\mathbf{x}_{k}\\}_{k=1}^{n}$:

  * •

The kernels
$\varphi_{\mathbf{x}_{1}},\varphi_{\mathbf{x}_{2}},\dots,\varphi_{\mathbf{x}_{n}}$
are linearly independent (for example $\varphi(r)=r^{2}$ in $V=\mathbb{R}$ is
not a radial basis function).

  * •

The kernels
$\varphi_{\mathbf{x}_{1}},\varphi_{\mathbf{x}_{2}},\dots,\varphi_{\mathbf{x}_{n}}$
form a basis for a Haar Space, meaning that the interpolation matrix

$\begin{bmatrix}\varphi(\|\mathbf{x}_{1}-\mathbf{x}_{1}\|)&\varphi(\|\mathbf{x}_{2}-\mathbf{x}_{1}\|)&\dots&\varphi(\|\mathbf{x}_{n}-\mathbf{x}_{1}\|)\\\
\varphi(\|\mathbf{x}_{1}-\mathbf{x}_{2}\|)&\varphi(\|\mathbf{x}_{2}-\mathbf{x}_{2}\|)&\dots&\varphi(\|\mathbf{x}_{n}-\mathbf{x}_{2}\|)\\\
\vdots&\vdots&\ddots&\vdots\\\
\varphi(\|\mathbf{x}_{1}-\mathbf{x}_{n}\|)&\varphi(\|\mathbf{x}_{2}-\mathbf{x}_{n}\|)&\dots&\varphi(\|\mathbf{x}_{n}-\mathbf{x}_{n}\|)\\\
\end{bmatrix}$ (1)

is non-singular [2] [3].

Frequently utilized varieties of radial basis functions consist of:

  * •

Gaussian RBF:

$\varphi(r)=\exp\left(-\frac{r^{2}}{2\sigma^{2}}\right)$

where $r$ is the distance between the input point and the center, and $\sigma$
is a parameter controlling the width of the Gaussian.

  * •

Multiquadric RBF:

$\varphi(r)=\sqrt{1+\left(\frac{r}{\sigma}\right)^{2}}$

where $r$ is the distance between the input point and the center, and $\sigma$
is a parameter controlling the shape of the function.

  * •

Inverse Multiquadric RBF:

$\varphi(r)=\frac{1}{\sqrt{1+\left(\frac{r}{\sigma}\right)^{2}}}$

where $r$ is the distance between the input point and the center, and $\sigma$
is a parameter controlling the shape of the function.

  * •

Thin Plate Spline RBF:

$\varphi(r)=r^{2}\log(r)$

where $r$ is the distance between the input point and the center.

Imagine a ball rolling around a landscape with hills and valleys. An attractor
acts like the bottom of a valley. Regardless of where you place the ball on
the landscape (starting conditions), if it rolls downhill long enough, it will
eventually settle at the valley’s bottom (the attractor). This signifies that
the system (the ball) tends towards a specific set of values (the valley’s
position) over time. Thus, formally defining an attractor involves identifying
a group of numeric values that a system naturally gravitates towards,
irrespective of its initial parameters.

Mathematical defintion of an attractor:
Let $t$ represent time and let $f(t,\cdot)$ be a function specifying the
dynamics of the system. If $a$ is a point in an $n$-dimensional phase space,
representing the initial state of the system, then $f(0,a)=a$, and for a
positive value of $t$, $f(t,a)$ is the result of the evolution of this state
after $t$ units of time. For example, if the system describes the evolution of
a free particle in one dimension, then the phase space is the plane
$\mathbb{R}^{2}$ with coordinates $(x,v)$, where $x$ is the position of the
particle, $v$ is its velocity, $a=(x,v)$, and the evolution is given by

$f(t,(x,v))=(x+tv,v).$

An attractor is a subset $A$ of the phase space characterized by the following
three conditions:

  1. 1.

$A$ is forward invariant under $f$: if $a$ is an element of $A$, then so is
$f(t,a)$ for all $t>0$.

  2. 2.

There exists a neighborhood of $A$, called the basin of attraction for $A$ and
denoted $B(A)$, which consists of all points $b$ that "enter" $A$ in the limit
$t\to\infty$. More formally, $B(A)$ is the set of all points $b$ in the phase
space with the following property: For any open neighborhood $N$ of $A$, there
is a positive constant $T$ such that $f(t,b)\in N$ for all real $t>T$.

  3. 3.

There is no proper (non-empty) subset of $A$ having the first two properties.

Since the basin of attraction contains an open set containing $A$, every point
that is sufficiently close to $A$ is attracted to $A$. The definition of an
attractor uses a metric on the phase space, but the resulting notion usually
depends only on the topology of the phase space [4]. In the case of
$\mathbb{R}^{n}$, the Euclidean norm [5] is typically used, which is defined
as

$\|\textbf{x}\|=\sqrt{x_{1}^{2}+x_{2}^{2}+\dots+x_{n}^{2}}$

. Using these concepts we propose implementing a multi-branched radial basis
neural network to help predict the chaotic and random behaviours of an
attractor.

##  2 Related Work

Radial Basis networks have been extensively studied and proven effective in
various classification tasks [6][7]. They offer a versatile framework for
pattern recognition and data analysis, leveraging the flexibility of radial
basis functions to model complex relationships within datasets. By capturing
the intricate dynamics and nonlinear interactions inherent in real-world
phenomena, Radial Basis networks contribute to advancing our understanding of
complex systems and facilitating informed decision-making in fields ranging
from communication systems [8][9] to computational biology [10][11].
While RBF layers offer valuable capabilities in certain modeling tasks, they
alone may not be sufficient for capturing the rich dynamics and predicting
chaotic and random behaviors in attractors. To address the complexities
inherent in chaotic systems, more sophisticated and adaptable modeling
approaches are required, which may involve combining RBF layers with other
architectural components and techniques tailored to the specific
characteristics of chaotic dynamics.
Attention mechanisms [12] have emerged as powerful tools in the realm of
neural networks, offering sophisticated mechanisms for selectively focusing on
relevant parts of input data while suppressing irrelevant information.
Originally inspired by human cognitive processes, attention mechanisms have
found widespread applications in various domains, including natural language
processing, computer vision, and sequential data modeling.

##  3 Dataset

We use a pre exisitng kaggle dataset [13]. This dataset comprises time series
data originating from an unidentified physics attractor, synthesized through
undisclosed governing rules. Manifesting intricate and chaotic dynamics, the
attractor presents a challenge for analysis.
The dataset encompasses 36,700 data points, each delineating the positions of
two points in a two-dimensional space at distinct time intervals. Collected
over approximately 28 minutes, the dataset offers insights into the
attractor’s behavior over time. Notably, the system undergoes periodic resets,
typically occurring upon reentry into a recurring loop.

Variable | Type | Definition
---|---|---
time | Float | The time in seconds since the start of the simulation
distance | Float | Distance between both objects
angle1 | Angle | Angle of the first object
pos1x | Float | X position of the first object
pos1y | Float | Y position of the first object
angle2 | Angle | Angle of the second object
pos2x | Float | X position of the second object
pos2y | Float | Y position of the second object
Table 1: Description of variables present in the dataset.

##  4 Methodology

It defines the network architecture consisting of several components:

  1. 1.

Branches: Three separate branches are utilized, each focusing on learning the
relationship between a specific pair of input columns.

     * •

An RBFLayer: Performs the Radial Basis Function transformation on the input
data; customizable parameters include the number of kernels $(K)$, output
features $(F_{\text{o}})$, radial function $(\varphi)$, norm function
$\|\cdot\|$, and normalization option (normalize). We use the
inverse_multiquadric function • ‣ 1 and Euclidean norm 1.

     * •

Dropout layer: Introduced with a probability of $0.3$ to mitigate overfitting.

     * •

AttentionLayer: Focuses on significant portions of the transformed data within
the branch.

     * •

Linear layers with $\operatorname{ReLU}(x)=\max\\{0,x\\}$ and
$\operatorname{tanh}(x)=\displaystyle\frac{e^{x}-e^{-x}}{e^{x}+e^{-x}}$
activation functions for additional feature extraction and transformation.

  2. 2.

Merging Layer: Following the processing of each pair of columns within their
respective branches, the outputs are concatenated. A linear layer with a ReLU
activation function integrates the combined information.

  3. 3.

Output Layer: A final linear layer with an output size of $3$ projects the
merged features onto the desired three-dimensional prediction.

Denote $x\in\mathbb{R}^{3}$ as the input vector having three features. Let
${\hat{y}\in\mathbb{R}^{3}}$ denote the output of the model. Each branch
accepts a pair of input features represented as $x_{i},x_{j}\in\\{1,2,3\\}$,
satisfying ${i\neq j}$.

The forward function governs the data flow through the network:

  1. 1.

Input Splitting: Separation of the input data $x$ into three distinct columns,
representing the features: $x=(x_{1},x_{2},x_{3})$.

  2. 2.

Branch Processing: Feeding each pair of columns into the assigned branch
(branch1, branch2, or branch3); subsequently processed through their
constituent layers, yielding an output per pair.

  3. 3.

Output Concatenation: The individual branch outputs, namely
$(out_{1},out_{2},out_{3})$, undergo concatenation along the feature
dimension.

  4. 4.

Merging: Transmission of the concatenated outputs through the merging layer
produces a unified representation.

  5. 5.

Prediction: Applying the merged features to the ultimate output layer yields
the three-dimensional prediction $(\hat{y}_{1},\hat{y}_{2},\hat{y}_{3})$.

Figure 1: Proposed multi layer architecture Figure 2: Single Sequential
Proposed Layer

This design enables the model to discern specific relationships amongst
diverse input feature pairs while combining the learned features through the
attention mechanism and merging stages for delivering the final prediction.

##  5 Training

We train the model on a singular NVIDIA A30 GPU. It takes 2 hours to train the
model for 2000 epochs with a batch_size=512. We use the Mean Squared Error
(MSE) loss function as our criterion:

$\mathrm{MSE}(\hat{\boldsymbol{y}},\boldsymbol{y})=\frac{1}{N}\sum_{i=1}^{N}\left(\hat{y}_{i}-y_{i}\right)^{2}$

Where $\hat{\boldsymbol{y}}$ represents the predicted values, $\boldsymbol{y}$
represents the actual target values, and $N$ is the total number of samples.
The MSE computes the average of the squared differences between predicted and
actual values, providing a measure of the model’s performance in minimizing
prediction errors. We utilize Adam [14] optimizer for our model:

$\displaystyle m_{t}$ $\displaystyle=\beta_{1}\cdot m_{t-1}+(1-\beta_{1})\cdot
g_{t},$ $\displaystyle v_{t}$ $\displaystyle=\beta_{2}\cdot
v_{t-1}+(1-\beta_{2})\cdot g_{t}^{2},$ $\displaystyle\hat{m}_{t}$
$\displaystyle=\frac{m_{t}}{1-\beta_{1}^{t}},$ $\displaystyle\hat{v}_{t}$
$\displaystyle=\frac{v_{t}}{1-\beta_{2}^{t}},$ $\displaystyle\theta_{t+1}$
$\displaystyle=\theta_{t}-\frac{\eta}{\sqrt{\hat{v}_{t}}+\epsilon}\cdot\hat{m}_{t},$

Where $m_{t}$ and $v_{t}$ are the first and second moment estimates, $g_{t}$
is the gradient, $\beta_{1}$ and $\beta_{2}$ are the exponential decay rates
for the moment estimates, $\hat{m}_{t}$ and $\hat{v}_{t}$ are bias-corrected
estimates, $\theta_{t}$ is the parameter at iteration $t$, $\eta$ is the
learning rate, and $\epsilon$ is a small constant to prevent division by zero.

Finally for comparison purposes we train with the same hyper parameters two
models a singular branch 2 and the proposed model 1. All implementations were
done in PyTorch[15].

##  6 Results

Figure 3: Loss over iterations of the Single Sequential Network Figure 4: Loss
over iterations of the Multi-Branched Network

  * •

Loss over iterations of the Single Sequential Network 3

    * –

The training loss for Object 1 (blue) starts high, sharply decreases, and then
fluctuates around a lower level with some spikes.

    * –

The training loss for Object 2 (orange) follows a similar pattern but
maintains a higher overall loss throughout the training process.

    * –

There are large spikes in the loss for both objects early in training,
indicating potential instability or difficulty in the initial learning phase.

    * –

The loss seems to stabilize and flatten out more towards the end of the
training iterations shown.

  * •

Loss over iterations of the Multi-Branched Network 4

    * –

The overall pattern is similar to Image 1, with Object 2’s loss (orange) being
consistently higher than Object 1’s loss (blue).

    * –

However, the initial large spikes in loss are more prominent and last longer
compared to Image 1.

    * –

The loss curves appear to flatten out and stabilize at a later point in the
training process compared to Image 1.

    * –

There are fewer small fluctuations and spikes in the loss curves once they
stabilize, suggesting potentially smoother convergence.

In summary, while the overall trend of Object 2 having higher training loss is
consistent across both images, the single sequential network exhibits more
pronounced initial instability and takes longer to stabilize compared to the
multi-branched architecture.

The results were saved in the form of GIFs. To view those GIFs kindly go to
this Google Drive Folder

##  7 Conclusion

In conclusion, this paper has explored the application of Radial Basis
Function Neural Networks (RBFNNs) in predicting chaotic and random behaviors
in attractors. Through a comprehensive review of related work, we have
highlighted the strengths and limitations of RBFNNs in capturing the complex
dynamics of chaotic systems. Leveraging insights from chaos theory and neural
network architecture, we have proposed novel approaches for enhancing the
predictive capabilities of RBFNNs with attention mechanisms. Furthermore, we
have presented a detailed description of our experimental setup, including the
dataset used, model architecture, and training procedure. Our results
demonstrate the effectiveness of our proposed methods in predicting chaotic
and random behaviors in attractors, surpassing previous approaches in terms of
predictive accuracy and computational efficiency. Overall, this paper
contributes to advancing our understanding of chaotic systems and lays the
groundwork for future research in utilizing RBFNNs for predictive modeling in
complex dynamical systems.

##  8 Reproducibility

Results can be reproduced from the code present in my GitHub repository.

##  9 Acknowledgement

We acknowledge the work of Alessio Russo who originally implemented RBFNNs in
PyTorch. His work is available on his GitHub[16].

## References

  * [1] Contributors to Wikimedia projects.  Radial basis function - Wikipedia, 2024.
  * [2] Gregory E. Fasshauer.  Meshfree Approximation Methods with MATLAB.  World Scientific Publishing Co. Pte. Ltd., Singapore, 2007.
  * [3] Holger Wendland.  Scattered Data Approximation.  Cambridge University Press, Cambridge, 2005.
  * [4] John Milnor.  On the concept of attractor.  Communications in Mathematical Physics, 99(2):177–195, Jun 1985.
  * [5] M. Emre Celebi, Fatih Celiker, and Hassan A. Kingravi.  On euclidean norm approximations, 2010.
  * [6] Yue Wu, Hui Wang, Biaobiao Zhang, and K.-L. Du.  Using Radial Basis Function Networks for Function Approximation and Classification.  International Scholarly Research Notices, 2012, March 2012.
  * [7] James A Leonard and Mark A Kramer.  Radial basis function networks for classifying process faults.  IEEE Control Systems Magazine, 11(3):31–38, 1991.
  * [8] Deng Jianping, Narasimhan Sundararajan, and P Saratchandran.  Communication channel equalization using complex-valued minimal radial basis function neural networks.  IEEE Transactions on neural networks, 13(3):687–696, 2002.
  * [9] Hao Yu, Tiantian Xie, Stanisław Paszczynski, and Bogdan M Wilamowski.  Advantages of radial basis function networks for dynamic system design.  IEEE Transactions on Industrial Electronics, 58(12):5438–5450, 2011.
  * [10] A Vande Wouwer, Christine Renotte, and Ph Bogaerts.  Biological reaction modeling using radial basis function networks.  Computers & chemical engineering, 28(11):2157–2164, 2004.
  * [11] Yu-Yen Ou et al.  Identifying the molecular functions of electron transport proteins using radial basis function networks and biochemical properties.  Journal of Molecular Graphics and Modelling, 73:166–178, 2017.
  * [12] Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N. Gomez, Lukasz Kaiser, and Illia Polosukhin.  Attention is all you need, 2023.
  * [13] NIKITRICKY.  Physics attractor time series dataset, 2023.
  * [14] Diederik P. Kingma and Jimmy Ba.  Adam: A method for stochastic optimization, 2017.
  * [15] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Köpf, Edward Yang, Zach DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala.  Pytorch: An imperative style, high-performance deep learning library, 2019.
  * [16] Alessio Russo.  Pytorch rbf layer, 2021.


# Characterizing Quantum Evolutions via a Recommender System

Priya Batra<EMAIL_ADDRESS>Anukriti Singh
<EMAIL_ADDRESS>T. S. Mahesh<EMAIL_ADDRESS>Department of
Physics and NMR Research Center,
Indian Institute of Science Education and Research, Pune 411008, India

###### Abstract

We explore the possibility of predicting the properties of quantum evolutions
via matrix factorization algorithm, a particular type of the recommender
system (RS). A system undergoing a quantum evolution can be characterized in
several ways. Here we choose (i) quantum correlations characterized by
measures such as entropy, negativity, or discord, and (ii) state-fidelity.
Using quantum registers with up to 10 qubits, we demonstrate that an RS can
efficiently characterize both unitary and nonunitary evolutions. After
carrying out a detailed performance-analysis of the RS in two-qubits, we show
that it can be used to distinguish a clean database of quantum correlations
from a noisy one. Moreover, we find that the RS may bring about a significant
computational advantage for estimating the quantum discord of a general
quantum state, for which no simple closed-form expression exists. Finally, we
show that the RS can efficiently characterize systems undergoing nonunitary
evolutions in terms of both reduced quantum correlations and state-fidelity.

Machine learning, recommender system, matrix factorization algorithm,
artificial intelligence, nonclassical correlations

##  I Introduction

Machine learning is increasingly featuring in almost every aspect of our
understanding of the world. A popular class of machine learning, namely
Recommender system (RS), is often used to orient consumers towards certain
products according to individual preferences Linden  _et al._ (2003); Funk
(2006); Schafer  _et al._ (2001). There are many approaches to build an RS,
such as the content based system, where recommendations are based on users
past experience Lieberman (1997); Maes (1995), knowledge based system, where
one uses the knowledge of users and items to render the recommendations Burke
(2000), and the widely used collaborative filtering, which exploits the user-
item correlation, i.e., the interconnection between an user’s preferences
among products with the recommendations provided by other users Goldberg  _et
al._ (1992); Schafer  _et al._ (2007).

Collaborative filtering is further classified as neighborhood methods and
latent-factor modeling. Matrix factorization algorithm (MFA) is a popular tool
for an RS via latent-factor modeling. It can predict rating of a specific item
by particular user based on the self-rating of other items as well as others’
ratings of various items. Koren  _et al._ (2009).

Recently there has been a substantial progress in using machine learning for
quantum information problems. Examples include applications in quantum
tomography Torlai  _et al._ (2018); Kieferová and Wiebe (2017), quantum error
correction Torlai and Melko (2017), quantum control August and Ni (2017),
understanding quantum phase transitions Carrasquilla and Melko (2017);
Canabarro  _et al._ (2019), and studying quantum-many-body problems Carleo and
Troyer (2017); Gao and Duan (2017). It has been shown recently that machine
learning techniques can work as state classifiers too. Sirui Lu and co-workers
have showed the separability criteria of entangled state using convex hull
approximation and supervised learning Lu  _et al._ (2018). Y. C. Ma and M. H.
Yung showed that it is possible to classify the separable and entangled states
using artificial neural networks Ma and Yung (2018). The work has been further
extended to experimental data Gao  _et al._ (2018) and later has been applied
to simultaneous learning of multiple nonclassical correlations as well Yang
_et al._ (2019). Valeria Cimini and co-workers proposed an artificial neural
network to calculate the negativity of Wigner function for multi-mode quantum
states Cimini  _et al._ (2020). There has also been a recent work to find the
connection of geometric and entropy discord using machine learning Li  _et
al._ (2019).

In this work, we first carry out a detailed performance study of an RS on a
two-qubit system. By studying the dependence of prediction efficiency on
systematically introduced noise in the input database, we find that an RS can
identify a noisy database or distinguish a fake database from a genuine one.
We then compare the computational efficiency of the RS prediction of discord
in two and three-qubit systems with the conventional method involving
subtraction of the maximum possible classical correlation from the total
correlation. We find that, within certain tolerance limits, an RS can be
significantly faster. We also demonstrate the feasibility of the RS
application in larger systems, by predicting changes in entropy, changes in
discord, and fidelity of quantum states in registers with up to 10 qubits.
Finally, we demonstrate the RS ratings of nonunitary evolutions by
successfully predicting discord changes and fidelity values of a two-qubit
system subjected to independent single-qubit decoherence channels.

The paper is organized as follows. We briefly introduce the RS via MFA in Sec.
II and then in Sec. III we describe adapting it to characterize unitary
quantum evolutions in terms of changes in quantum-correlation. In Sec. IV, we
describe using an RS to characterize unitary evolutions in terms of state-
fidelity. Finally, after describing the RS characterization of nonunitary
evolutions in Sec. V, we conclude in Sec. VI.

##  II Recommender System via Matrix Factorization Algorithm

MFA represents the user-item interaction in a lower dimensional latent space
Wikipedia contributors (2020). User-item interaction matrix can be decomposed
in two lower dimension matrices with a latent-factor dimensionality $f$.
Consider a set of $m$ users and a set of $n$ items. Each user $i$ is
represented using a parameter vector $\Theta^{(i)}\in\mathbb{R}^{f}$ and each
item $j$ is represented using a feature vector $X^{(j)}\in\mathbb{R}^{f}$ (see
Fig. 1(a)). Here $\mathbb{R}^{f}$ is the coordinate space of dimension $f$
over real numbers. The interaction between a user $i$ and an item $j$ is
modeled by the scalar product

$r_{i,j}=\Theta^{(i)}\cdot X^{(j)}=\sum_{l=1}^{f}\Theta^{(i)}_{l}X^{(j)}_{l}$
(1)

and can be conceived as the predicted rating. Thus the task reduces to
finding, for all users and all items, the vectors $\Theta$ and $X$, which are
consistent with known ratings and thereby make the best predictions about the
unknown elements.

Initially one starts with random parameter as well as feature vectors and
evaluate rating elements $r_{i,j}$. Let $\kappa=\\{(i,j)\\}$ be the set of
user-item pairs for which the ratings ${\cal R}_{i,j}$ are known. The mismatch
between the actual and evaluated ratings is given by

$J_{0}=\sum_{(i,j)\in\kappa}(r_{i,j}-{\cal R}_{i,j})^{2}.$ (2)

The next step is to optimize $\Theta$ and $X$ vectors by minimizing $J_{0}$.
To avoid over-fitting, two regularization terms are added, such that the
objective function is of the form

$J=\frac{J_{0}}{2}+\frac{\lambda}{2}\sum_{i=1}^{m}\|\Theta^{(i)}\|+\frac{\lambda}{2}\sum_{j=1}^{n}\|X^{(j)}\|,$
(3)

where $\lambda$ is the regularization parameter and $\|\cdot\|$ denotes the
norm of the vector. Here we use first order gradient descent algorithm for the
minimization task. The gradients in the $k$th iteration can be casted as

$\displaystyle G_{\Theta^{(i,k)}}=\frac{\partial
J^{(k)}}{\partial\Theta^{(i,k)}}=\sum_{j\in\kappa_{(i,.)}}(r_{i,j}^{(k)}-{\cal
R}_{i,j})X^{(j,k)}+\lambda\Theta^{(i,k)},$ $\displaystyle
G_{X^{(j,k)}}=\frac{\partial J^{(k)}}{\partial
X^{(j,k)}}=\sum_{i\in\kappa_{(.,j)}}(r_{i,j}^{(k)}-{\cal
R}_{i,j})\Theta^{(i,k)}+\lambda X^{(j,k)},~{}~{}~{}$ (4)

where $\kappa_{(i,.)},\kappa_{(.,j)}\subset\kappa$ with the respective indeces
being fixed. Once the gradients are calculated, the parameter and feature
vectors for $(k+1)$th iteration are updated according to

$\displaystyle\Theta^{(i,k+1)}$ $\displaystyle:=$
$\displaystyle\Theta^{(i,k)}-\alpha G_{\Theta^{(i,k)}}$ $\displaystyle
X^{(j,k+1)}$ $\displaystyle:=$ $\displaystyle X^{(j,k)}-\alpha G_{X^{(j,k)}},$
(5)

where $\alpha$ is the suitable step size. After reaching a desired number $K$
of iterations, or desired value of the objective function, one can calculate
the final rating $r_{i,j}^{(K)}=\Theta^{(i,K)}\cdot X^{(j,K)}$ for all unknown
elements $(i,j)\notin\kappa$ Wikipedia contributors (2020); Andrew Ng,
Stanford Machine Learning Course .

##  III Characterizing Unitary Evolutions via correlation changes

###  III.1 Adopting MFA recommender system:

Here the RS is based on an input database containing user recommendations of
various items. It is insightful to consider the example of machine
recommendations depending on viewers’ ratings of movies (Fig. 1a). Every
viewer rates a few of the movies, which also leaves an imprint of the viewer’s
tastes or preferences. Also considering the recommendations provided by other
viewers, the RS then predicts ratings of a particular viewer for movies not
yet seen/rated by the viewer. In our RS, quantum states are viewers and
quantum evolutions are movies (Fig. 1b). Quantum evolutions transform the
states and thereby possibly change their internal quantum correlations. After
undergoing a quantum evolution, a state can rate the evolution according to
the change in its quantum correlation. The task now is to predict the unknown
ratings in the database (Fig. 1b).

Figure 1: (a) Movie database based on viewers’ rating. The RS assigns an
$f$-dimensional parameter vector $\Theta^{(i)}$ for each viewer $i$ and a
$f$-dimensional feature vector $X^{(j)}$ for each movie $j$ as illustrated.
(b) An analogous database of change $\Delta C_{i,j}$ of quantum correlations
in states $\rho_{i}$ introduced by quantum evolutions $U_{j}$.

In our analysis, we first setup a rating database ${\cal R}$ with a set
$\\{U_{1},U_{2},\cdots,U_{n_{u}}\\}$ of $n_{u}$ randomly generated unitary
operators and a set $\\{\rho_{1},\rho_{2},\cdots,\rho_{n_{s}}\\}$ of randomly
generated $n_{s}$ quantum states. A random unitary operator is generated by
matrix exponentiation of a random anti-Hermitian generator. Random pure states
are generated simply by normalizing a random vector of complex elements. We
use Bures method Al Osipov  _et al._ (2010); Maziero (2015) to generate a
random mixed state

$\rho=A_{U}A_{U}^{\dagger}/\mbox{Tr}[RR^{\dagger}]~{}\mbox{with}~{}R=(\mathbbm{1}+U)A.$

Here $\mathbbm{1}$ is the identity matrix, $U$ is a random unitary operator
and $A$ is a random complex matrix.

An element of the database matrix ${\cal R}$, corresponding to $i$th state and
$j$th unitary operator is the change in a measure $C$ of quantum correlation,

${\cal R}_{i,j}=\Delta C_{i,j}=C(U_{j}\rho_{i}U_{j}^{\dagger})-C(\rho_{i}).$
(6)

We consider the following three measures of quantum correlation:

(i) Entropy ($S$): The Von Neumann entropy

$S(\rho^{A}_{i})=-\mbox{Tr}[\rho^{A}_{i}\log\rho^{A}_{i}]=-\mbox{Tr}[\rho^{B}_{i}\log\rho^{B}_{i}]=S(\rho^{B}_{i})$
(7)

is considered as the entropy of entanglement for a pure state
$\rho^{AB}_{i}=|{\psi^{AB}_{i}}\rangle\langle{\psi^{AB}_{i}}|$ Plenio and
Virmani (2014). Here $\rho^{A}_{i}=\mbox{Tr}_{B}[\rho^{AB}_{i}]$ and
$\rho^{B}_{i}=\mbox{Tr}_{A}[\rho^{AB}_{i}]$ are the reduced density operators
obtained after tracing out qubit $B$ or qubit $A$ respectively. The entropy
change $\Delta S_{i,j}=S(U_{j}\rho^{A}_{i}U^{\dagger}_{j})-S(\rho^{A}_{i})$ is
listed in the database as entropy rating.

(ii) Negativity ($N$): For a general quantum state $\rho^{AB}_{i}$, pure or
mixed, negativity is a convex entanglement monotone that is defined as Plenio
and Virmani (2014)

$N(\rho^{AB}_{i})=\frac{\|R_{A}\|-1}{2},$ (8)

where $R_{A}$ is the partial transpose of $\rho^{AB}$ with respect to the
subsystem $A$, and $\|R_{A}\|=\mbox{Tr}\sqrt{R_{A}^{\dagger}R_{A}}$ is the
trace norm of $R_{A}$. The change in negativity brought about by a unitary
operator $U_{j}$ is $\Delta
N_{i,j}=N(U_{j}\rho^{AB}_{i}U_{j}^{\dagger})-N(\rho^{AB}_{i})$. Since the two-
qubit negativity is bounded between 0 and 1/2, it is convenient to compare
twice of negativity with other measures.

(iii) Quantum discord ($D$): It is based on the notion that if one subtracts
maximum possible classical correlation from the total correlation, the residue
must be from quantum correlation. Discord is defined in terms of difference
between two different ways of estimating mutual information Ollivier and Zurek
(2001); Katiyar  _et al._ (2012):
$I(\rho^{AB}_{i})=S(\rho^{A}_{i})+S(\rho^{B}_{i})-S(\rho^{AB}_{i})$,
$J(\rho^{AB}_{i})=S(\rho^{B}_{i})-S(\rho^{B}_{i}|\rho^{A}_{i})$, and

$D(\rho^{AB}_{i})=I(\rho^{AB}_{i})-\mbox{max}_{\\{\Pi^{A}\\}}J(\rho^{AB}_{i}),$
(9)

wherein the maximization is carried out over all possible measurement bases
$\\{\Pi^{A}\\}$ in subsystem A.

Unlike negativity, discord may be nonzero even for certain separable states,
which therefore qualify to be labeled as quantum correlated. Although this
fact supports discord as a stronger measure of quantum correlation, unlike
entropy and negativity, there is no simple analytical expression for
estimating discord. Moreover, the maximization of mutual information over all
possible measurement bases becomes prohibitively expensive for larger system
sizes. In this regard, it will be interesting to see how machine learning
works for predicting the change in quantum discord $\Delta
D_{i,j}=D(U_{j}\rho^{AB}_{i}U^{\dagger}_{j})-D(\rho^{AB}_{i})$. If an input
state is initially uncorrelated or of zero discord, then the RS is essentially
estimating the discord content in the transformed state itself.

###  III.2 Characterizing two-qubit unitaries

We first consider a million-element rating database with $n_{u}=1,000$
randomly generated two-qubit unitaries and $n_{s}=1,000$ randomly generated
two-qubit quantum states, and setup the rating matrix ${\cal R}$ with
$n_{s}\times n_{u}$ elements. As described before, $\Delta S$ ratings are for
a database with random pure states, while $\Delta N$ and $\Delta D$ values are
both for a common database with random pure/mixed states. From the complete
database, we randomly remove a set $\bar{\kappa}=\\{(i,j)\\}$ of $n_{r}$
elements, which are to be predicted by the RS. The prediction error is
measured by the root-mean-square (RMS) deviation

$\delta_{C}=\left[\sum_{(i,j)\in\bar{\kappa}}(\Delta C_{0_{i,j}}-\Delta
C_{i,j})^{2}\right]^{1/2},$ (10)

between the actual correlation changes $\\{\Delta C_{0_{i,j}}\\}$ and their
predicted values $\\{\Delta C_{i,j}\\}$.

Figure 2: Red dots represent the RS predictions of changes in two-qubit
correlations: entropy change $\Delta S$ (a), negativity change $2\Delta N$
(b), and discord change $\Delta D$ (c), plotted against the actual values
($\Delta S_{0}$, $2\Delta N_{0}$, $\Delta D_{0}$) and the number $n_{r}$ of
unknown elements in a database of 1000 states and 1000 unitary operators. Blue
lines represent the ideal case $\Delta C_{0}=\Delta C$. RMS deviations
$\delta_{C}$ between predicted and actual ratings are also indicated.

Results of the RS predictions for a two-qubit register are displayed in Fig.
2. Here, the predicted values $\Delta C$ are plotted versus the actual values
$\Delta C_{0}$ for various numbers $n_{r}$ of unknown ratings. Entropy ratings
(Fig. 2(a)) are for pure states, while the negativity and discord ratings
(Fig. 2(b)) are for mixed states. It is clear that the RS is quite successful
in predicting the changes in all the correlation measures. The RMS deviations
remained less than 0.05, except for $n_{r}=9\times 10^{5}$ which corresponds
to 90% the million-element database being unknown. The systematic dependence
of RMS deviations on $n_{r}$ is expected, since as the database becomes more
and more sparse, the minimum in the latent parameter space gets shallower, and
harder it gets to predict the ratings.

###  III.3 Dependence on the dimensions of database

We now fix the number of unknown ratings $n_{r}=100$, and vary the dimension
of the database. We constructed a set of databases with number of states
$n_{s}$ varying from 10 to 10,000, and the number of unitaries $n_{u}$ varying
from 10 to 1000. Fig. 3(a-c) display the RMS deviations between the predicted
and actual values of entropy, negativity, and discord values respectively. As
expected, for a fixed number $n_{r}$ of unknown ratings, the RS efficiency in
terms of RMS deviations between the actual and predicted correlations improves
as the size of the database grows. We observe an interesting diagonal symmetry
indicating better efficiency with roughly equal number of states and
evolutions. It is also interesting to note that while the entropy rating is
the most sensitive to the database size, the discord rating is least
sensitive.

Figure 3: RMS deviations between the predicted and actual changes of entropy
(a), negativity (b), and discord (c), for various sizes of the database, but
with fixed number $n_{r}=100$ of unknown ratings.

###  III.4 Rating the Werner state

To observe the RS in action with a concrete example, we replace one of the
random evolutions with the two-qubit controlled-NOT operation

$U_{\mbox{CNOT}}=|{0}\rangle\langle{0}|\otimes\mathbbm{1}_{2}+|{1}\rangle\langle{1}|\otimes\sigma_{x},$
(11)

where $\mathbbm{1}_{2}$ is the $2\times 2$ identity operation and $\sigma_{x}$
is the NOT gate, i.e., the Pauli x-operator. Further we choose one of the
input states to be the separable state,

$\rho_{-}=(1-\epsilon)\mathbbm{1}_{4}/4+\epsilon|{-}\rangle\langle{-}|\otimes|{1}\rangle\langle{1}|,$
(12)

where $|{-}\rangle=(|{0}\rangle-|{1}\rangle)/\sqrt{2}$ is the single-qubit
superposition state and the scalar quantity $\epsilon$ is the purity of state.
While the state $\rho_{-}$ is a pure state for $\epsilon=1$ and becomes
completely mixed for $\epsilon=0$, it is separable for all values of
$\epsilon$, and hence has zero negativity and discord. Upon acted by the CNOT
operation, $\rho_{-}$ transforms into the Werner state

$\rho_{W}=(1-\epsilon)\mathbbm{1}_{4}/4+\epsilon|{S_{0}}\rangle\langle{S_{0}}|,$
(13)

which is the convex sum of the maximally mixed state and the singlet state
$S_{0}=(|{01}\rangle-|{10}\rangle)/\sqrt{2}$. It is known that Werner state is
quantum correlated and thus has nonvanishing discord for all values
$\epsilon>0$ Katiyar  _et al._ (2012). However, it is entangled and hence has
nonzero negativity only for $\epsilon>1/3$. Since the initial state $\rho_{-}$
is quantum uncorrelated, the change in negativity or negativity of upon
transforming into $\rho_{W}$ by the CNOT operation correspond to the amount of
discord or negativity in $\rho_{W}$ itself. Thus it is interesting to see the
RS prediction of negativity and discord in this particular case. Fig. 4
plotting the predicted values $\Delta N$ and $\Delta D$ versus the purity
factor $\epsilon$ shows an excellent agreement between actual and predicted
values. RMS deviations also remained less than 0.05 indicating the high
quality of predictions.

Figure 4: Symbols indicate the predicted values of $\Delta N$ and $\Delta D$
values for the quantum uncorrelated input state $\rho_{-}$ transformed into
the Werner state $\rho_{W}$ by the CNOT operation. Solid lines represent
actual values and symbols represent predicted values as indicated by the
legend box. The RMS deviations are $\delta_{2N}=0.026$ and $\delta_{D}=0.038$.

###  III.5 Identifying noisy database

Since the RS takes an input database, it is interesting to ask the following
questions. What would be the dependence of rating efficiency, if the input
database itself is noisy, unreliable, or even fake? Can we use this dependency
for asserting if a database, with or without a large number of unknown
elements, is clean or not? In order to answer these questions, we choose
discord rating in a two-qubit register. We first construct a $1000\times 1000$
noise-free discord database ${\cal R}^{(0)}$, from which we setup a set of
databases ${\cal R}^{(\eta)}$, with varying noise parameter $0<\eta<1$. The
noisy entries are of the form

${\cal R}_{i,j}^{(\eta)}=\eta s_{ij}+(1-\eta){\cal R}_{i,j}^{(0)},$ (14)

where $s_{ij}$ is a random number such that $-1<s_{ij}<1$. Thus, as $\eta$
discretely varies from 0 to 1, the corresponding database goes from a clean
discord database to one with completely random entries. In each database,
corresponding to a fixed value of $\eta$, we randomly remove 500,000 entries,
which constitutes 50% unknown elements. We now attempt to complete the ratings
of the unknown elements starting from each of these databases. The results
shown in Fig. 5 indicates that for small noise parameter $\eta\leq 0.01$,
there is no significant effect on the RMS error. However, as the noise builds
up, the RS finds it increasingly difficult to rate the entries, and the RMS
error increases by about five times for 10% noise. Finally, as expected, the
RS completely fails for the completely random database corresponding to
$\eta=1$. Thus, if a reference database is provided, the RS can help identity
if a given database is noisy or fake.

Figure 5: Rating two-qubit quantum discord based on a noise-free database (a)
or noisy databases (b-d) with noise parameter $\eta$ as shown in each case. In
subplot (e), the RMS error $\delta_{D}$ is plotted versus the noise parameter
$\eta$. In all the cases, the number of unknown elements $n_{r}=500,000$,
which corresponds to 50% of the 1000$\times$1000 database. Figure 6: Rating
quantum discord for two- (a,b) and three-qubit (c,d) systems, with 50% (a,c)
and 75% (b,d) of unknown elements. The RMS error and ratio
$\tau_{\mathrm{cal}}/\tau_{\mathrm{MR}}$ of computational time for standard
calculation to that of the RS are shown in each case. Figure 7: Predicted
changes in quantum correlations (vertical axes) versus actual values
(horizontal axes). The columns correspond to various sizes of registers as
indicated by the number $n_{q}$ of qubits. Top and bottom rows correspond to
change in entropy $\Delta S$ and change in discord $\Delta D$ respectively.
The scaled RMS deviation $\delta^{\prime}$ between the predicted and actual
ratings is shown in each case.

###  III.6 Computational Time

Since a simple analytical expression for calculating quantum discord does not
exist, it is interesting to compare the computational time
$\tau_{\mathrm{RS}}$ for the RS prediction of quantum discord with the
calculation time $\tau_{\mathrm{cal}}$ by the standard approach that involves
subtracting maximum possible classical correlation from total correlation
Ollivier and Zurek (2001); Katiyar  _et al._ (2012). As described earlier, the
maximization is carried out over all possible measurement bases. Although
there exists, for a two-qubit system, an optimal set of measurements that
substantially reduces computational time Lu  _et al._ (2011), no such
shortcuts are known for larger systems. In the case of three-qubit system
$ABC$, we measured qubit $A$ along 40 projectors uniformly distributed over
the single-qubit Bloch sphere, and evaluated $D(\rho_{\rho_{A|BC}})$. Fig. 6
displays the discord ratings of two- and three-qubit systems. In the case of
two-qubits, the ratio $\tau_{\mathrm{cal}}/\tau_{\mathrm{RS}}$ was 6 and 10
respectively for 50% as well as 75% unknown elements. Similarly, in the case
of three-qubits, the ratio $\tau_{\mathrm{cal}}/\tau_{\mathrm{RS}}$ was 83 and
123 again for 50% as well as 75% unknown elements respectively. The RMS error
is also displayed in each case. Thus, depending the desired tolerance levels,
one can choose to adopt an RS for faster evaluation of the discord database
starting from a sparse database. The speedups are up to one and two orders of
magnitude in two- and three-qubit systems respectively.

###  III.7 Discord rating of multi-qubit unitaries

We now study the efficiency of predicting the correlation changes in larger
quantum registers, with sizes up to 10 qubits. For each register, we
constructed a database with 100 random unitary operators and 100 random
quantum states. The number of unknown ratings was fixed at $n_{r}=100$. In
larger qubit systems, there are multiple ways of partitioning the system for
estimating quantum correlations. For the sake of numerical simplicity, we
first partition the whole system of $n_{q}$ qubits into two parts: a two-qubit
part $P$ and $n_{q}-2$ qubit part $Q$. After calculating the reduced density
operator $\rho^{P}=\mbox{Tr}_{Q}[\rho^{PQ}]$, we estimated the quantum
correlation in the two-qubit part $P$ as described in the previous section. In
Fig. 7, the predicted values of changes in quantum correlations $\Delta S$
(top row) or $\Delta D$ (bottom row) are plotted versus actual values $\Delta
S_{0}$ (top row) or $\Delta D_{0}$ (bottom row) for various sizes of quantum
registers with upto 10 qubits. As expected, the absolute range of correlation
changes decreased steadily with the number of qubits. This effect was more
severe in the case of negativity, which mostly vanished in registers with more
than three qubits, and hence we now omit it from the discussions. Furthermore,
to compare the RMS errors with respective ranges, we rescale it as
$\delta_{C}^{\prime}=\delta_{C}/\sqrt{m}$, where $m=\max\\{\Delta
C_{0_{i,j}},\Delta C_{i,j}\\}$ where maximum is taken over all the $n_{r}$
rated elements. Despite the smaller size of database, the predictions were
largely in agreement with the actual values, and the scaled RMS errors mostly
remained below 0.3, thus confirming the feasibility of the RS predictions for
various system sizes.

##  IV Characterizing Unitary Evolutions via State-Fidelity

Instead of rating quantum correlations, one can also rate the state-fidelity.
For a set of pure input states $\\{|{\psi_{i}}\rangle\\}$ and unitaries
$\\{U_{i}\\}$, the database elements

${\cal
R}_{i,j}=F_{i,j}(\psi_{i},U_{j})=|\langle\psi_{i}|U_{j}|\psi_{i}\rangle|^{2}$
(15)

describe the fidelity of output state $U_{j}|{\psi_{i}}\rangle$ with the input
state $|{\psi_{i}}\rangle$. Similarly, for a set of mixed input states
$\\{\rho_{i}\\}$ we may use the Uhlmann trace distance

${\cal
R}_{i,j}=F_{i,j}(\rho_{i},U_{j})=\left|\mbox{Tr}\sqrt{\rho_{i}^{1/2}U_{j}\rho_{i}U_{j}^{\dagger}\rho_{i}^{1/2}}\right|^{2}.$
(16)

Of course, instead of the state-fidelity w.r.t. the input state, one could
have also chosen any other target state.

Again using a database of 1000 random states and 1000 unitaries, we can rate
the fidelities of $n_{r}$ elements out of a total of million elements. The
results are shown in Fig. 8. Interestingly, the rating is extremely successful
for pure states (Fig. 8(a)), for which up to 50% of unknown elements can be
predicted with a low RMS error of $\delta_{F}\leq 0.003$ ($n_{r}=5\times
10^{5}$). For mixed states the RMS errors are relatively higher, but still
$\delta_{F}\leq 0.02$ for predicting up to 50% of unknown elements. However,
in both cases, the rating fails for 90% of unknown elements.

Figure 8: Predicted fidelity $F$ versus actual fidelity $F_{0}$ and the number
$n_{r}$ of unknown elements for pure (a) and mixed (b) states. Blue lines
represent the ideal case $F=F_{0}$. RMS errors $\delta_{F}$ are show in each
case.

We have also studied the fidelity rating in multi-qubit registers. Fig. 9
displays the results of the fidelity rating in registers with up to 10 qubits.
In each case, we used a database of 100 states and 100 unitaries, while the
number of unknown entries was fixed to 100. It is clear that the fidelity
rating is largely successful in all these cases.

Figure 9: Predicted fidelity $F$ versus actual fidelity $F_{0}$ for pure (top
row) and mixed states (bottom row) of $n_{q}$-qubit registers as indicated. In
each case, the size of the database is $100\times 100$ and the number $n_{r}$
of unknown elements is fixed to $100$. Figure 10: The RS predictions of
nonunitary evolutions. Predicted values (red points) for changes in discord
(a) and fidelity (b) versus actual values and the number $n_{r}$ of unknown
elements. In each case, the size of the database is $1000\times 1000$. Blue
lines represent the ideal curves. RMS deviations between the predicted and
actual values are shown in each case.

##  V Rating nonunitary evolutions

In order to study rating of nonunitary evolutions we chose independent single-
qubit decoherence channels acting on a two-qubit system Nielsen and Chuang
(2002). Each channel has an operator-sum representation in terms of Kraus
operators $\\{M_{k}\\}$, which preserve the trace of the density operator by
obeying $\sum_{k}M_{k}^{\dagger}M_{k}=\mathbbm{1}$. A two-qubit input state
$\rho_{i}^{AB}$ is transformed according to

$\displaystyle{\cal E}_{j}(\rho_{i}^{AB})$ $\displaystyle=$
$\displaystyle\frac{1}{2}\sum_{k}(M_{k}^{j_{A}}\otimes\mathbbm{1})\rho_{i}^{AB}(M_{k}^{j_{A}^{\dagger}}\otimes\mathbbm{1})$
(17) $\displaystyle+\frac{1}{2}\sum_{l}(\mathbbm{1}\otimes
M_{l}^{j_{B}})\rho_{i}^{AB}(\mathbbm{1}\otimes M_{l}^{j_{B}^{\dagger}}),$

wherein channels $j_{A}$ and $j_{B}$ have acted on qubits $A$ and $B$
respectively. We consider the scenario wherein
$j_{A},j_{B}\in\\{\mathbbm{X},\mathbbm{Y},\mathbbm{Z},\mathbbm{D},\mathbbm{A}\\}$,
i.e., the two qubits undergo a random combinations of the following single-
qubit channels:

  * (i)

Bit-flip:
$M_{1}^{\mathbbm{X}}=\sqrt{x}~{}\mathbbm{1},~{}M_{2}^{\mathbbm{X}}=\sqrt{1-x}~{}\sigma_{x},$

  * (ii)

Phase-flip:
$M_{1}^{\mathbbm{Z}}=\sqrt{z}~{}\mathbbm{1},~{}M_{2}^{\mathbbm{Z}}=\sqrt{1-z}~{}\sigma_{z},$

  * (iii)

Bit & phase-flip:
$M_{1}^{\mathbbm{Y}}=\sqrt{y}~{}\mathbbm{1},~{}M_{2}^{\mathbbm{Y}}=\sqrt{1-y}~{}\sigma_{y},$

  * (iv)

Depolarization:
$M_{1}^{\mathbbm{D}}=\sqrt{1-\frac{3d}{4}}~{}\mathbbm{1},~{}M_{2}^{\mathbbm{D}}=\frac{\sqrt{d}}{2}~{}\sigma_{x},\\\
M_{3}^{\mathbbm{D}}=\frac{\sqrt{d}}{2}~{}\sigma_{y},~{}M_{4}^{\mathbbm{D}}=\frac{\sqrt{d}}{2}~{}\sigma_{z},$
and

  * (v)

Amplitude damping:
$M^{\mathbbm{A}}_{1}=\left[\begin{array}[]{cc}1&0\\\
0&\sqrt{1-a}\end{array}\right],M^{\mathbbm{A}}_{2}=\left[\begin{array}[]{cc}0&\sqrt{a}\\\
0&0\end{array}\right]$.

Here $x,y,z,a,d$ are randomly chosen channel probabilities and $\sigma_{x}$,
$\sigma_{y}$, and $\sigma_{z}$ are Pauli operators.

Again we form a database of 1000 mixed two-qubit states and 1000 random
nonunitary evolutions, each of which corresponds to a pair of channels from
the above, with a particular pair of channel probabilities, acting
independently on the individual qubits. Fig. 10 displays the predicted values
of change in discord

$\Delta D_{{i,j}}=D({\cal E}_{j}(\rho_{i}))-D(\rho_{i})$ (18)

(Fig. 10(a)) and of Fidelity

$F_{i,j}(\rho_{i},{\cal E}_{j})=\left|\sqrt{\sqrt{\rho_{i}}{\cal
E}_{j}(\rho_{i})\sqrt{\rho_{i}}}\right|^{2}$ (19)

(Fig. 10(b)), plotted versus the actual values for various numbers $n_{r}$ of
unknown elements. Since the decoherence channels only reduce quantum
correlations, $\Delta D$ is essentially negative. It is evident that the
prediction is highly successful, with the RMS deviation from the actual values
being 0.01 or lesser even at 90% of unknown elements. Thus, in summary, we
find that an RS can efficiently characterize both unitary and nonunitary
quantum evolutions in terms of either correlation changes or state-fidelity.

##  VI Conclusions

We have adopted a class of recommender systems, called the matrix
factorization algorithm, to characterize quantum evolutions in terms of
quantum correlations as well as state-fidelity. First, using two-qubit
databases, we have carried out a detailed analysis of the recommender system
predicting three types of correlations, namely entropy, negativity, and
discord. We found that the recommender system was able to efficiently predict
well over 50% of unknown elements in the database. By systematically varying
the database sizes and keeping the number of unknown ratings fixed, we confirm
that the prediction efficiency improves with the size of the database. As a
particular example, we show that predictions of negativity and discord of
Werner state prepared with different purities match very well with the
expected trend that entanglement exists only for purity above 1/3 and discord
at all purity $>0$.

By introducing noise into the database in a systematic way, we observed the
prediction efficiency deteriorating with noise. We proposed that this fact,
along with a reference database, could be used to distinguish a genuine or
clean database form a noisy or a fake database.

Predicting quantum discord of a general state, for which a closed-form
expression is not known, is most interesting. The analysis of computational
time hints that the machine prediction, within a fixed tolerance limit, can be
far more efficient. In million-element databases of two- and three-qubit
registers, we observed a speedup of one and two orders of magnitude
respectively, with less than 0.1 root-mean-square deviations.

We have also demonstrated the capability of the recommender system in
characterizing nonunitary evolutions by subjecting it to predict correlation
and fidelity ratings of states undergoing certain decoherence channels. The
feasibility of using the recommender system in larger systems is demonstrated
by successfully characterizing both unitary and nonunitary evolutions in
registers with up to ten qubits.

Although here we only discussed databases of random states and random
evolutions, it is straightforward to adapt it to other interesting scenarios
such as in predicting quantum correlations in a time-evolving many-body system
initialized with various possible states. It would also be interesting to
predict the behavior of quantum systems under evolutions dictated by more
sophisticated master equations. Further improvisations in the precision and
speed of machine predictions will go a long way as well. We anticipate that
such machine learning techniques not only pave the way for realizing quantum
technologies but also are helpful to gain deeper insights into the quantum
theory itself.

## Acknowledgments

Authors thank Krithika and Soham for useful discussions. P. B. gratefully
acknowledges support from the Prime Minister’s Research Fellowship and the
discussion meeting Statistical Physics of Machine Learning (ICTS/SPMML2020/01)
at International Centre for Theoretical Sciences Bangalore for all useful
discussions. A. S. gratefully acknowledges support from the Summer
Program-2019 at IISER Pune.

## References

  * Linden  _et al._ (2003) G. Linden, B. Smith, and J. York, IEEE Internet computing 7, 76 (2003).
  * Funk (2006) S. Funk, “Netflix update: Try this at home,” (2006).
  * Schafer  _et al._ (2001) J. B. Schafer, J. A. Konstan, and J. Riedl, Data mining and knowledge discovery 5, 115 (2001).
  * Lieberman (1997) H. Lieberman, in  _Proceedings of the ACM SIGCHI Conference on Human factors in computing systems_ (1997) pp. 67–74.
  * Maes (1995) P. Maes, in  _Readings in human–computer interaction_ (Elsevier, 1995) pp. 811–821.
  * Burke (2000) R. Burke, Encyclopedia of library and information systems 69, 175 (2000).
  * Goldberg  _et al._ (1992) D. Goldberg, D. Nichols, B. M. Oki, and D. Terry, Communications of the ACM 35, 61 (1992).
  * Schafer  _et al._ (2007) J. B. Schafer, D. Frankowski, J. Herlocker, and S. Sen, in  _The adaptive web_ (Springer, 2007) pp. 291–324.
  * Koren  _et al._ (2009) Y. Koren, R. Bell, and C. Volinsky, Computer 42, 30 (2009).
  * Torlai  _et al._ (2018) G. Torlai, G. Mazzola, J. Carrasquilla, M. Troyer, R. Melko, and G. Carleo, Nature Physics 14, 447 (2018).
  * Kieferová and Wiebe (2017) M. Kieferová and N. Wiebe, Phys. Rev. A 96, 062327 (2017).
  * Torlai and Melko (2017) G. Torlai and R. G. Melko, Phys. Rev. Lett. 119, 030501 (2017).
  * August and Ni (2017) M. August and X. Ni, Physical Review A 95, 012335 (2017).
  * Carrasquilla and Melko (2017) J. Carrasquilla and R. G. Melko, Nature Physics 13, 431 (2017).
  * Canabarro  _et al._ (2019) A. Canabarro, F. F. Fanchini, A. L. Malvezzi, R. Pereira, and R. Chaves, Phys. Rev. B 100, 045129 (2019).
  * Carleo and Troyer (2017) G. Carleo and M. Troyer, Science 355, 602 (2017).
  * Gao and Duan (2017) X. Gao and L.-M. Duan, Nature communications 8, 1 (2017).
  * Lu  _et al._ (2018) S. Lu, S. Huang, K. Li, J. Li, J. Chen, D. Lu, Z. Ji, Y. Shen, D. Zhou, and B. Zeng, Physical Review A 98, 012315 (2018).
  * Ma and Yung (2018) Y.-C. Ma and M.-H. Yung, npj Quantum Information 4, 1 (2018).
  * Gao  _et al._ (2018) J. Gao, L.-F. Qiao, Z.-Q. Jiao, Y.-C. Ma, C.-Q. Hu, R.-J. Ren, A.-L. Yang, H. Tang, M.-H. Yung, and X.-M. Jin, Physical review letters 120, 240501 (2018).
  * Yang  _et al._ (2019) M. Yang, C.-l. Ren, Y.-c. Ma, Y. Xiao, X.-J. Ye, L.-L. Song, J.-S. Xu, M.-H. Yung, C.-F. Li, and G.-C. Guo, Physical review letters 123, 190401 (2019).
  * Cimini  _et al._ (2020) V. Cimini, M. Barbieri, N. Treps, M. Walschaers, and V. Parigi, arXiv preprint arXiv:2003.03343 (2020).
  * Li  _et al._ (2019) X.-Y. Li, Q.-S. Zhu, M.-Z. Zhu, Y.-M. Huang, H. Wu, and S.-Y. Wu, EPL (Europhysics Letters) 127, 20009 (2019).
  * Wikipedia contributors (2020) Wikipedia contributors, “Matrix factorization (recommender systems) — Wikipedia, the free encyclopedia,”  (2020).
  * (25) Andrew Ng, Stanford Machine Learning Course, “https://www.coursera.org/learn/machine-learning,” .
  * Al Osipov  _et al._ (2010) V. Al Osipov, H.-J. Sommers, and K. Życzkowski, Journal of Physics A: Mathematical and Theoretical 43, 055302 (2010).
  * Maziero (2015) J. Maziero, Brazilian Journal of Physics 45, 575 (2015).
  * Plenio and Virmani (2014) M. B. Plenio and S. S. Virmani, in  _Quantum Information and Coherence_ (Springer, 2014) pp. 173–209.
  * Ollivier and Zurek (2001) H. Ollivier and W. H. Zurek, Physical review letters 88, 017901 (2001).
  * Katiyar  _et al._ (2012) H. Katiyar, S. S. Roy, T. Mahesh, and A. Patel, Physical Review A 86, 012309 (2012).
  * Lu  _et al._ (2011) X.-M. Lu, J. Ma, Z. Xi, and X. Wang, Physical Review A 83, 012327 (2011).
  * Nielsen and Chuang (2002) M. A. Nielsen and I. Chuang, “Quantum computation and quantum information,” (2002).

## Appendix

Identifying local and nonlocal operators: Suppose the objective is not to
obtain quantitative values for the changes in quantum correlations introduced
by unitary operators, but rather to identify local and nonlocal unitary
operators. For this purpose, we generate a database with 1000 randomly
generated unitaries, 500 of which are by construction local operators of the
form $U_{A}\times U_{B}$ that do not change quantum correlations. The
remaining $4\times 4$ dimensional random unitaries $U_{AB}$ are mostly
nonlocal and introduce changes in quantum correlations. Fig. 11 shows the
three-dimensional simultaneous plots of predicted entropy, negativity, and
discord values for various numbers $n_{r}$ of unknown ratings. It is clear
that most of the points fall into one of the two bins - negligible correlation
changes resulting from local unitaries, or significant correlation changes
brought about by nonlocal unitaries.

Figure 11: Binning unitaries into local or nonlocal by predicting the changes
in quantum correlations. Subplots correspond to different numbers $n_{r}$ of
predicted ratings as indicated in each case.


# Invariants of magnetic lines for Yang-Mills solutions

Akhmet’ev P.M.111National Research University “Higher School of Economics”,
Moscow, Russia; IZMIRAN, Troitsk, Moscow, Russia and Dvornikov M.S.222IZMIRAN,
Troitsk, Moscow, Russia

###### Abstract

We construct a new Yang-Mills $3D$-solution on the space of negative scalar
curvarure. We discuss a problem of non-abelian gauge symmetry is broken with
the assumption that a scalar curvature of the domain is a negative small
parameter. In this case we use the following fact: a geometrical scale related
with Vassiliev’s discriminant of magnetic lines coincids with a physical
Kolmogorov scale. This gives an estimation of $\alpha$-effect by the
dispersion of the asymptotic ergodic Hopf invariant in the limite with a
negative scalar curvature parameter.

##  1 Introduction

Since the theoretical discovery of the Yang–Mills theory in 1954 in Ref. [1],
the importance of non-abelian gauge fields for modern physics is enormous. It
is mainly because of the fact that the experimentally confirmed standard model
of elementary particles is based on $U(1)\times SU(2)\times SU(3)$ gauge
group. The main problem in the study of non-abelian gauge fields is that the
corresponding wave equations are generically nonlinear. It leads to the
appearance of various non-perturbative solutions, like instantons, which play
an important role in the construction of the models for the quantum
chromodynamics (QCD) vacuum [2]. Numerous non-perturbative methods for the
study of Yang–Mills fields are reviewed in Ref. [3].

Despite the fact that the non-abelian gauge symmetry of the standard model is
broken presently owing to the Higgs mechanism, this symmetry was restored in
the early universe before the electroweak phase transition (EWPT) which
happened at $T=100\,\text{GeV}$ or $t\sim 10^{-10}\,\text{s}$ since the Big
Bang [4]. All particles, including gauge bosons, were massless before EWPT.
The behavior of gauge bosons before EWPT can influence the evolution of
Maxwell electromagnetic fields in the present universe. It is applied
especially to the hypercharge field which is the progenitor of photons [5]. If
a nonzero hypermagnetic field existed in the early universe, it can result in
the leptogenesis and, hence, the baryogenesis owing to the analogue of the
Adler anomaly [6].

Non-abelian gauge fields in the early universe evolve in the non-Minkowski
spacetime. The spacetime of the expanding universe has the
Friedmann–Robertson–Walker (FRW) metric. Despite the fact that FRW metric is
conformally flat, the influence of nontrivial geometry on the evolution of a
Yang–Mills field has to be studied. Various quantum effects involving
Yang–Mills fields in curved spacetime are described in Ref. [7]. The
quantization of Yang–Mills fields in curved spacetime was reviewed in Ref.
[8].

The powerful tool for the study of the electromagnetic fields evolution in
highly conducting medium involves the magnetic hydrodynamics (MHD)
approximation. In this approach, one can exclude the electric field and deal
with a magnetic field only. As a result one obtains the induction equation
[9]. If P-symmetry is violated in a system, the induction equation acquires
the $\alpha$-term responsible for the magnetic field instability. The
mechanism for the $\alpha$-term generation in solar physics is described in
Ref. [10].

However, if one deals with the electromagnetic field in a curved spacetime,
the induction equation is difficult to formulate since one cannot separate the
electromagnetic field tensor into electric and magnetic fields in an invariant
manner. One can make it in the FRW metric using comoving variables (see, e.g.,
Ref. [11]). However, in this situation, one deals with the spectra of the
magnetic energy and the magnetic helicity. The same problem arises for non-
abelian fields. In this case, the magnetic field is a matrix valued function.
In the present work, we consider one of the ways to introduce the
$\alpha$-parameter for Yang–Mills fields in the invariant manner in a curved
spacetime.

It is clear that properties of $SU(2)$ cannot be totally described by means of
helicity integral and higher invariant of magnetic lines are into account. A
starting point for the present paper is a work [12], where the authors shows
that Chern-Simons integral is related with higher Massey product of magnetic
lines of Yang-Milles solution. However, the Massey product are not universal
invariants in MHD this invariant do not satisfy ergodic asymptotic properties,
as the helicity integral. In the textbook [13], based on an Arnold’s paper
[14], the asymptotic ergodic Hopf invariant is introduced, this invariant is
an expression of the helicity density of magnetic lines in an invariant form.
Several authors [16] introduced a notion of field line helicity, which
coincides with the convolution of asymptotic Hopf invariant in the case when
magnetic line are concentrated in a closed domain and are tangent to the
boundary of the domain. But in a general case, when magnetic lines are across
to a component of the boundary of the domain, an extra terms in the
construction of the field line helicity is required. At present a various
modifications of helicity integrals are studied: relative magnetic helicity
[15], winding magnetic helicity [17]. All this constructions are applied to
Maxwellian magnetic field, but allow no satisfactory description of properties
of $SU(2)$-fields (see Section 2 for helicity term in CS-integral).

In the paper by the first author [18] a construction of higher invariant of
magnetic lines for Maxwell field are introduced in an asymptotic and ergodic
form and in [19] are related with the $\beta$-term in the mean field theory.
From mathematical point of view there are difficulties in the definition. Even
in the case of helicity, the conception of the asymptotic ergodic Hopf
invariant in the case magnetic fields is not generic gives difficulties, see
[13] Ch.III, Sec. 4.3. In the case magnetic lines are closed (this case is
also non-generic, one speaks about magnetic field concentrated in magnetic
tubes) we get well-defined notion of density of higher magnetic invariant
($M$-invariant), but not of higher Massey invariant originally was introduced
in MHD in [20].

As we have shown linking number of field lines, or, the magnetic helicity,
allows no description of magnetic fields in a curved domain, higher invariants
of magnetic lines are required. Based on Kolmogorov’s approach on
hydrodynamics, which is known in MHD as the Kolmogorov’s magnetic spectrum, we
study a contribution of the scalar curvature parameter of the domain to
$\alpha$-term in the mean field equation [21]. The construction is based on
the idea that a scale of the space, which a priori give contribution to higher
invariant of magnetic lines, has to lead invariant gauge. This is observed in
the case the density of magnetic lines is distributed as a magnetic flow of
Kolmogorov spectrum [19]. If we assume that $M_{5}$-invariant of magnetic
lines distributes in the domain uniformly with the scale, we get the
Kolmogorov’s $\beta$-term scale distribution. But, the distribution of
magnetic lines with uniformly $M_{5}$-density depends on the scalar curvature
parameter and its asymptotic for small curvature parameter is easily
calculated. With an additional assumption that $\alpha$-term is distributed
with the same scale, we get the required estimation of the asymptotic of the
$\alpha$-parameter with respect to the scalar curvature of the domain.

This work is organized in the following way. In Sec. 1, we construct a Yang-
Milles solution on the standard 3D-sphere $S^{3}$ with positive scalar
curvature and give a short history of the question. In Sec. 2, we generalize
the solution for the domain with negative scalar curvature. We use
calculations from the textbook [13] and give a Yang-Milles solution for this
domain. The Yang-Milles solution is bracken, and we get a Maxwell field with
the Kolmogorov distribution, which is closed to a component of Yang-Milles
field. We have a non-compact domain in this case, and we assume that a field
is quasi-periodic. This assumption means, that the Chern-Simons wave, which is
described the large scale component of magnetic field, is closed to geodesic
flow in the domain. In Sec. 3 we recall properties of $M_{5}^{sim}$ invariant
[19] and calculate this invariant for a component of Yang-Milles field first,
and for closed Maxwell magnetic field, in geometrical terms of Vassiliev
discriminant theory [22]. This geometrical approach gives an integral
contribution of asymptotic of linking number of closed magnetic lines.
Physically this interpreted as linking numbers of magnetic lines in the
magnetic flow in mean-field equations, and this gives an asymptotic of the
$\alpha$-term by the scalar curvature parameter. In Sec. 4 we give
conclusions.

##  2 Yang-Mills field and Chern-Simons invariant

In this section a connection $\vec{\bf A}^{R}$ in $SU(2)$-bundle $\eta^{R}$
over the standard 3D sphere $S^{3}$ is studied. The bundle $\eta^{R}$ is
isomorphic to the trivial bundle $S^{3}\times{\bf H}$. Let us consider $S^{3}$
as the unite sphere in ${\bf H}$ with the base vectors $\\{+1,+{\bf i},+{\bf
j},+{\bf k}\\}$. We identify the bundle $\eta^{R}$ by the isomorphism:

$\phi^{R}:{\bf H}\to\eta^{R},\quad\phi^{R}(\xi)=\xi\circ{\bf x},\quad{\bf
x}\in S^{3},\quad\xi\in{\bf H}.$

In this formula $\circ$ is the quaternion product. This formula determines
$\vec{\bf A}^{R}$ as $\phi^{R}$-image of imaginary quaternions base vectors
$\\{{\bf i},{\bf j},{\bf k}\\}$. By this definition $\vec{\bf A}^{R}$ is the
Levi-Civita connection over the bundle $\eta^{R}$ (as the real bundle this is
4D bundle). In the construction we use the metric, and we write-down
connections and its curvatures as vectors in the standard coordinate system.
There is an isomorphism $T(S^{3})=Im(\eta^{R})$, where $T(S^{3})$ is the
tangent bundle, $Im(\eta^{R})\subset\eta^{R}$ is a 3D subbundle generated by
imaginary quaternions. Also this formula determines the right-polarized frame
on $T(S^{3})$ by the $\phi^{R}$-image of the imaginary base quaternions.

The right-action by the base quaternions in the fiber of $\eta^{R}$ is an
isometry of the fiber, this isometry is described by the corresponding matrix,
which is called Pauli matrix.

Yang-Milles field is a $2$-form ${\vec{\bf F}}$, which is the curvature of the
connection $\vec{\bf A}$ in $SU(2)$-bundle $\eta$. This form satisfies the
equations:

$\displaystyle D{\vec{\bf F}}=0;\quad D^{*}({\vec{\bf F}})=0.$ (1)

The curvature ${\vec{\bf F}}^{R}$ of the connection $\vec{\bf A}^{R}$
satisfies the equation (1(1)), but not the equation (1(2)). To get a Yang-
Milles field, we replace $\vec{\bf A}^{R}\mapsto\vec{\bf A}^{L}$ as following.

Recall, the fiber $\eta_{{\bf x}}$ over a point ${\bf x}\in S^{3}$ is defined
by the right action on the quaternion ${\bf x}$ on the standard fiber ${\bf
H}$. Let us define the left action by the conjugated quaternion $\bar{{\bf
x}}$, namely:

$\phi^{L}:{\bf H}\to\eta^{L},\quad\phi^{L}(\xi)=\bar{{\bf
x}}\circ\xi,\quad{\bf x}\in S^{3},\quad\xi\in{\bf H};$ $\bar{{\bf
x}}=\overline{a+b{\bf i}+c{\bf j}+d{\bf k}}=a-b{\bf i}-c{\bf j}-d{\bf k}.$

Let us define the connection $\vec{\bf A}^{L}$ in $\eta^{L}$. Let us define
the Riemannian metric on $\eta^{L}$ analogously with the case $\eta^{R}$
above. The fiber $\eta^{L}_{{\bf x}}$, ${\bf x}=+1\in S^{3},$ is isometric to
${\bf H}$. The triple

$\vec{\bf A}^{L}=(\vec{\bf A}_{i}^{L},\vec{\bf A}_{j}^{L},\vec{\bf
A}_{k}^{L})$

is defined as the result of the left translation of the base quaternions ${\bf
i}_{\eta},{\bf j}_{\eta},{\bf k}_{\eta}\\}$ from the fiber $\eta^{L}_{+1}$
into the space of the bundle. Because $Im(\eta^{L})$ is isomorphic to
$T(S^{3})$, we may assume that $\vec{\bf A}^{L}$ is the frame in $T(S^{3})$,
which is a $1$-form of the connection. In the prescribed coordinate system we
may consider $\vec{\bf A}^{L}$ as vectors, by definition $(\vec{\bf
A}_{i}^{L},\vec{\bf A}_{j}^{L},\vec{\bf A}_{k}^{L})=(\bar{{\bf i}},\bar{{\bf
j}},\bar{{\bf k}})$, where in the right hand side we get the conjugated
quaternion frame on $T(S^{3})$.

Let us define a $2$-form ${\vec{\bf F}}^{L}$, the Yang-Milles field, using the
Riemannian metric, this form is defined as the vector-field by the formula:

${\vec{\bf F}}^{L}=(-\vec{\bf A}_{i}^{L},-\vec{\bf A}_{j}^{L},-\vec{\bf
A}_{k}^{L}).$

In this case the equation (1;1) is written as:

${\vec{\bf F}}^{L}={\operatorname{rot}}\ \vec{\bf A}_{l}^{L}+[\vec{\bf
A}^{L}_{l+1},\vec{\bf A}^{L}_{l+2}],$

where $[,]$ is the commutator $[\vec{\bf A}_{l}^{L},\vec{\bf
A}_{l+1}^{L}]=\vec{\bf A}_{l+2}^{L}$, $l={\bf i},{\bf j},{\bf k}$,
${\operatorname{rot}}$ is the standard vorticity operator on $S^{3}$. Because

$D^{\ast}d(\vec{\bf A}^{L})={\operatorname{rot}}\
{\operatorname{rot}}(\vec{\bf A}_{l}^{L})+[{\operatorname{rot}}(\vec{\bf
A}_{l+1}^{L}),\vec{\bf A}_{l+2}^{L}]+[\vec{\bf
A}_{l+1}^{L},{\operatorname{rot}}(\vec{\bf A}_{l+2}^{L})]=0,$

the equation (1(2)) is written in the form:

${\operatorname{rot}}[\vec{\bf A}_{l+1}^{L},\vec{\bf A}_{l+2}^{L}]+[\vec{\bf
A}_{l+1}^{L},[\vec{\bf A}_{l}^{L},\vec{\bf A}_{l+1}^{L}]]+[\vec{\bf
A}_{l+2}^{L},[\vec{\bf A}_{l}^{L},\vec{\bf A}_{l+2}^{L}]]=0.$

The equations (1(1))(1(2)) are satisfied because of the equation:

${\operatorname{rot}}\vec{\bf A}_{l}^{L}=-2\vec{\bf A}_{l}^{L},\quad[\vec{\bf
A}_{l}^{L},\vec{\bf A}_{l+1}^{L}]=\vec{\bf A}_{l+2}^{L}.$

The following integral:

$\displaystyle CS({\vec{\bf F}}^{L})=\int_{S^{3}}tr(\vec{\bf A}^{L}\wedge
d\vec{\bf A}^{L}+\vec{\bf A}^{L}\wedge\vec{\bf A}^{L}\wedge\vec{\bf A}^{L})$
(2)

is called the Chern-Simons functional.

The first term of the integral (2) is defined as the average linking number of
magnetic lines of the fields:

$\displaystyle\vec{\bf B}_{{\bf i}}^{L}={\operatorname{rot}}\vec{\bf A}_{{\bf
i}}^{L},\vec{\bf B}_{{\bf j}}^{L}={\operatorname{rot}}\ \vec{\bf A}_{{\bf
j}}^{L},\vec{\bf B}_{{\bf k}}^{L}={\operatorname{rot}}\vec{\bf A}_{{\bf
k}}^{L}.$ (3)

The second term is an average number of higher linking numbers of magnetic
lines of the triple: $\\{\vec{\bf B}_{{\bf i}}^{L},\vec{\bf B}_{{\bf
j}}^{L},\vec{\bf B}_{{\bf k}}^{L}\\}$.

### An history of the question

In a paper by Hornig and Mayer [12] for the right-hand fields the equation is
satisfied:

${\operatorname{rot}}\vec{\bf A}_{l}^{R}=+2\vec{\bf A}_{l}^{R}.$

In this case Yang-Milles equation (1(2)) is not satisfied.

By the quaternion conjugation $\vec{\bf A}^{R}\mapsto\vec{\bf A}^{L}$ the
proper number of the operator ${\operatorname{rot}}$ changes the sign. The
equation:

$[\vec{\bf A}_{l}^{L},\vec{\bf A}_{l+1}^{L}]=\vec{\bf A}_{l+2}^{L}$

is satisfied, because the base $\eta^{L}_{+1}$ remains standard.

The space $S^{3}$ from the holomorphic point of view is the right-polarized
space. One may considered this space alternatively as the left-polarized space
by a different representation. In this construction $S^{3}$ is the
$Spin$-covering for the geodesic flow on the standard $2D$-sphere $S^{2}$. Let
us shows that in this caser Yang-Mills equations are generalized.

##  3 Geodesic flows of positive and negative curvature

In this section we use definition and calculations from [13] Ch.5; 4.4. The
previous calculation for Yang-Mills field on $S^{3}$ are reformulated for the
$Spin$-covering over the space

$SO(3)=\\{\xi\in T_{\ast}(S^{2}):||\xi||=1\\}.$

All the calculations remains analogous, if we replace $S^{2}$ into the
Lobachevski plane $\Lambda$.

Let us consider the space

$T\Lambda=\\{\xi\in T_{\ast}(\Lambda):||\xi||=1\\},$

where $\Lambda$ is the Lobachevski plane with the unite negative scalar
curvature $-1$. Let us define the right triple of vector fields:

$\\{f^{R},e^{R},\tilde{e}^{R}\\}.$

The field $f^{R}$ is the isometric rotation trough a positive angle in the
prescribed fiber, $e^{R}$ is the geodesic flow, $\tilde{e}^{R}$ is the
conjugated geodesic flow on the Lobachevski plane. The sectional curvature on
the plane $(e^{R},\tilde{e}^{R})$ equals to $-1$, sectional curvatures in the
planes $(f^{R},e^{R})$, $(f^{R},\tilde{e}^{R})$ equal to $-2$. The last two
equations follow is proved as following. Replace $\Lambda$ by $S^{2}$, the
scale in the vertical fiber $S^{1}\subset S^{3}$ over $\tilde{T}(S^{2})$ is
normalized by the standard metric in the base of the double covering $\times
2:S^{1}\to S^{1}$. The calculation for $S^{2}\mapsto\Lambda$ is analogous, the
only difference that the scalar positive curvature in all sections changes the
sign.

Take the parameter $\lambda>0$, the scale of the metric of the vertical
coordinate $f^{R}$ such that the following equations are satisfied:

$(e^{R},e^{R})=(\tilde{e}^{R},\tilde{e}^{R})=1;\quad(f^{R},f^{R})=\lambda^{2}.$

In this case the following equations are satisfied:

${\operatorname{rot}}\ e^{R}=-\frac{e^{R}}{\lambda};\quad{\operatorname{rot}}\
\tilde{e^{R}}=-\frac{\tilde{e^{R}}}{\lambda};\quad{\operatorname{rot}}\
f^{R}=f^{R}.$

The volume form, generated by this vector fields equals to
$\mu(e^{R},\tilde{e}^{R},f^{R})=\lambda$. The parameter $\lambda$, by this
equation, changed the phase volume into $\lambda$ times.

Let us pass into the left base, which is given by the following triple of the
vector fields:

$\\{f^{L}=-\frac{f^{R}}{\lambda},\ e^{L}=e^{R},\
\tilde{e}^{L}=\tilde{e}^{R}\\}.$

For this base we get the relation:

${\operatorname{rot}}\ e^{L}=-\frac{e^{L}}{\lambda};\quad{\operatorname{rot}}\
\tilde{e^{L}}=-\frac{\tilde{e^{L}}}{\lambda};\quad{\operatorname{rot}}\
f^{L}=-\frac{f^{L}}{\lambda}.$

Let us denote by $\tilde{T}(\Lambda)$ the double covering over $T(\Lambda)$.
The Riemannian metric in this space is taken analogously with the case
$\tilde{T}(S^{2})$: in the horizontal plane this metric relates with the
metric on $T(\Lambda)$, in the vertical circle the length of the vector is
twice greater then the angle length on the base. For $\lambda=1$ we get a
manifold with a constant scalar curvature, because all sectional curvatures
equal $-1$. In a general case the sectional curvature in the horizontal plane
$(e^{L},\tilde{e}^{L})$ equals to $-1$, and sectional curvatures in the planes
$(f^{L},e^{L})$, $(f^{L},\tilde{e}^{L})$ equal to $\lambda^{-1}$.

In the Riemannian manifold $\tilde{T}(\Lambda)$ the base fields is re-denoted
by the following gauge:

$f^{L}\mapsto-2f^{L},\quad e^{L}\mapsto
e^{L},\quad\tilde{e}^{L}\mapsto\tilde{e}^{L},$

the resulting fields are not re-denoted.

For the fields on $\tilde{T}(\Lambda)$ the following relations are satisfied:

${\operatorname{rot}}\
e^{L}=-2\frac{e^{L}}{\lambda};\quad{\operatorname{rot}}\
\tilde{e}^{L}=-2\frac{\tilde{e}^{L}}{\lambda};\quad{\operatorname{rot}}\
f^{L}=-2\frac{f^{L}}{\lambda}.$

Let us consider the case $\lambda=1$. As in the previous section, the triple
$\\{f^{L},e^{L},\tilde{e}^{L}\\}$ determines the Yang-Milles solution.

Let us consider the case of an arbitrary $\lambda>0$. The scale of the space
is normalized by its volume. By this scale the fields

$\displaystyle\vec{\bf A}_{{\bf
i}}^{\lambda}=\lambda^{\frac{1}{3}}f^{L};\quad\vec{\bf A}_{{\bf
j}}^{\lambda}=\lambda^{\frac{1}{3}}e^{L},\quad\vec{\bf A}_{{\bf
k}}^{\lambda}=\lambda^{\frac{1}{3}}\tilde{e}^{L}$ (4)

are normalized base fields. The sectional curvature in the plane $({\bf
j},{\bf k})$ (the horizontal plane) equals to $\lambda^{-\frac{2}{3}}$. The
sectional curvature in the plane: $({\bf i},{\bf j})$, $({\bf i},{\bf k})$
contains $\lambda^{-\frac{5}{3}}$. The following equations are satisfied:

${\operatorname{rot}}\vec{\bf A}^{\lambda}_{{\bf
i}}=-\frac{2}{\lambda}\vec{\bf A}^{\lambda}_{{\bf
i}},\quad{\operatorname{rot}}\vec{\bf A}^{\lambda}_{{\bf
j}}=-\frac{2}{\lambda}\vec{\bf A}^{\lambda}_{{\bf
j}},\quad{\operatorname{rot}}\vec{\bf A}^{\lambda}_{{\bf
k}}=-\frac{2}{\lambda}\vec{\bf A}^{\lambda}_{{\bf k}}.$

In the case $\lambda\to+\infty$ the space become flat.

Because the scale of the space determines a proportional scale of the
brackets, the relation

$[\vec{\bf A}_{{\bf i}}^{\lambda},\vec{\bf A}_{{\bf j}}^{\lambda}]=\vec{\bf
A}_{k}^{\lambda},\quad[\vec{\bf A}_{{\bf k}}^{\lambda},\vec{\bf A}_{{\bf
i}}^{\lambda}]=\vec{\bf A}_{{\bf j}}^{\lambda},\quad[\vec{\bf A}_{{\bf
j}}^{\lambda},\vec{\bf A}_{{\bf k}}^{\lambda}]=\vec{\bf A}_{{\bf
i}}^{\lambda}$

remains as in the case $\lambda=1$.

Let us denote $\vec{\bf B}_{{\bf i}}^{\lambda}=-\frac{2}{\lambda}\vec{\bf
A}_{{\bf i}}^{\lambda}$, $\vec{\bf B}_{{\bf
j}}^{\lambda}=-\frac{2}{\lambda}\vec{\bf A}_{{\bf j}}^{\lambda}$, $\vec{\bf
B}_{{\bf k}}^{\lambda}=-\frac{2}{\lambda}\vec{\bf A}_{{\bf k}}^{\lambda}$. Let
us calculate densities of magnetic helicities. We get:

$(\vec{\bf A}_{{\bf i}}^{\lambda},\vec{\bf B}_{{\bf i}}^{\lambda})=(\vec{\bf
A}_{{\bf j}}^{\lambda},\vec{\bf B}_{{\bf j}}^{\lambda})=(\vec{\bf A}_{{\bf
k}}^{\lambda},\vec{\bf B}_{{\bf k}}^{\lambda})=-2.$

We see that a sub-integral kernel of the each field is constant and depends
not of the parameter $\lambda$. As the result, the density of magnetic
helicities, the first term in the integral, depends not on sectional
curvatures.

Let us consider the second term of the integral (2) as a function of the
parameter $\lambda$. We get that the sub-integral 3-form

$\displaystyle\vec{\bf A}_{{\bf i}}^{\lambda}\wedge\vec{\bf A}_{{\bf
j}}^{\lambda}\wedge\vec{\bf A}_{{\bf k}}^{\lambda}$ (5)

is proportional to $\lambda^{-1}$. Therefore the density of the second term in
2 is proportional to $\lambda^{-1}$.

Let us present an easy geometrical reason of this calculation. Accordingly
Hornig and Mayer [12], an analogous case of magnetic fields on
$\tilde{T}(S^{2})$), the density $\rho$ of the second term in the integral (2)
is given by the triples magnetic lines distribution, which are detected by its
starting points and its tangent vectors, and which are pairwise linked in the
phase space.

Let us consider a point $x\in\tilde{T}\Lambda$, consider its projection
$\bar{x}$ on $\Lambda$ and calculate the density of magnetic lines, starting
in a neighborhood of $\bar{x}$, which are pairwise intersected, as a function
of the parameter $\lambda$.

Let us consider $3$ concentric circles of a great radius $R_{1},R_{2},R_{3}$
in the central point $\bar{x}$. The ratio of the parallelism angle from
$x_{1}\in R_{1}$ to the line $l_{2}$ and the perimeter of the circle $R_{1}$
equals to $\lambda^{-\frac{1}{3}}$, when $R_{1}\to+\infty$. The integral
asymptotic average of the geodesic triangles when radii tends to $+\infty$ is
of the order $\lambda^{-1}$, because pairwise intersections are independently.
As the result, we get $\rho\sim\lambda^{-1}$, this gives the asymptotic of the
second term in (2), as a function of the scalar curvature of the horizontal
section.

##  4 $M$-invariant of magnetic lines

We assume that Yang-Mills equation for (4), is broken. Let us consider the
only potential, say, the potential $\vec{\bf A}_{j}$, which generates the
magnetic field $\vec{\bf B}_{{\bf j}}$, which is closed to the Yang-Milles
component. Let us estimate the square of the density distribution of triangles
of magnetic lines of the field $\vec{\bf B}_{{\bf j}}$ on the horizontal
plane. This value is of the order $\sim\lambda^{-\frac{2}{3}}$. An analogous
calculation is presented in [23], Lemma 5.1, with the assumption that there
exist an additional symmetry of magnetic lines. (Note, that in [23] a scale of
the angle and, therefore, the Riemannian metric in $\tilde{T}(\Lambda)$ is
taken an alternative way, this leads to non-suitable choice of the parameter
$a$ in the statement of Main Theorem in Section $3$).

To estimate the square of the density, let us take all triangles, with
intersections of the sides on angles $\geq\varepsilon$, where $\varepsilon>0$
is a small positive parameter. Then, let us take the limit $\varepsilon\to
0+$. This limit is of the order $\sim\lambda^{-\frac{1}{3}}$, the square is of
the order $\sim\lambda^{-\frac{2}{3}}$

Let us consider the rescaling into $l$-times: ${\bf x}\mapsto l{\bf x}$. In
this case the volume form is changed in $l^{3}$ times. The density of the
second term in the integral (2) remains unchanged, because each potential (4),
which is included in (5) is changed by the factor $l^{-1}$. In the mean field
theory [21], the magnetic field $\vec{\bf B}_{{\bf j}}$ can be consider as a
large scale magnetic field. From physical assumption that a small-scale
magnetic field is Kolmogorov distributed, with the exponent $-\frac{5}{3}$. By
this assumption, the turbulent diffusion generate a magnetic flow with the
exponent $-\frac{7}{6}$ (the $\beta$-term in the mean field equation).
Moreover, because the $\alpha$-term, the small-scale component twists a little
magnetic lines of the large-scale magnetic field $\vec{\bf B}_{{\bf j}}$. A
origin of the $\alpha$-effect and its amplitude is outside of our
investigation. We only assume that $\alpha$ and $\beta$-terms in the mean
field equation corresponds each other, and the $\alpha$-term is small. With
this assumption we estimate the $\alpha$-term with respect to the curvature
parameter $\lambda$.

Accordingly to [18], a higher asymptotic ergodic invariant $M^{sim}_{5}$ of
magnetic lines has a density exponent, which coincides with the exponent of
the magnetic flow of the Kolmogorov spectrum in $k$-space. In the density of
$M^{sim}_{5}$-invariant for the magnetic field $\vec{\bf B}_{{\bf j}}$ is
presented an only non-trivial term. This is the main term described in Theorem
17 [19] (an auxiliary formula), with respect to the contribution of linking
numbers of magnetic lines, which is presented in Theorem 17 [18]. This term is
quadratic with respect to the magnetic flow, it calculates the square of the
density of geodesic triangles of the projections of Chern-Simons magnetic
line. Each triangle is taken proportional to the product of squares of
magnetic flows. Accordingly to the combinatorial formula of
$M_{5}^{sim}$-invariant[18], in a random quintuple of magnetic lines all
triples of lines are considered. In the case a triple represents a triangle,
its contribution is taken proportionally to the product of $10$ pairwise
linking coefficients in the considered quintuple of magnetic lines.

We assume that the parameter $\lambda^{-1}$ is small. Let us divide quintuples
of magnetic lines into two groups: quintuples of pairwise far lines and all
other quintuples, which contain at least two closed magnetic lines. A pair of
far magnetic lines, which represents by a pair of curves in the horizontal
plane, closed to geodesics, are linked with the coefficient $-1$ (the only
intersection points of this pair), or, unliked, when the projection are
disjoint. Below it was shown that pairwise coefficients are distributed with a
density which is not depended on $\lambda$, the average values of the square
of triangles is of the order $\lambda^{-\frac{2}{3}}$. Because we assume that
the integral over closed magnetic lines is small, we may assume that a
contribution of triangles with at least two closed sides into
$M^{sim}_{5}$-integral is proportional to $\alpha^{2}$, because this term
determines a contribution of the second moment of the density of Hopf
asymptotic ergodic invariant. The dispersion of the Hopf asymptotic ergodic
invariant is investigated in [24]. From this we get that an assumption that
$\alpha$ and $\beta$-terms in the mean field equation are balanced, we get the
asymptotic $\alpha\sim\lambda^{-\frac{1}{3}}$.

##  5 Discussion

The well-known Gauss invariant (and its modern generalization suggested by
Arnold) plays an important role in magnetic field generation in celestial
bodies. The Arnold form of the Gauss integral admits a gauge invariant large-
scale density. One may assume that the second-order time-derivative of this
density is homogenous in the domain by a small-scale velocity, this gives an
additional property for large-scale magnetic flow. In the present work, using
this idea we have studied the following question: if the $\alpha$-parameter in
dynamo mean field equation depends on the scalar curvature of the space.

We assume that the large-scale magnetic field, Chern-Simons magnetic wave, is
closed to a component of a Yang-Mills solution and a small scale magnetic
field satisfies the Kolmogorov distribution in $k$-space. Because of this
assumption, the $\beta$-term is predicted by singularities of field lines of
large scale magnetic field, the $\alpha$-term and $\beta$-term have to get a
common distribution. This gives an asymptotic of $\alpha$-term with respect to
the scalar curvature parameter.

After a short introduction in Sec. 1, we have constructed a new solution of
Yang-Mills equation on the standard 3D sphere, based on an ansatz proposed by
Hornig and Mayer [12], in Sec. 2. In Sec. 3, we have translated Yang-Mills
solutions for the space of a negative scalar curvature, Chern-Simons wave is
constructed as geodesic flow on the Lobachevsky plane. We have used definition
and calculations from Arnold-Khesin textbook [13].

In Sec. 4 we have proved the main result, assuming that the Chern-Simons wave,
which is closed to a component of Yang-Milles solution, represents a large
scale magnetic Maxwell field and admits the Kolmogorov’s distribution. This
assumption with the main result of the paper [18] by the first author proves
that the magnetic flow is visualized using magnetic lines singularities (using
geodesic triangles) in the framework of Vassiliev theory. Because calculations
of singularities of fields lines depend on the curvature in the domain, we get
a relationship of physical parameters in mean fields equation with the
curvature. This gives an asymptotic for a small curvature parameter.

Our results in Sec. 4 require that the curvature of a spacetime to be nonzero,
namely, negative. According to the measurements present in Ref. [25], our
universe is flat with the high level of accuracy. However, one can consider
perturbations of a metric which can result in a nonzero curvature.

The contribution of the gravity to the dynamo $\alpha$-parameter can result
from a nonzero torsion terms (see, e.g., Ref. [26]). In the present work, we
have studied another mechanism for the contribution of gravity to the
$\alpha$-parameter in torsionless spaces. It is valid for non-abelian Yang-
Mills fields. This contribution can be important, e.g., for the study of
gravity perturbations on the evolution of a hypercharge field in the early
universe.

The perturbations of gravitational field are known to cause the deviation of
the beam of light propagation from the geodesics. For the first time, this
phenomenon was predicted in Ref. [27]. This subject was reviewed in Ref. [28].
The recent mathematical development of this problem is given in Ref. [29].
Thus, observing the propagation of photons, emitted by distant sources, one
can study the spectral distribution of density perturbations in the universe
[30].

According to the Einstein equations, matter density perturbations result in
fluctuations of gravitational fields. Relying on our results, one can study
the contribution of these density fluctuations to the $\alpha$-parameter. In
this vein, we can explore the dynamics of stochastic primordial
(hyper-)magnetic fields. Note that the research of stochastic cosmic magnetic
fields was also carried out in Ref. [31]. The proposed study is important for
observational constraints on primordial magnetic fields [32].

## Acknowledgments

We are thankful to V.B.Semikoz for useful discussions.

## References

  * [1] C.-N. Yang and R. L. Mills, Conservation of Isotopic Spin and Isotopic Gauge Invariance, Phys. Rev. 96, 191–195 (1954).
  * [2] E. Shuryak and T. Schäfer, The QCD vacuum as an instanton liquid, Ann. Rev. Nucl. Part. Sci. 47, 359–394 (1997).
  * [3] M. Q. Huber, Nonperturbative properties of Yang–Mills theories, Phys. Rept. 879, 1–92 (2020) [arxiv:1808.05227].
  * [4] D. S. Gorbunov and V. A. Rubakov, Introduction to the Theory of the Early Universe: Hot Big Bang Theory (World Scientific, Singapore, 2011), pp. 217–218.
  * [5] M. Dvornikov and V. B. Semikoz, Lepton asymmetry growth in the symmetric phase of an electroweak plasma with hypermagnetic fields versus its washing out by sphalerons, Phys. Rev. D 87, 025023 (2013) [arxiv:1212.1416].
  * [6] M. Giovannini and M. E. Shaposhnikov, Primordial hypermagnetic fields and triangle anomaly, Phys. Rev. D 57, 2186–2206 (1998) [hep-ph/9710234].
  * [7] A. A. Grib, S. G. Mamaev, and V. M. Mostepanenko, Vacuum Quantum Effects in Strong Fields (Atomizdat, Moscow, 1980).
  * [8] S. Hollands, Renormalized Quantum Yang–Mills Fields in Curved Spacetime Rev. Mod. Phys. 20, 1033–1172 (2008) [arxiv:0705.3340].
  * [9] L. D. Landau and E. M. Lifshitz, Electrodynamics of Continuous Media (Pergamon, Oxford, 1984), p. 228.
  * [10] M. Stix, The Sun: An Introduction (Springer, Berlin, 2004), 2nd ed., pp. 368–372.
  * [11] M. Dvornikov and V. B. Semikoz, Instability of magnetic fields in electroweak plasma driven by neutrino asymmetries, J. Cosmol. Astropart. Phys. 05, 002 (2014) [arxiv:1311.5267].
  * [12] G. Hornig and Ch. Mayer, Towards a third-order topological invariant for magnetic fields, J. Phys. A: Math. Gen. 35 3945 arXiv:physics/02 (2002).
  * [13] Arnol’d V. I. and Khesin B. A., Topological Methods in Hydrodynamics Applied Mathematical Sciences vol 125 (1998) (New York: Springer-Verlag)
  * [14] Arnol’d V. I. (1974) Proc. Summer School in Differential Equations (Erevan) Armenian SSR Acad. Sci. [English translation: 1986 Sel. Math. Sov. 5 327–45].
  * [15] MacTaggart, David; Valli, Alberto, Relative magnetic helicity in multiply connected domains, J. Phys. A: Math. Theor. 58 (2023) 435701 (15pp).
  * [16] Russell A.J., Yeates A.R., Hornig G. and Wilmont-Smith A.L. (2015) Evolution of field line helicity during magnetic reconnection Phys.Plasmas 22 032106.
  * [17] Xiao, Daining; Prior, Christopher B.; Yeates, Anthony R., Spherical winding and helicity, J. Phys. A: Math. Theor. 56 (2023) 205201 (31pp).
  * [18] Peter M. Akhmet’ev, On a higher integral invariant for closed magnetic lines, revisited, Journal of Geometry and Physics 170 (2021) 104379.
  * [19] P.M.Akhmet’ev, Topological meaning of the slope of the Kolmogorov spectrum of magnetic turbulence: $M_{5}$-invariant of magnetic lines and its combinatorial formula, Journal of Geometry and Physics 178(2022)104583.
  * [20] M.V. Monastyrsky, V.S. Retakh, Topology of linked defects in condensed matter, Commun. Math. Phys. 103 (1986) 445–459.
  * [21] H.K. Moffatt, Magnetic Field Generation in Electrically Conducting Fluids, Cambridge University Press, (1978).
  * [22] Arnold, V.I. (1994). The Vassiliev Theory of Discriminants and Knots. In: Joseph, A., Mignot, F., Murat, F., Prum, B., Rentschler, R. (eds) First European Congress of Mathematics . Progress in Mathematics, vol 3. Birkhäuser Basel. $https://doi.org/10.1007/978-3-0348-9110-3_{1}$
  * [23] P. M. Akhmet’ev, Knot Invariants in Geodesic Flows, Differential equations and dynamical systems, Collected papers, Trudy Mat. Inst. Steklova, 308, Steklov Math. Inst. RAS, Moscow, 2020, 50–64; Proc. Steklov Inst. Math., 308 (2020), 42–55.
  * [24] P. Akhmet’ev and I. V. Vyugin, Dispersion of the Arnold’s Asymptotic Ergodic Hopf Invariant and a Formula for its Calculation, Arnold Math. J. 6, 199–211 (2020) [arXiv:1906.12131].
  * [25] N. Aghanim et al. (Planck Collaboration), Planck 2018 results. VI. Cosmological parameters, Astron. Astrophys. 641, A6 (2020) [arXiv:1807.06209].
  * [26] K. Bamba, C.-Q. Geng, and L.-W. Luo, Generation of large-scale magnetic fields from inflation in teleparallelism, J. Cosmol. Astropart. Phys. 10, 058 (2012) [arxiv:1208.0665].
  * [27] Ya. B. Zel’dovich, Observations in a Universe Homogeneous in the Mean, Sov. Astron. 8, 13 (1964).
  * [28] S. R. Dolan, Geometrical optics for scalar, electromagnetic and gravitational waves on curved spacetime, Int. J. Mod. Phys. D 27, 1843010 (2017) [arxiv:1806.08617].
  * [29] E. A. Illarionov and D. D. Sokoloff, Relative efficiency of three mechanisms of vector field growth in a random media, Phys. Rev. E 107, 044110 (2023) [arXiv:2303.14388].
  * [30] P. Fleury, J. Larena, and J.-P. Uzan, Weak lensing distortions beyond shear, Phys. Rev. D 99, 023526 (2019) [arxiv:1809.03924].
  * [31] A. Mack, T. Kahniashvili, and A. Kosowsky, Microwave Background Signatures of a Primordial Stochastic Magnetic Field, Phys. Rev. D 65, 123004 (2002) [astro-ph/0105504].
  * [32] P. A. R. Ade et al. (Planck Collaboration), Planck 2015 results. XIX. Constraints on primordial magnetic fields, Astron. Astrophys. 594, A19 (2016) [arxiv:1502.01594].