Analog Nature of Quantum Adiabatic Unstructured Search

TL;DR

This paper investigates the robustness of the quantum adiabatic unstructured search algorithm, revealing that its quadratic speedup is generally fragile against various imperfections such as control limitations, Hamiltonian errors, and thermal interactions.

Contribution

It provides a comprehensive analysis of how different non-idealities affect the performance and speedup of the quantum adiabatic unstructured search algorithm.

Findings

Quadratic speedup is not robust to control imperfections.

Hamiltonian misspecification undermines the algorithm's advantage.

Thermal interactions can eliminate the quantum speedup.

Abstract

The quantum adiabatic unstructured search algorithm is one of only a handful of quantum adiabatic optimization algorithms to exhibit provable speedups over their classical counterparts. With no fault tolerance theorems to guarantee the resilience of such algorithms against errors, understanding the impact of imperfections on their performance is of both scientific and practical significance. We study the robustness of the algorithm against various types of imperfections: limited control over the interpolating schedule, Hamiltonian misspecification, and interactions with a thermal environment. We find that the unstructured search algorithm's quadratic speedup is generally not robust to the presence of any one of the above non-idealities, and in some cases we find that it imposes unrealistic conditions on how the strength of these noise sources must scale to maintain the quadratic speedup.

Equations

$$H_p=\text{openone}-|m⟩⟨m|,$$

$$H(s(t))=(1-s(t))H_b+s(t)H_p,$$

$$|{\epsilon }_0(s)⟩$$

$$|{\epsilon }_1(s)⟩$$

$$\delta (s)$$

$$\cos \theta (s)$$

$$\sin \theta (s)$$

$$|{\epsilon }_{k+1}(s)⟩=\frac{1}{\sqrt{2}}\left(|f(k)⟩-|\overline{f(k)}⟩\right),1\le k\le \frac{N}{2}-1,$$

$$|{\epsilon }_2^{\prime }(s)⟩=\sqrt{\frac{N-2}{N-1}}(|\overline{m}⟩$$

$$-\frac{1}{N-2}\sum _{j=1}^{N/2-1}\left(|f(j)⟩+|\overline{f(j)}⟩\right)),$$

$$|{\epsilon }_{k+1}^{\prime }(s)⟩=\sqrt{\frac{2(k-1)}{k}}(\frac{1}{2}\left(|f(k)⟩+|\overline{f(k)}⟩\right)$$

$$-\frac{1}{2(k-1)}\sum _{j=1}^{k-1}\left(|f(j)⟩+|\overline{f(j)}⟩\right)),2\le k\le \frac{N}{2}-1,$$

$$f(j)=\left\{\begin{array}[]{@{}ll@{}}j-1,&\text{if}\ j-1<\min(m,\overline{m})\\ j,&\text{otherwise}\end{array}\right..$$

$$|⟨{\epsilon }_0(s)|{\sigma }_i^z|{\epsilon }_{k+1}(s)⟩|$$

$$1\le k\le \frac{N}{2}-1,$$

$$|⟨{\epsilon }_0(s)|{\sigma }_i^z|{\epsilon }_2^{\prime }(s)⟩|$$

$$|⟨{\epsilon }_0(s)|{\sigma }_i^z|{\epsilon }_{k+1}^{\prime }(s)⟩|$$

$$\frac{\text{d}s}{\text{d}t}=ϵ{\delta }^2(s),$$

$$s(t)=\frac{1}{2}+\frac{1}{2\sqrt{N-1}}\tan \left[\left(2\frac{t}{\mathcal{T}}-1\right){\tan }^{-1}\sqrt{N-1}\right],$$

$$\mathcal{T}=\frac{N}{ϵ\sqrt{N-1}}{\tan }^{-1}\sqrt{N-1}\approx \frac{\pi }{2ϵ}\sqrt{N},$$

$$H(s)=(1-s)H_b+s(1+\chi )H_p$$

$$s_{\ast }=\frac{N(\chi +2)-2(\chi +1)}{N(\chi +2{)}^2-4(\chi +1)},$$

$$P_s\approx \{\begin{array}{rl}\exp (-2.11N{\sigma }^2),& \text{mbox}for\sigma <\frac{1}{\sqrt{7N}}\\ \frac{1}{N},& \text{mbox}for\sigma >\sqrt{\frac{3}{N}}\end{array}.$$

$$R(s)$$

$$P_{s}\big{|}_{\mathrm{thermal}}=\frac{1}{1+e^{-\beta\delta}+(N-2)e^{-\beta(1+\delta)/2}}\ .$$

Full text

I Introduction

Quantum adiabatic optimization (QAO) Finnila et al. (1994); Brooke et al. (1999); Kadowaki and Nishimori (1998); Farhi et al. (2001); Santoro et al. (2002); van Dam et al. (2001); Reichardt (2004) is a paradigm of computing in which a slowly time-evolving Hamiltonian that uses continuously decreasing quantum fluctuations is employed in order to prepare the ground state of a target Hamiltonian in an analog, rather than digital, manner Young et al. (2008, 2010); Hen and Young (2011); Hen (2012); Farhi et al. (2012); Schützhold and Schaller (2006); Žnidarič and Horvat (2006). As such, it is expected by many to be a simpler way of carrying out quantum-assisted calculations experimentally Vandersypen et al. (2001); Bian et al. (2013); Babbush et al. (2014); Rønnow et al. (2014); Johnson et al. (2011); Hen et al. (2015); Boixo et al. (2013); Albash et al. (2015); Johnson et al. (2011); Tolpygo et al. (2015a, b); Jin et al. (2015); Lanting et al. (2017); Perdomo-Ortiz et al. (2011).

To date, there is only a handful of quantum adiabatic optimization algorithms whose runtime is provably superior to their classical counterparts Roland and Cerf (2002); Hen (2014a, b, c); Somma et al. (2012)111We exclude here algorithms derived via the polynomial equivalence theorem between the circuit and adiabatic models of quantum computing Aharonov et al. (2007) since the ‘final’ Hamiltonians are not necessarily diagonal in the computational basis.. First and foremost among these is the quantum adiabatic unstructured search (QAUS) algorithm — an oracular algorithm for identifying a marked state in an unstructured list. Originally devised by Roland and Cerf Roland and Cerf (2002) (but see also Refs. Farhi and Gutmann (1996); van Dam et al. (2001) for earlier variants), the algorithm consists of encoding the search space in a ‘problem Hamiltonian,’ $H_{p}$, that is constant across the entire search space except for one ‘marked’ configuration $|m\rangle=|m_{1}m_{2}\ldots m_{n}\rangle$ whose cost is lower than the rest. Here, $m_{i}\in\{0,1\}$ are the bits of the $n$-bit marked configuration (the number of elements in the search space is thus $N=2^{n}$ where $n$ is the number of elements). Similar to its gate-based counterpart, Grover’s unstructured search algorithm Grover (1997), the runtime of the Roland and Cerf algorithm scales as $O(\sqrt{N})$, which is to be contrasted with the linear scaling with $N$ of the number of queries required classically for finding the marked item.

While the asymptotic scaling of the runtime of quantum adiabatic algorithms such as QAUS give an accounting of the ‘time resources’ used by the algorithm, one should be careful about not accounting for other resources, particularly precision, due to the analog nature of the algorithm Aaronson (2005); Vergis et al. (1986); Jackson (1960). Failing to do so has practical ramifications since any physical implementation of an analog algorithm is expected to have some fixed precision.

In this study, we examine the robustness of the QAUS algorithm to finite precision as exhibited by several noise models Roland and Cerf (2005); Mandrà et al. (2015). For completeness we also revisit the thermal robustness of the algorithm Wild et al. (2016); Amin et al. (2009); Tiersch and Schützhold (2007); Åberg et al. (2005a, b); de Vega et al. (2010); Mostame et al. (2010) using a specific decoherence model. While these forms of imperfection are expected to appear together in any physical implementation, we treat each type separately here. We find that the quadratic speedup of the QAUS algorithm is sensitive to both finite precision and thermal effects, requiring both precision and temperature to scale in physically unreasonable ways to maintain the quantum speedup. For the former, we do this using two forms of Hamiltonian implementation errors that shift the position of the minimum gap, and only with a precision that scales exponentially with the system size can the quadratic speedup be maintained.

While it is well accepted that scalable quantum computing is not possible without fault tolerance Aharonov et al. (2006), there is as of yet no known accuracy-threshold theorem for the adiabatic paradigm of quantum computing. Therefore, while fault tolerance schemes can be applied to the gate-based approach for solving unstructured search Grover (1997), no equivalent schemes exist to date for the adiabatic approach. Our study therefore calls into question the practical significance of the QAUS asymptotic speedup in the absence of physically meaningful schemes to mitigate and correct for these errors. Specifically, if we are to rely on such speedups to give rise to a significant separation between the computational costs of quantum and classical algorithms at some maximum size, there is a significant engineering challenge to realize the necessary quantum system with a sufficiently high precision, a feat that becomes increasingly harder with growing system size.

We begin with a brief overview of the algorithm and then move on to discuss the various types of imperfections considered and their impact on performance. In the concluding section we discuss the meaning and implications of our results.

II The quantum adiabatic unstructured search algorithm

The unstructured search problem Hamiltonian is a one-dimensional projection onto the marked state:

[TABLE]

where $|m\rangle\langle m|$ is the projection onto the marked state, which belongs to the computational basis $\{|0\rangle,|1\rangle,\dots,|N-1\rangle\}$, and $\openone$ is the identity operator. To achieve the quadratic speedup, a carefully tailored variable-rate annealing schedule $s(t)$ is chosen that interpolates the Hamiltonian between a ‘beginning’ Hamiltonian $H_{b}$ and the problem Hamiltonian $H_{p}$, varying slowly as a function of time $t\in[0,\mathcal{T}]$ in the vicinity of the minimum energy gap between the ground state and first excited state and more rapidly in places where the energy gap is large Jansen et al. (2007); Lidar et al. (2009); Kato (1950); Roland and Cerf (2002). Here, $H_{b}$ is a one-dimensional projection onto the equal superposition of all computational basis states, i.e., $H_{b}=\openone-|+\rangle\langle+|$, where $|+\rangle=\bigotimes_{i=1}^{n}|+\rangle_{i}$, and $|+\rangle_{i}=\frac{1}{\sqrt{2}}\left(|0\rangle_{i}+|1\rangle_{i}\right)$, and the total Hamiltonian is given by

[TABLE]

where we have assumed the boundary conditions $s(0)=0$ and $s(\mathcal{T})=1$ at the beginning and end of the interpolation respectively. While the efficiency of a generic adiabatic algorithm may depend sensitively on the form of $H_{b}$ Žnidarič and Horvat (2006); Farhi et al. (2008), the one-dimensional projection above gives rise to the optimal scaling performance Roland and Cerf (2002).

If the initial state is taken to be the ground state of $H(0)$, i.e., $|+\rangle$, then the evolution according to $H(s)$ is restricted to the two-dimensional subspace spanned by $|m\rangle$ and $|m^{\perp}\rangle=\frac{1}{\sqrt{N-1}}\sum_{i\neq m}^{N}|i\rangle$. The ground state and first excited state of the system are in this subspace throughout the interpolation and can be written as:

[TABLE]

with eigenvalues $\varepsilon_{0}(s)=\frac{1}{2}\left(1-\delta(s)\right)\ ,\quad\varepsilon_{1}(s)=\frac{1}{2}\left(1+\delta(s)\right)$ respectively and

[TABLE]

The remaining $N-2$ energy eigenstates are outside the aforementioned two-dimensional subspace and have energy $1$ throughout the interpolation. For later convenience, we write them as:

[TABLE]

where $\overline{f(j)}=N-1-f(j)$ is the integer associated with the negation of the bit-representation of the integer $f(j)$ and

[TABLE]

This particular form of the excited states is useful because, $\sigma_{i}^{z}\left(|f\rangle+|\overline{f}\rangle\right)=\pm\left(|f\rangle-|\overline{f}\rangle\right)$ and $\sigma_{i}^{z}\left(|f\rangle-|\overline{f}\rangle\right)=\pm\left(|f\rangle+|\overline{f}\rangle\right)$, irrespective of the qubit index $i$ and the state $|f\rangle$. This then means that we have the following relations:

[TABLE]

The optimized annealing schedule of Roland and Cerf Roland and Cerf (2002) that defines the QAUS algorithm satisfies a ‘local’ adiabatic condition Hen (2014a, b):

[TABLE]

where $\epsilon$ is a small constant. The optimized annealing schedule satisfying the interpolation boundary conditions is given by

[TABLE]

with the optimal runtime being Roland and Cerf (2002)

[TABLE]

i.e., it is proportional to the square root of the dimension of the Hilbert space, similarly to its gate-based counterpart Grover (1997). For a sufficiently small $\epsilon$, this choice guarantees that a system prepared in the ground state of $H(0)$ remains close to the instantaneous ground state throughout the evolution using $H(s)$.

III Finite schedule precision

The QAUS algorithm is analog in nature, in that it requires continuously varying the strengths of $H_{b}$ and $H_{p}$ throughout the evolution Roland and Cerf (2002, 2003); Hen (2014a, b). For the local adiabatic interpolation, Eq. (7), the annealing schedule $s(t)$ changes exponentially slowly around the minimum gap, which is on the order of $1/\sqrt{N}$, in a region of width $1/\sqrt{N}$ Childs and Goldstone (2004); Hen (2017). In any conceivable physical setting however, we expect only a limited control over the interpolating schedule, and here we ask whether this restriction adversely affect the performance of the QAUS algorithm.

We begin our exploration by specifying the schedule $s(t)$ as a piecewise linear schedule between equally spaced time points $0,t_{1},t_{2},\ldots,\mathcal{T}$ with $t_{j}=j\Delta t$ for different spacings $\Delta t$ such that the schedule at $s(t_{j})$ coincides with the original QAUS schedule given by Eq. (7) [see Fig. 1]. A numerical investigation reveals that a piecewise linear schedule with only two intermediate points (3-piece interpolation) suffices to achieve the quadratic speedup. This is demonstrated in Fig. 1, which depicts the probability of success $P_{s}$, the probability of measuring the marked state at the end of the evolution, as a function of problem size $n$ for three- and four-piece interpolations. The results show that already with a 3-piece schedule and a total time given by Eq. (9), a constant (with system size) probability of success is achieved. Higher-piece interpolations give, as expected, higher success probabilities.

We thus find that the smooth $s(t)$ schedule in Eq. (8) is not necessary to obtain a quadratic speedup for as long as the linear slope at the minimum gap, $s=1/2$, scales as $1/\sqrt{N}$. Since the region of the minimum gap shrinks exponentially as $1/\sqrt{N}$, this requires ‘hitting’ the location of the minimum gap with increasing precision as the problem size grows.

To illustrate the above point, we consider the scenario of a slightly shifted schedule $s(t)$ that ‘misses’ the location of the minimum gap by a small but fixed amount. This is equivalent to the case where the Hamiltonian itself is slightly misspecified:

[TABLE]

where $\chi$ is a small fixed constant and the schedule $s(t)$ is taken to be the unperturbed one [Eq. (8)]. For the above Hamiltonian, the gap is minimal at

[TABLE]

which, in the limit of $N\to\infty$ becomes $s_{*}=(\chi+2)^{-1}$. By employing the original QAUS annealing schedule, it is easy to see that there will always be a problem size $n_{*}$ beyond which the schedule does not sufficiently slow-down in the vicinity of the minimum gap. We confirm these expectations in Fig. 2 with simulation results for different values of displacements $\chi$ corresponding to displaced minimal gaps. Any nonzero value of $\chi$ (equivalently, any nonzero displacement of the minimum gap) eventually leads to an exponentially decreasing probability of success, with the transition to exponential behavior occurring at larger values of $n$ for smaller displacements $\chi$.

To have the quadratic speedup, a fixed success probability must be maintained for growing system size. The results in Fig. 2 show that to achieve this, the distance of $s_{*}$ from $1/2$ must be decreased accordingly. We can ask how big a perturbation is allowed, or equivalently how many bits of precision are required, for the schedule to achieve this. Since the gap is small to within a width of $1/\sqrt{N}$, we expect to require approximately $n/2$ bits of precision. Thus the schedule must be precise to within $O(n)$ bits of precision in order to maintain the quadratic speedup. This is confirmed by the numerical data in Fig. 2, where we see that for $n=11,13,15,17$, we require approximately a factor of 2 decrease in the distance of $s_{*}$ from $1/2$. We further discuss the feasibility of the increasing precision requirement in the concluding section.

IV Noisy Hamiltonian

The noise model in the previous section still restricted the unitary evolution to the two-dimensional subspace spanned by $|m\rangle$ and $|m^{\perp}\rangle$. We now extend our analysis by considering noise that prevents the evolution from being restricted to this subspace. Specially, we consider the QAUS algorithm perturbed by a noise Hamiltonian $\tilde{H}$ such that the total Hamiltonian is now given by $H^{\prime}(s)=H(s)+\tilde{H}$, where the noise Hamiltonian $\tilde{H}$ has matrix elements in the $|0\rangle,\dots,|N-1\rangle$ basis that are drawn randomly from a Gaussian distribution with mean zero and standard deviation $\sigma$. Our model of noise has the elements fixed throughout the evolution, which is different from the time-dependent noise model studied in Ref. Roland and Cerf (2005). The adaptive step Runge-Kutta-Fehlberg algorithm was used for the efficient numerical solution of the time dependent Schrödinger equation Press et al. (1992); Cash and Karp (1990).

Figure 3 shows the dependence of the probability of success $P_{s}$ on $\sigma$ for various $N$ values. The data can be fitted by

[TABLE]

The probability of success approaches $1/N$ in the large noise limit. In this limit the Hamiltonian is random so measuring the marked state occurs with probability $1/N$. The initial exponential decay of $P_{s}$ is a function of $N\sigma^{2}$. This means that for a constant noise strength $\sigma$, the probability of success decays exponentially with $N$, and the only way to mitigate it is to require that $\sigma$, the noise strength, scales as $1/\sqrt{N}$.

V Interaction with a thermal bath

The robustness of the QAUS algorithm in the presence of interactions with an external environment has been extensively studied Wild et al. (2016); Amin et al. (2009); Tiersch and Schützhold (2007); Åberg et al. (2005a, b); de Vega et al. (2010). A generic interaction breaks the symmetry that restricts the system evolution to the lowest two eigenstates, and for completeness here we show how the exponential number of excited states within a constant energy gap places (unrealistic) requirements on the temperature (or overall energy scale of the Hamiltonian) to maintain performance even for possibly the most innocuous noise model Albash and Lidar (2015).

We consider a model of decoherence between a quantum annealing system of qubits and a thermal environment described by the Markovian adiabatic master equation Albash et al. (2012). (We assume that this model holds throughout the anneal, even though we expect the validity conditions of the microscopic derivation of the model to break down near the minimum gap.) We focus on the case where each qubit is connected to its own independent bath of bosonic harmonic oscillators. The excitation rate from the ground state to an excited state $|\varepsilon_{i}(s)\rangle$ at any point $s$ is generically given by $R_{0\to i}(s)=\sum_{\alpha=1}^{n}\gamma(\Delta_{i})e^{-\beta\Delta_{i}(s)}|\langle\varepsilon_{0}(s)|A_{\alpha}|\varepsilon_{i}(s)\rangle|^{2}$, where $\Delta_{i}(s)$ is the energy gap from the ground state to the excited state and $\beta$ is the inverse-temperature of the bosonic bath. $\gamma(\Delta_{i})$ encodes the spectral density of the bosonic bath, the bath correlations, and the system-bath coupling strength $g$, and $A_{\alpha}$ is the system operator part of the system-bath interaction. We consider $A_{\alpha}=\sigma_{\alpha}^{z}$ corresponding to a ‘dephasing’ bath. For concreteness, we can consider an Ohmic spectral density, such that $\gamma(\Delta)=2\pi g^{2}\Delta/(1-e^{-\beta\Delta})$ Albash et al. (2012).

Of relevance to us is the excitation rate during the anneal from the instantaneous ground state to all the excited states outside the two-dimensional subspace, which is given by

[TABLE]

where we have used the relations in Eq. (6).

If follows from Eqs. (4) that $\sin^{2}\frac{\theta}{2}$ is a monotonically decreasing function of $s$ for $0\leq s\leq\frac{1}{2}$ and $N>1$, such that $\sin^{2}\frac{\theta}{2}=1-\frac{1}{N}$ for $s=0$ and $\sin^{2}\frac{\theta}{2}=\frac{1}{2}-\frac{1}{2\sqrt{N}}$ for $s=\frac{1}{2}$. Hence, we expect that the excitation rate to the excited states outside the two-dimensional subspace for the first half of the anneal to scale as $\sim ng^{2}/(e^{\beta(1+\delta)/2}-1)$ for large $n$. In conjunction with a total annealing time that scales as $\sqrt{N}$ (Eq. (9)), we can expect the open-system dynamics to not be restricted to the two-dimensional subspace for a constant temperature and system-bath coupling.

We note that if the system thermalizes on the instantaneous Hamiltonian, the probability of success at any point in the anneal is given by

[TABLE]

For any fixed nonzero temperature, this gives a probability of being in the instantaneous ground state that scales as $1/N$ for any point along the interpolation.

In order to suppress the excitations to outside the two-dimensional subspace during the anneal, we can scale the inverse temperature $\beta$ linearly with $n$ for a constant $g$, which will exponentially (in $n$) suppress thermal excitations out of the ground state and will ensure that the instantaneous thermal state always has a finite population on the instantaneous ground state. These results are consistent with the analysis of Ref. Wild et al. (2016), although in that work the overall energy scale of the Hamiltonian $E_{0}$ was scaled linearly, such that $\beta E_{0}\sim n$. Alternatively the system-bath coupling $g$ can be scaled down at least as $N^{-1/4}$ for a constant $\beta$ in order to ensure that thermal excitations are suppressed during the entire evolution.

VI Conclusions and discussion

We studied the robustness of the quantum adiabatic unstructured search algorithm against various types of imperfections from limited control over the adiabatic schedule to Hamiltonian misspecification to an interaction with a decohering bath. Our findings can be summarized as follows. In the presence of finite perturbations to the Hamiltonian, the probability of hitting the marked state decreases exponentially with system size $n$ if the interpolating schedule is not adjusted accordingly. This results in the loss of the quadratic speedup of the error-free algorithm. The scaling is similar when we consider a noise model that introduces Gaussian noise to the matrix elements, which now does not restrict the evolution to a two-dimensional subspace: the probability of hitting the marked state now decreases exponentially in $N$ for a fixed standard deviation of the noise. Our results indicate that the standard deviation must be scaled down as $1/\sqrt{N}$, which can also be derived from the analysis of Ref. Jarret et al. (2018)222We thank Michael Jarret for pointing this out.. While neither of the above noise models have been constructed with a physical mechanism in mind, these noise models reproduce effects we expect to generically occur. We expect generic noise to break the symmetries of the Hamiltonian that restrict the evolution to a particular subspace, and we expect generic noise to shift the position of the minimum gap in a noise-instance dependent way. Our results show that without the interpolation schedule slowing down precisely at the noise-shifted minimum gap, the quadratic speedup of the QAUS algorithm will be lost.

We emphasize that even if the Hamiltonians $H_{b}$ and $H_{p}$ can be implemented precisely, the annealing schedule $s(t)$ still needs to be controlled with exponential precision around the minimum gap, even if we use piece-wise linear interpolations. This need for growing precision must inevitably translate to the need of additional resources, without which the QAUS algorithm cannot retain its quadratic speedup. This is the signature of analog computing, and our results illustrate the need for alternative methods that would combat the exponentially growing precision requirement.

Our work also has some implications for algorithms that require access to continuous-query Hamiltonian oracles or query other properties of the Hamiltonian (e.g., Refs. Jarret et al. (2018); Cleve et al. (2009)) wherein a sub-routine returning, e.g., the value of the gap, is called. Our work suggests that the value of the gap needs to be returned with growing precision as a function of system size and hence requires growing space resources that needs to be accounted for.

We also revisited the thermal stability of the algorithm, studying it in the framework of the weak-coupling Markovian adiabatic master equation. Here, in the absence of specific fine tuning, the presence of an exponential number of excited states at a fixed energy gap away from the ground state already imposes serious constraints on the temperature and/or the system-bath coupling just to ensure the evolution is restricted to the two-dimensional subspace. The former needs to be scaled down at least inversely proportional to the system size, or the latter must be scaled down at least exponentially with the system size.

We finally point out that our analyses above is an asymptotic one. Since any physical device will have a finite fixed size, one could imagine noise strengths that are sufficiently reduced to make the QAUS algorithm successful. Such a device may still have practical uses if the computational costs of the quantum and classical algorithms are well seperated, and our results do not exclude such a possibility.

Acknowledgements.

We thank Michael Jarret for useful discussions. The research is based upon work (partially) supported by the Office of the Director of National Intelligence (ODNI), Intelligence Advanced Research Projects Activity (IARPA), via the U.S. Army Research Office contract W911NF-17-C-0050. The views and conclusions contained herein are those of the authors and should not be interpreted as necessarily representing the official policies or endorsements, either expressed or implied, of the ODNI, IARPA, or the U.S. Government. The U.S. Government is authorized to reproduce and distribute reprints for Governmental purposes notwithstanding any copyright annotation thereon. L. B. was supported within the framework of State Assignment No. 0033-2019-0007 of Russian Ministry of Science and Higher Education.