# Pattern generation and symbolic dynamics in a nanocontact vortex oscillator

Myoung-Woo Yoo, Damien Rontani, Jérémy Létang, S. Petit-Watelot, Thibaut Devolder, Marc Sciamanna, Karim Bouzehouane, Vincent Cros, Joo-Von Kim

► To cite this version:

Myoung-Woo Yoo, Damien Rontani, Jérémy Létang, S. Petit-Watelot, Thibaut Devolder, et al. Pattern generation and symbolic dynamics in a nanocontact vortex oscillator. Nature Communications, 2020, 11 (1), pp.601. 10.1038/s41467-020-14328-7. hal-03005535

HAL Id: hal-03005535

https://hal.science/hal-03005535

Submitted on 14 Nov 2020

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

ARTICLE

https://doi.org/10.1038/s41467-020-14328-7

OPEN

# Pattern generation and symbolic dynamics in a nanocontact vortex oscillator

Myoung-Woo Yoo ID 1*, Damien Rontani ID 2, Jérémy Létang 1, Sébastien Petit-Watelot ID 3, Thibaut Devolder 1, Marc Sciamanna 2, Karim Bouzehouane ID 4, Vincent Cros ID 4 & Joo-Von Kim ID 1*

Harnessing chaos or intrinsic nonlinear behaviours of dynamical systems is a promising avenue toward unconventional information processing technologies. In this light, spintronic devices are promising because of the inherent nonlinearity of magnetization dynamics. Here, we demonstrate experimentally the potential for chaos-based schemes using nanocontact vortex oscillators by unveiling and characterizing their waveform patterns and symbolic dynamics using time-resolved electrical measurements. We dissociate nonlinear deterministic patterns from thermal fluctuations and show that the emergence of chaos results in the unpredictable alternation of well-defined patterns. With phase-space reconstruction techniques, we perform symbolic analyses of the time series and show that the oscillator exhibits maximal entropy and complexity at the centre of its incommensurate region. This suggests that such vortex-based systems are promising nanoscale sources of entropy that could be exploited for information processing.

<footnote n="1">¹Centre de Nanosciences et de Nanotechnologies, CNRS, Université Paris-Saclay, 10 boulevard Thomas Gobert, 91120 Palaiseau, France. ²Chaire Photonique, LMOPS EA 4423 Laboratory, Université de Lorraine & CentraleSupélec, 2 rue Edouard Belin, F-57070 Metz, France. ³Institut Jean Lamour, CNRS, Université de Lorraine, Campus Artem, 2 allée André Guinier, 54011 Nancy, France. ⁴Unité Mixte de Physique, CNRS, Thales, Université Paris-Saclay, 1 avenue Augustin Fresnel, 91767 Palaiseau, France. *email: myoung-woo.yoo@c2n.upsaclay.fr; joo-von.kim@c2n.upsaclay.fr</footnote>
Nonlinear dynamics and chaos are powerful frameworks with which many phenomena in physics, biology, and engineering can be understood¹. Chaos refers to the high sensitivity of a nonlinear dynamical system to perturbations in its initial conditions, where the temporal evolution is unpredictable on the long term. From the perspective of applications in information processing, chaos has attracted much attention over the last two decades² because chaotic waveforms are random-like, yet deterministic and potentially controllable. They have found various applications in information technologies, such as encrypted communications at the physical layer³, ultrafast random number generation⁴,⁵, data processing⁶, computing⁷, secure-key exchange⁸, radar applications⁹, precision sensing¹⁰, and encoding information via symbolic dynamics¹¹.

Spintronic devices based on magnetic multilayers are good candidates for chaos-based applications, because magnetization dynamics in ferromagnets is intrinsically nonlinear¹²⁻²⁰. Moreover, such dynamics can be driven and detected by spin-dependent transport phenomena, such as spin-transfer torques, magnetoresistive effects, and (inverse) spin-Hall effects²¹,²², giving rise to devices such as spin-torque and spin-Hall nanos oscillators²³,²⁴ that can be integrated into conventional semiconductor electronics²⁵. Because magnetization dynamics can occur at the nanoscale at microwave frequencies, spintronic devices hold much promise for highly compact, GHz-rate information processing using chaos.

One example of chaos in a nanoscale spintronic device can be found in nanocontact vortex oscillators (NCVOs)²⁶,²⁷. In the NCVO, the gyration and switching of the vortex core can be induced by spin-transfer torques and oscillating output signals can be detected by the magnetoresistance. In contrast to vortex oscillators based on nanopillars²⁸⁻³⁰, the NCVO can exhibit nontrivial dynamics that involves a self-phase-locking phenomenon between the core gyration and core switching²⁶. If the ratio between the frequencies of these two processes is irrational, the behaviour is chaotic²⁷.

Here, we demonstrate experimentally that the chaotic regime of the NCVO involves simple aperiodic waveform patterns. These can be encoded into bit sequences, which are correlated with the core-polarity state of the magnetic vortex. First, we describe time-resolved signals from the NCVO at 77 K and validate their chaotic characteristics from sensitivity to initial conditions and correlation dimension analysis. Then, we show that the time traces are in fact only composed by a few waveform patterns which are ordered aperiodically in the chaotic regime. By reconstructing attractor geometries from the measured time series, we reveal the symbolic dynamics of chaotic NCVOs, which is in good agreement with the patterns observed in simulation. We extract bit sequences based on this symbolic analysis and show that the generated bits can achieve maximal values of the Shannon block entropy and Lempel–Ziv complexity.

## Results

**Chaos in an NCVO.** The NCVO comprises an extended spinvalve multilayer with a metallic point contact (approximately 20 nm in diameter) on the top of the surface (Fig. 1a)²⁶. When an electric current is applied through the contact, the component of the current flow perpendicular to the film generates an Oersted field (blue arrow in Fig. 1a) which promotes a magnetic vortex in the free layer and generates a Zeeman energy potential for it that is centred on the nanocontact³¹. The current component in the film plane (orange arrows in Fig. 1a) pushes the vortex core out from the centre by exerting spin-transfer torques³². The competition between the two effects results in a stable gyration²⁶,³¹,³³,³⁴ and switching dynamics²⁶ of the vortex core around the

nanocontact. In general, the shape of the core trajectory around the nanocontact is not circular³⁴, as evidenced by a rich harmonic content in the power spectrum. This results from the presence of an antivortex or domain walls in the extended film that appear during the nucleation process²⁶.

Figure 1b shows a map of the power spectral density of the magnetoresistance oscillations at 77 K as a function of the applied current $I_{dc}$, measured with a spectrum analyser. The NCVO exhibits three dynamical regimes: pure-gyration, commensurate, and incommensurate states²⁶. When $I_{dc}$ is lower than a threshold for core reversal (~10.3 mA here), only the gyration frequency, $f_0$, is observed (red circles in Fig. 1b) with its harmonics ($2f_0, 3f_0, ...$) because of the noncircular core trajectories. If $I_{dc}$ is larger than the threshold, core reversal appears in addition to the gyration. This dynamical state is accompanied by additional sidebands at $f_0 \pm f_{mod}$ (yellow crosses in Fig. 1b), where $f_{mod}$ is a modulation frequency that is related to the periodicity of the core reversal. The ratio $f_{mod}/f_0$ as a function of $I_{dc}$ is shown in Fig. 1c, which is distinguished by two plateaus with monotonic increases elsewhere. At the plateaus (yellow regions in Fig. 1c), the frequency spectrum shows clear peaks and $f_{mod}/f_0$ remains constant with $I_{dc}$. In this case, the ratio can be expressed as simple integer fractions (1/3 and 1/4 in this experiment), because the core reversal process is phase-locked to the gyration²⁶. The relation between the core dynamics and $f_{mod}/f_0$ is shown in Supplementary Fig. 1 and discussed in more detail in Supplementary Note 1. In the incommensurate state (red regions in Fig. 1c), in contrast, the frequency spectrum becomes more complex, where $f_{mod}/f_0$ varies with $I_{dc}$ from one plateau to another, which indicates that no phase locking occurs between the gyration and core reversal.

The incommensurate state represents chaotic behaviour²⁷. Sensitivity to initial conditions, a hallmark of chaos, can be seen in the time traces. These were obtained at 77 K with a single-shot oscilloscope, which were then filtered using a pattern matching technique to reduce the measurement noise (Supplementary Fig. 2, Supplementary Note 2). By overlaying several segments of the time traces with very similar initial conditions, we can obtain a visual measure of the sensitivity to initial conditions in the commensurate (Fig. 1d) and incommensurate states (Fig. 1e). In the commensurate state (Fig. 1d), the waveforms remain coherent over tens of nanoseconds, with evidence of jitter setting in at around 50 ns. In the incommensurate state, however, the coherence is lost below 10 ns (Fig. 1e), which is due to the sensitivity to initial conditions. We can further verify the presence of chaos by analysing the fractal geometry of the reconstructed attractor³⁵. To this end, we compute the correlation dimension $D_c$ from the filtered time series using the correlation sum $C_m(\epsilon)$ in Fig. 1f and its derivative with respect to a geometric scaling $\epsilon$ (see Methods). We estimate the geometric dimension in Fig. 1g by looking at constant values of the derivative. For the commensurate case at $I_{dc} = 14$ mA (yellow lines in Fig. 1g) the dimension is found to be $D_c = 1.05 \pm 0.01$, which is very close to 1 and consistent with limit-cycle dynamics of a NCVO presenting a small amount of jitter in the position of the core-polarity switching. For the commensurate case at $I_{dc} = 12.6$ mA (not shown in Fig. 1f-g for clarity), the correlation dimension is $D_c = 1.15 \pm 0.01$, which is further from 1 compared with $I_{dc} = 14$ mA, because of a higher level of jitter in the position of the core-polarization switching. For the incommensurate case (red lines in Fig. 1g), however, the dimension is found to be $D_c = 1.70 \pm 0.04$, which is consistent with a fractal geometry associated with temporal chaos. In Fig. 1h, we present the analysis based on correlation sum to determine a lower bound of the Kolmogorov–Sinai (KS) entropy $h_{KS}$ by computing the $K_2$-entropy metric (see Methods section). In Fig. 1i, by
<figure><img src="image_3.png" /><figcaption>**Fig. 1** Chaotic characteristics of the output time traces. **a** A schematic of a nanocontact vortex oscillator. **b** A map of the power spectra as a function of input current amplitudes, *I*<sub>dc</sub>. The red circles and yellow cross marks indicate a fundamental frequency, *f*<sub>0</sub>, and its upper sideband at *f*<sub>0</sub> + *f*<sub>mod</sub>, respectively. *f*<sub>mod</sub> is a modulation frequency. **c** *f*<sub>mod</sub>/*f*<sub>0</sub> as a function of *I*<sub>dc</sub>. The yellow and red regions represent the commensurate and the incommensurate states, respectively. The dotted horizontal lines indicate plateaus in which the self-phase-locking occurs. **d** Eighteen different time traces which have identical initial conditions in the commensurate state (*I*<sub>dc</sub> = 14.0 mA). **e** As in **d** but in the incommensurate state at *I*<sub>dc</sub> = 13.2 mA. **f** Correlation integrals *C*<sub>*m*</sub>(*e*) as a function of a geometric scaling *ϵ* and **g** its derivatives of ∂ln*C*<sub>*m*</sub>(*e*)/∂ln*ϵ* for embedding dimensions from *m* = 10 (bottom curve) to *m* = 24 (top curve) by increment of 2. The red and yellow curves are for *I*<sub>dc</sub> = 13.2 and 14.0 mA, respectively. The horizontal lines indicated by red and yellow arrows correspond to estimate of *D*<sub>*c*</sub> based on the flats in the scaling *ϵ*-interval [1, 10<sup>0.3</sup>]. **h** Metric K<sub>2,m</sub>(*ϵ*) as a function of *ϵ* for embedding dimensions from *m* = 10 (top curve) to *m* = 24 (bottom curve) by increment of 2 at *I*<sub>dc</sub> = 13.2 mA and **i** 〈K<sub>2,m</sub>(*ϵ*)〉<sub>*ϵ*</sub> in the scaling range [1, 10<sup>0.3</sup>] for the asymptotic determination of the K<sub>2</sub>-entropy.</figcaption></figure>

extrapolating the tendency of ⟨K'<sub>2,m</sub>(e)⟩<sub>e</sub> with increasing values of *m* with an exponential fit, we estimate that *K*₂ = 0.12 ± 0.02 ns⁻¹ for *I*<sub>dc</sub> = 13.2 mA. These results are consistent with a previous study using the titration of chaos with added noise³⁶ to identify the presence of chaos in the NCVO²⁷. Note that the responses in the incommensurate state can be reproduced by micromagnetic simulations even at 0 K as shown previously²⁶. We contend therefore that the measured chaotic characteristics are mostly deterministic, rather than stochastic as driven by thermal fluctuations.

**Pattern generation.** A feature of the chaos generated by the NCVO involves distinct waveforms that repeat aperiodically. Representative experimental time series are shown in Fig. 2a–c. Here we use a time axis that is normalized with respect to the core gyration period, 1/f₀, such that tf₀ represents the number of core gyrations. Both commensurate and incommensurate states show similar features where distinct oscillatory patterns are delimited by cusps. These patterns correspond to a number

<figure><img src="image_4.png" /></figure>

of orbits of the vortex core around the nanocontact, with the cusps representing a core reversal event; the position of these cusps are indicated by the dots and dotted lines in Fig. 2a–c. Note that core reversal results in the change in the sense of gyration (i.e. clockwise to counterclockwise, and vice versa). These features in the measured time series are reproduced in micromagnetic simulations (see Supplementary Fig. 3 and Supplementary Note 3).

From the intervals between the core switching events in Fig. 2a–c, we can define a gyration number for the core switching, *n* = Δtf₀, which is shown in Fig. 2d–f. In the commensurate state, *n* exhibits a simple time evolution. At *I*<sub>dc</sub> = 12.6 mA, the switching always occurs every two gyrations (Fig. 2d), so *n* remains constant at 2. Similarly, at *I*<sub>dc</sub> = 14.0 mA, core reversal occurs after one and two gyrations successively, a process which repeats periodically; in this case *n* oscillates between 1 and 2 as shown in Fig. 2f. In the incommensurate state (*I*<sub>dc</sub> = 13.2 mA), however, *n* switches between 1 and 2 in an aperiodic fashion (Fig. 2e), which is consistent with the chaotic dynamics expected in this regime.
<figure><img src="image_5.png" /><figcaption>**Fig. 2** Pattern generation from time series. **a–c** Representatives of experimentally measured time series at *I*<sub>dc</sub> = 12.6, 13.2, and 14.0 mA. The black dots with dotted lines indicate the core-polarity switching events. We normalized the time axes *tf*<sub>0</sub> which is identical with a number of core gyrations. The yellow, blue, and red regions indicate waveform patterns denoted by *pn* = −2, −1, and +2, respectively, where *p* is a core polarity and *n* is a required gyration number for the core-polarity switching. **d–f** Time evolutions of a required gyration number for the core switching, *n* = Δ*tf*<sub>0</sub>, obtained from **a–c** by calculating intervals between the black dots. **g** Schematics of the possible core-polarity switching scenarios. The core trajectories (top panels) and expected output waveforms (bottom panels) are shown for different *p* and *n* combinations. The colours of the oscillatory patterns (red, blue, and yellow) correspond with those in **a–c**.</figcaption></figure>

The required gyration number for core switching, *n*, is always approximately integer (typically 1 or 2 in this experiment) in both commensurate and incommensurate cases as shown in Fig. 2d–f. This is consistent with simulation results, in which core-polarity switching occurs only in a restricted region of the film plane close to the nanocontact, where conditions for core reversal are met<sup>26</sup>. In addition the core polarity, *p*, can only have two values, +1 and −1, so we hypothesize that in general there are only four possible patterns for the commensurate or incommensurate states, *pn* = +1, +2, −1, and −2. In other words, the nonlinear physical properties of the NCVO result in a sequence that represents a combination of these four patterns.

We plot schematic core trajectories of the possible switching scenarios in Fig. 2g, along with the expected time series for different *pn*. Without loss of generality, we assume here that the vortex has clockwise chirality and the reference layer is saturated in the +*y* direction. Based on the schematic waveforms, we can identify the corresponding oscillatory patterns from the time series, as indicated by the background colours in Fig. 2a–c. In the commensurate state (Fig. 2a, c), the time series are composed of two *pn* patterns, which repeat periodically. They involve *pn* = {−2, +2} for *I*<sub>dc</sub> = 12.6 mA (Fig. 2a) and {−1, +2} for 14.0 mA (Fig. 2c), respectively. In the incommensurate state (Fig. 2b), there exist three *pn* patterns, *pn* = {−2, −1, +2}, which appear without a well-defined periodicity. This shows

that the NCVO generates simple oscillatory patterns even in the chaotic state.

**Symbolic dynamics and bit sequences.** We further analyse the pattern generation from the perspective of symbolic dynamics<sup>37</sup>. The principle of symbolic dynamics is to find an adequate partition of the system's Poincaré section in its phase space (see Methods), such that every time there is a transition from one region of the section to another, a symbol is emitted. As a result, the nonlinear dynamics of the system can be reduced to a sequence of symbols. However, generating symbols with the proper partitions and surfaces from model-free, experimental, and scalar time series is a challenging problem in general. First, we reconstruct attractor geometries from the measured time series using a three-dimensional delay embedding, which is sufficiently large to completely unfold the reconstructed attractors when considering the commensurate and incommensurate states shown in Fig. 3a–c<sup>35</sup>. In the commensurate state, the NCVO exhibits a limit cycle, because the trace is periodic (Fig. 3a, c). On the contrary, the attractor for the incommensurate state is more intricate (Fig. 3b). As explained in the Methods section, we set a proper Poincaré surface (white surfaces in Fig. 3a–c) allowing for a potential generating partition. Note that the surface is not unique, and we can choose any different plane for the symbolic analysis. On the surface, we obtain the Poincaré maps from the
<figure><img src="image_6.png" /><figcaption>**Fig. 3** Reconstructed attractor geometries and symbolic dynamics. **a–c** Reconstructed attractor geometries by a method of delay from the measured time traces at $I_{dc}$ = 12.6, 13.2, and 14.0 mA, respectively. The white plane is arbitrarily chosen Poincaré surface of section (see Methods). (Inset) Poincaré maps at the surfaces. The red dashed lines indicate a simple partition to divide the plane into two regions, $R_A$ and $R_B$, for encoding symbols, A and B. The Poincaré surface of section and partitions are identical for all $I_{dc}$ in these figures. **d–f** Dynamics of symbols defined from the partition on the Poincaré maps. Above the graphs, corresponding $pn$ patterns and generated bit sequences are represented. The bits are defined as $0 \equiv [A,A,A,B,B,B]$ and $1 \equiv [A,B,B,B,B]$. **g–i** Rules of the symbolic dynamics at $I_{dc}$ = 12.6, 13.2, and 14.0 mA. **j** Generated bit sequences for long term at $I_{dc}$ = 12.6, 13.2, and 14.0 mA.</figcaption></figure>

intersection (insets of Fig. 3a–c) in which we find several clusters of the points. Note that here, for the Poincaré map, we do not consider the transverse orientation. We use a simple partition to divide the map into two different regions, $R_A$ and $R_B$ (red dashed lines in the insets of Fig. 3a–c), then record the symbols, A or B, when the attractor cross the surface either in $R_A$ or $R_B$, respectively. The encoded symbolic sequences are shown in Fig. 3d–f, where the corresponding $pn$ patterns are shown above the graphs. By comparing the $pn$ patterns and symbolic dynamics, we can see that $pn = -2, -1$, and +2 correspond to [A,A,A], [A], and [B,B, B], respectively. This result shows that the choice of the partition for the determination of the symbolic dynamics is in good agreement with the $pn$ sequences in both the commensurate and incommensurate states. Other partition choices are also possible but may render identifying the symbolic sequences more difficult (see Supplementary Fig. 4 and Supplementary Note 4).

We can find simple rules in the symbolic sequences (Fig. 3d–f). In the commensurate case, the sequences show only one repeated cycle: [A,A,A,B,B,B] and [A, B, B, B] for $I_{dc}$ = 12.6 and 14.0 mA, respectively (Fig. 3g, i). In the incommensurate state, however, two possible cycles coexist in the sequence and appears erratically over

time (Fig. 3h). To simplify the analysis of complexity, we define binary symbols attributed to the two different patterns accessible by the NCVO: $0 \equiv [A,A,A,B,B,B]$ and $1 \equiv [A,B,B,B]$. Then, we extract the bit sequences as represented in Fig. 2d–f above the graphs. In the commensurate states (Fig. 2d, f), the NCVO generates only one type of bits: 0 for $I_{dc}$ = 12.6 mA and 1 for $I_{dc}$ = 14.0 mA. On the other hand, in the incommensurate state (Fig. 2e), the NCVO generates bits in no apparent order. To better illustrate this, we plot a bit sequence generated in the chaotic regime for a longer duration (Fig. 3j). In the NCVO system, the bits are generated at an average rate of ~ 131 Mbit s⁻¹, which is much faster than random number generation using spintronic devices such as stochastic magnetic tunnel junctions³⁸,³⁹. The bits are produced by chaotic polarity switching of the core after a few revolutions around the nanocontact, so the generation rate is mainly proportional to the gyration frequency, $f_0$, which is in the range of 0.1–1 GHz in typical soft magnetic materials. Here, the gyration frequency can reach up to ~ 600 MHz, so the potential bit generation could reach up to ~ 300Mbit s⁻¹.

Interestingly, the possible cycles in the incommensurate state correspond to those the NCVO already exhibits in its
<figure><img src="image_7.png" /><figcaption>**Fig. 4** Complexity and control of bit sequences. a Probability of 0 and 1 in the generated bit sequences as a function of $I_{dc}$. b Probability of moving from the current state, *i*, to the next state, *j*, $P_{ij}$. The inset shows the Markov chain for a 1 bit information. c Shannon block entropy rate, *h*, as a function of $I_{dc}$ for different block lengths, 3 bits and 6 bits. d Normalized Lempel-Ziv complexity, $C_{LZ}$, as a function of $I_{dc}$.</figcaption></figure>

commensurate state at $I_{dc}$ = 12.6 and 14.0 mA, respectively (Fig. 3h). This shows that not only the complexity is driven by the core switching but also that the type of patterns generated in the incommensurate (chaotic) state is fundamentally restricted to accessible patterns associated to the two neighbouring commensurate states. Hence, within the incommensurate region, we anticipate that the probability of appearance of one among the two accessible patterns can be controlled by $I_{dc}$.

**Complexity of the generated bit sequences.** To assess the complexity of the bit sequence extracted from the symbolic dynamics, we compute the probability of each bit as a function of $I_{dc}$ (Fig. 4a). We also estimate the probability of transition from one bit to the next, $P_{ij}$, where *i* and *j* are the current and the next binary states (Fig. 4b), while assuming a Markovian property, i.e., a one-bit memory depth associated to the bit sequence (see Supplementary Fig. 5 and Supplementary Note 5). In Fig. 4a, as $I_{dc}$ increases and the NCVO transitions from a commensurate to an incommensurate state, we observe that the probability of symbol 0 (resp. symbol 1) decreases smoothly from $P_0$ = 1 (resp. increases from $P_1$ = 0) monotonically and reaches the point where $P_0 = P_1 = 0.5$ at about $I_{dc} \approx 13.1$ mA before continuing to decrease to $P_0 = 0$ (resp. to increase to $P_1 = 1$). Similarly, in Fig. 4b, we observe similar behaviour for the transition probabilities and specifically at $I_{dc} \approx 13.1$ mA, we have $P_{0,0} = P_{0,1} \approx 0.5$, $P_{1,0} \approx 0.75$, and $P_{1,1} \approx 0.25$. With this particular transition probability matrix, the entropy rate of the Markov chain describing the NCVO's pattern dynamics can be computed and reads $h_\infty \approx 0.93$ bit per binary symbol for a maximum entropy rate of 1 bit per binary symbol (see Supplementary Note 6). These results imply that for certain DC currents it becomes difficult to predict the next generated bit from the current one in the incommensurate state. But as the current is varied from this operating point, the entropy rate decreases because the NCVO dynamics approaches one of the two neighbouring commensurate states. As such, one of the two accessible

temporal patterns will start to dominate the other and its corresponding bit (0 or 1) will become more probable.

To further assess the unpredictability of the generated bits, we evaluate the complexity of the bit sequences generated by using two metrics from information theory: the Shannon block entropy rate, *h*, and the normalized Lempel-Ziv complexity, $C_{LZ}$. The Shannon block entropy rate measures the amount of uncertainty carried on average by each bit: for a binary source of information, as the one generated by a chaotic NCVO, we have maximum entropy $h_{max}$ = 1 bit per binary symbol. We compute the entropy rate as a function of $I_{dc}$ (Fig. 4c) and consider binary blocks of length *n* = 3 and *n* = 6 to obtain more robust estimates, while ensuring the estimation remains invariant with respect to these choices of block lengths (see Methods). We observe that the uncertainty from the bit stream generated by the NCVO is non-monotonic inside the incommensurate region; it gradually increases to reach a peak value of *h* ≈ 0.94 at $I_{dc}$ = 13.1 mA before decreasing again as the 1/3 commensurate state is approached. This is consistent with an asymmetric distribution for the probability mass function of the generated bit, which indicates that the bit sequence generated by the NCVO in its incommensurate state inherits complexity from the Kolmogorov-Sinai entropy created by the chaotic dynamics³⁵. The determination of the lower bound of KS entropy, $h_{KS} \approx 0.12$ ns⁻¹ at $I_{dc}$ = 13.2 mA, provides the maximum level of entropy extraction from the NCVO at ~100Mbit s⁻¹: This is consistent with the average generation rate of bits achieved by the symbolic encoding used in our study. We perform a similar analysis for the Lempel-Ziv complexity $C_{LZ}$, which measures the diversity (i.e., lack of redundancy) of binary patterns encountered in a binary sequence. We observe a similar trend with finite values in the incommensurate region, reaching almost the maximum value of $C_{LZ,max}$ = 1, while $C_{LZ} \approx 0$ in the commensurate state. This means the bits generated by the chaotic NCVO cannot be efficiently compressed because of its maximal complexity. Since *h* and $C_{LZ}$ almost attain their maximum values of ~1, the raw generated bit sequences have suitable statistical features to be considered as a physical source of entropy for information processing. The probability and complexity assessments have been performed on more than 9300 bit strings for each $I_{dc}$, which are obtained by the pattern recognition method shown in Supplementary Fig. 6 and described in Supplementary Note 7.

## Discussion

Based on our analysis, the chaotic magnetization dynamics of the NCVO has desirable properties for information processing applications. First, the complexity of the magnetization dynamics is not encoded in the amplitude of the waveform but in the alternation of regular patterns. In that sense, it natively exhibits very similar properties to the chaos generated by specifically engineered systems proposed in ref. ⁴⁰ and hence has high resilience to perturbations. The extraction of the bits, once the two patterns are identified and stored digitally, could be done in real time with a low computational cost using matched filters⁴¹ implemented on field programmable-gate arrays or digital signal processors. A recent application of chaotic dynamics similar to those of NCVOs was proposed for WIFI communications (i.e., for challenging environments with multi-path interference, jamming, distortion) using electronic circuits⁴²,⁴³, because of the capacity to recover more easily regular patterns in low signal-to-noise conditions compared with chaotic waveforms in amplitude and/or phase of the signal. Second, the NCVO generates directly chaotic bits with high level of entropy (e.g., entropy rate of $h_\infty \approx 0.93$ bit per binary symbol at $I_{dc}$ = 13.1 mA in the incommensurate state), hence making it a promising technology
for physical random number generation with minimal post-processing. Furthermore, with chaotic effects taking place at the nanoscale, dense integration of NCVOs on a single spintronic chip would allow for the parallel generation of random bits with aggregate rate in the tens or hundreds of Gbit s⁻¹ similar to the approach already used in ref. 44 with microelectronic circuits. Finally, the properties of the chaotic NCVO dynamics can be easily tailored with the DC current injected in the nanocontact. As a result, one could design advanced control strategies by injected small perturbation in the current to tune the chaos and encode data in pattern alternation, as described in refs. 11,45, and to create a robust chaos-based encryption at the physical layer.

In addition, knowing the underlying structure of the temporal chaotic dynamics and its connection to the timing of the core-switching events, one could also design experimental strategies to control electrically the core dynamics to encode information secretly. This paves the way for nanoscale chaos-based information processing using the nonlinear dynamics of spintronic devices.

## Methods

**Sample fabrication and measurements.** The nanocontacts are fabricated using the atomic force microscope nano-indentation method⁴⁶ on the top of the sputtered deposited multilayer with the composition SiO₂/Cu (40 nm)/Co (20 nm)/Cu (10 nm)/Ni₈₁Fe₁₉ (20 nm)/Au (6 nm)/photoresist (50 nm)/Au (nanocontact)²⁷. The diameter of the contact is ~ 20 nm. The vortex is first nucleated by reversing the free layer magnetization with an in-plane applied magnetic field in the presence of a $I_{dc}$ = 16 mA current applied through the nanocontact. The vortex gyration around the nanocontact results in magnetoresistance oscillations that are detected after amplification as voltage fluctuations in the frequency domain by a spectrum analyser and in the time domain by a single-shot oscilloscope. RF switches are used to connect either of these two equipments to the sample, hence allowing for both time- and frequency-domain measurements to be made sequentially. The experiments are conducted in a cryostat at liquid nitrogen temperature (77 K) to better isolate the chaotic dynamics, which is a deterministic process but can appear as a thermal noise, from thermal fluctuations which are true stochastic processes. Note that all the data in this study are measured using a single device. Further details of the experimental setup and measurement procedure are described elsewhere²⁶,²⁷.

**Thermal noise filtering from time traces.** To improve the signal-to-noise ratio of the experimental time traces, we used an averaging filter. We collected similar short-term waveforms (~ 7.5 ns) from full time series by calculating convolutions, then averaged over them. This method is applicable in our system because the output time traces are composed of only two or three patterns even in the chaotic regime. More details can be found in Supplementary Fig. 2 and Supplementary Note 2.

**Time-delay embedding and phase-space reconstruction.** A time-delay embedding procedure is used to form an *m*-dimensional vector space related to the original phase space by a diffeomorphism preserving topological invariants of the original attractor, if $m > 2d_A$ with $d_A$ dimension of the original attractor⁴⁷. The vectors in the reconstructed phase space are obtained from univariate time-resolved series as follows:

$$ \mathbf{v}_n^{(m)} = [V(t_n), V(t_n - \tau), \dots, V(t_n - (m-1)\tau)], \quad (1) $$

with the measured voltage, $V(t_n)$, sampled at discrete $t_n = n\Delta t$ with $\Delta t = 12.5$ ps the experimental sampling period. In this study, we choose the embedding dimension $m = 3$ and the time-delay embedding $\tau \approx 1/(4f_0)$ which is close to the first zero of the autocorrelation function of $V(t)$³⁵.

**Poincaré section and symbolic analysis.** To simplify the definition of the Poincaré section ($x_3 = 0$) in the symbolic analysis, we apply a unitary transformation to the lag coordinates (which does not affect the topological equivalence between the reconstructed and original phase spaces) and form $\mathbf{x}_n^{(m)} = \mathbb{U}\mathbf{v}_n^{(m)}$, where $\mathbb{U}$ is a rotation matrix,

$$ \mathbb{U} = \begin{bmatrix} 1 & 0 & 0 \\ 0 & \cos \theta_1 & -\sin \theta_1 \\ 0 & \sin \theta_1 & \cos \theta_1 \end{bmatrix} \begin{bmatrix} \cos \theta_3 & -\sin \theta_3 & 0 \\ \sin \theta_3 & \cos \theta_3 & 0 \\ 0 & 0 & 1 \end{bmatrix}, \quad (2) $$

with $\theta_1 = -20^\circ$ and $\theta_3 = 67^\circ$. The partition on the Poincaré section is set as $x_2 = 0.76x_1 - 0.005$ (Fig. 3a-c) delimiting both regions, $R_A$ and $R_B$. This choice, despite being arbitrary, allows us to capture both the incommensurate and commensurate regimes of the NCVO ($I_{dc}$ = 12.6 - 14.0 mA). A reconstructed orbit is encoded with symbolic sequences, whose length is determined by the number of times the reconstructed attractor intersects the two regions.

**Computation of correlation dimension and entropy.** The correlation dimension $D_c$ provides insight on the fractal dimension of an attractor and hence is used for the detection of chaos from the filtered time series. Its computation relies on the Grassberger–Procacia (GP) algorithm⁴⁸ involving the correlation sum $C_m(\epsilon)$, which gives the average number of neighbouring vectors $\mathbf{x}_j^{(m)}$ within the range $\epsilon > 0$ from any given vectors $\mathbf{x}_i^{(m)}$ of the attractor obtained from the time-delay embedding procedure. It is defined as

$$ C_m(\epsilon) = \frac{2}{(N - n_T)(N - 1 - n_T)} \sum_{j=1}^{N} \sum_{i=1}^{N} \Theta(\epsilon - \| \mathbf{x}_i^{(m)} - \mathbf{x}_j^{(m)} \|), \quad (3) $$

with $\epsilon$ representing the typical radius of the neighbourhood surrounding the vector $\mathbf{x}_i^{(m)}$, $\| \cdot \|$ the norm-2, and $\Theta$ the Heaviside function. To avoid the bias induced by finite-size effects of the time series and time-correlation of neighbouring vectors, we introduce the Theiler correction $|i-j| > n_T$ to select eligible neighbours for the computation⁴⁹. The presence of self-similarity imposes the correlation sum to approximately satisfy a linear growth in log-scale and hence as a constant value (plateau) for its derivative. Hence, the correlation dimension $D_c$ is given by

$$ D_c = \lim_{\epsilon \to 0} \lim_{N \to \infty} \frac{\partial \ln C_m(\epsilon)}{\partial \ln \epsilon}. \quad (4) $$

Using the correlation sum $C_m(\epsilon)$, it is also possible to compute a lower bound for the Kolmogorov–Sinai (KS) entropy $h_{KS}$ known as the $K_2$-entropy⁵⁰,⁵¹ and defined as

$$ K_2 = \lim_{m \to \infty} \lim_{\epsilon \to 0} \lim_{N \to \infty} K_{2,m}(\epsilon) \text{ with } K_{2,m}(\epsilon) = \frac{1}{\tau} \ln \left[ \frac{C_m(\epsilon)}{C_{m+1}(\epsilon)} \right], \quad (5) $$

with $\tau = n_T \Delta t$ the embedding time lag derived from the sampling time of the experimental time series³⁵. The computation of $K_2$ can be also achieved by taking the limit of $K_{2,m}(\epsilon) = 1/4\tau \ln[C_m(\epsilon)/C_{m+4}(\epsilon)]$, which is the quantity introduced in ref. ⁵⁰ and used here. In practice, the approximate value of $K_2$ is determined by extrapolating the average value $(K_{2,m}(\epsilon))_\epsilon$ on the scaling range (here $[1, 10^{0.3}]$) as $m$ takes larger values. In our analysis, we have used this approach after normalizing the experimental filtered time series as described in the previous section. We use $N = 1.2 \times 10^5$ samples and $n_T = 15$ for the Theiler correction, with the embedding dimension $m$ chosen between 10 and 24 for the correlation dimension and entropy calculations. The neighbourhood radii $\epsilon$ in the range $[10^{-1}, 10]$.

**Computation of the Shannon block entropy.** We consider a random source of *n*-bit words from the dictionary $\{s_i\}_{n_i}$ with $1 \le n_i \le 2^n$. The words are obtained by sliding a window of *n* bits in width along the bit stream of length $N_b$ resulting from the symbolic analysis of the NCVO dynamics. The mathematical definition of the Shannon block entropy of the *n*-bit word source is given by

$$ H_n = - \sum_{s_i \in \{s_1, \dots, s_n\}} p(s_i) \log_2 p(s_i), \quad (6) $$

with $p(s_i)$ the probability of appearance of symbol $s_i$. We can then determine the entropy per binary symbol (or entropy rate with a maximum value at 1 bit per symbol) using the limit $h = \lim_{n \to \infty} H_n/n$. In finite binary sequences of length $N_b$, the use of blocks of length *n* usually leads to more robust estimates of the entropy compared with the direct estimation from the bit sequence. The probability of each word is determined with the likelihood estimator $\hat{p}(s_i) \approx \#(s_i)/N_b$ with $N_b = N_b - n + 1$. Finally, we use the upper limit for the block size given by

$$ N_b h \ge n 2^n \ln 2 \quad (7) $$

for binary words as suggested in ref. ⁵². Due to the finite size $N_b = 9300$ bits for the bit stream obtained from the NCVO's symbolic dynamics, we use block lengths in the range $n \in \{3, 6\}$.

**Computation of the Lempel–Ziv complexity.** The Lempel–Ziv complexity $C_{LZ}$ of a binary sequence measures the number of patterns present and is the basis of LZ77 compression⁵³,⁵⁴. For a large binary sequence of size *n*, it can be shown that the number of patterns $c(n)$ behaves asymptotically like the ratio $n/\log_2 n$. In a sequence of length *n* ≫ 1, we can use this ratio as a normalization factor for the number of patterns in order to ensure that the complexity measure remains bounded in the range $C_{LZ} \in [0, 1]$. In order to compute $C_{LZ}$, we use the algorithmic procedure presented in ref. ⁵⁵. If the binary sequence is generated by a stationary and ergodic process, the Lempel–Ziv complexity coincides with the entropy rate $h^{52}$ in the limit of large *n*. Here, the complexity $C_{LZ}$ is computed from a sequence of $N_b = 9300$ bits.

**Reporting summary.** Further information on research design is available in the Nature Research Reporting Summary linked to this article.
**Data availability**

The data sets generated and/or analysed during the current study are available from the
corresponding author on reasonable request.

**Code availability**

All relevant codes are available from the corresponding author on reasonable request.

Received: 12 August 2019; Accepted: 21 December 2019;

Published online: 30 January 2020

**References**

1. Strogatz, S. H. *Nonlinear Dynamics and Chaos: With Applications to Physics, Biology, Chemistry, and Engineering.* 2nd edn (Westview Press, Philadelphia, 2015).

2. Sciamanna, M. & Shore, K. A. Physics and applications of laser diode chaos. *Nat. Photonics* **9**, 151–162 (2015).

3. Argyris, A. et al. Chaos-based communications at high bit rates using commercial fibre-optic links. *Nature* **438**, 343–346 (2005).

4. Uchida, A. et al. Fast physical random bit generation with chaotic semiconductor lasers. *Nat. Photonics* **2**, 728–732 (2008).

5. Li, W. et al. Fast physical random-number generation based on room-temperature chaotic oscillations in weakly coupled superlattices. *Phys. Rev. Lett.* **111**, 044102–5 (2013).

6. Rontani, D., Choi, D., Chang, C.-Y., Locquet, A. & Citrin, D. S. Compressive sensing with optical chaos. *Sci. Rep.* **6**, 35206 (2016).

7. Ditto, W. L. & Sinha, S. Exploiting chaos for applications. *Chaos* **25**, 097615 (2015).

8. Keuninckx, L. et al. Encryption key distribution via chaos synchronization. *Sci. Rep.* **7**, 43428 (2017).

9. Lin, F.-Y. & Liu, J.-M. Chaotic Lidar. *IEEE J. Sel. Top. Quant. Electron.* **10**, 991–997 (2004).

10. Myneni, K., Barr, T. A., Reed, B. R., Pethel, S. D. & Corron, N. J. High-precision ranging using a chaotic laser pulse train. *Appl. Phys. Lett.* **78**, 1496–1498 (2001).

11. Hayes, S., Grebogi, C. & Ott, E. Communicating with chaos. *Phys. Rev. Lett.* **70**, 3031–3034 (1993).

12. Wigen, P. E. (ed). *Nonlinear Phenomena and Chaos in Magnetic Materials* (World Scientific, Singapore, 1994).

13. Alvarez, L. F., Pla, O. & Chubykalo, O. Quasiperiodicity, bistability, and chaos in the Landau-Lifshitz equation. *Phys. Rev. B* **61**, 11613–11617 (2000).

14. Bertotti, G., Serpico, C. & Mayergoyz, I. D. Nonlinear magnetization dynamics under circularly polarized field. *Phys. Rev. Lett.* **86**, 724–727 (2001).

15. Lee, K.-J., Deac, A., Redon, O., Nozières, J.-P. & Dieny, B. Excitations of incoherent spin-waves due to spin-transfer torque. *Nat. Mater.* **3**, 877–881 (2004).

16. Yang, Z., Zhang, S. & Li, Y. C. Chaotic dynamics of spin-valve oscillators. *Phys. Rev. Lett.* **99**, 134101–4 (2007).

17. Slavin, A. & Tiberkevich, V. Nonlinear auto-oscillator theory of microwave generation by spin-polarized current. *IEEE Trans. Magn.* **45**, 1875–1918 (2009).

18. Pylypovskyi, O. V., Sheka, D. D., Kravchuk, V. P., Mertens, F. G. & Gaididei, Y. Regular and chaotic vortex core reversal by a resonant perpendicular magnetic field. *Phys. Rev. B* **88**, 014432 (2013).

19. Bondarenko, A. V., Holmgren, E., Li, Z. W., Ivanov, B. A. & Korenivski, V. Chaotic dynamics in spin-vortex pairs. *Phys. Rev. B* **99**, 054402 (2019).

20. Montoya, E. A. et al. Magnetization reversal driven by low dimensional chaos in a nanoscale ferromagnet. *Nat. Commun.* **10**, 543 (2019).

21. Ralph, D. C. & Stiles, M. D. Spin transfer torques. *J. Magn. Magn. Mater.* **320**, 1190–1216 (2008).

22. Hoffmann, A. Spin Hall effects in metals. *IEEE Trans. Magn.* **49**, 5172–5193 (2013).

23. Kim, J.-V. in *Solid State Physics* (eds. Camley, R. E. & Stamps, R. L.) Vol. 63, Cap. 4, 217–294 (Academic Press, San Diego, 2012).

24. Chen, T. et al. Spin-torque and spin-Hall nano-oscillators. *Proc. IEEE* **104**, 1919–1945 (2016).

25. Villard, P. et al. A GHz spintronic-based RF oscillator. *IEEE J. Solid-State Circuits* **45**, 214–223 (2010).

26. Petit-Watelot, S. et al. Commensurability and chaos in magnetic vortex oscillations. *Nat. Phys.* **8**, 682–687 (2012).

27. Devolder, T. et al. Experimental observation of chaos in nanocontact vortex oscillators. *Phys. Rev. Lett.* **123**, 147701 (2019).

28. Pribiag, V. S. et al. Magnetic vortex oscillator driven by d.c. spin-polarized current. *Nat. Phys.* **3**, 498–503 (2007).

29. Dussaux, A. et al. Large microwave generation from current-driven magnetic vortex oscillators in magnetic tunnel junctions. *Nat. Commun.* **1**, 8 (2010).

30. Locatelli, N. et al. Dynamics of two coupled vortices in a spin valve nanopillar excited by spin transfer torque. *Appl. Phys. Lett.* **98**, 062501 (2011).

31. Mistral, Q. et al. Current-driven vortex oscillations in metallic nanocontacts. *Phys. Rev. Lett.* **100**, 257201 (2008).

32. Zhang, S. & Li, Z. Roles of nonequilibrium conduction electrons on the magnetization dynamics of ferromagnets. *Phys. Rev. Lett.* **93**, 127204 (2004).

33. Pufall, M., Rippard, W., Schneider, M. & Russek, S. Low-field current-hysteretic oscillations in spin-transfer nanocontacts. *Phys. Rev. B* **75**, 140404 (2007).

34. Keatley, P. S. et al. Direct observation of magnetization dynamics generated by nanocontact spin-torque vortex oscillators. *Phys. Rev. B* **94**, 060402 (2016).

35. Kantz, H. & Schreiber, T. *Nonlinear Time Series Analysis* (Cambridge University Press, 2004).

36. Poon, C.-S. & Barahona, M. Titration of chaos with added noise. *Proc. Natl. Acad. Sci. USA* **98**, 7107 (2001).

37. Bollt, E. M. Review of chaos communication by feedback control of symbolic dynamics. *Int. J. Bifurcat. Chaos* **13**, 269–285 (2003).

38. Fukushima, A. et al. Spin dice: a scalable truly random number generator based on spintronics. *Appl. Phys. Express* **7**, 083001-5 (2014).

39. Vodenicarevic, D. et al. Low-energy truly random number generation with superparamagnetic tunnel junctions for unconventional computing. *Phys. Rev. Appl.* **8**, 9–9 (2017).

40. Corron, N. J., Blakely, J. N. & Stahl, M. T. A matched filter for chaos. *Chaos* **20**, 023123 (2010).

41. Proakis, J. G. & Salehi, M. *Digital Communications*, 5th edn (McGraw-Hill, 2007).

42. Ren, H., Baptista, M. & Grebogi, C. Wireless communication with chaos. *Phys. Rev. Lett.* **110**, 184101 (2013).

43. Ren, H., Bai, C., Liu, J., Baptista, M. & Grebogi, C. Experimental validation of wireless communication with chaos. *Chaos* **26**, 083117 (2016).

44. Rosin, D., Rontani, D. & Gauthier, D. J. Ultrafast physical generation of random numbers using hybrid Boolean networks. *Phys. Rev. E* **87**, 040902 (2013).

45. Hayes, S., Grebogi, C. & Ott, E. Experimental control of chaos for communication. *Phys. Rev. Lett.* **73**, 1781–1784 (1994).

46. Bouzhououane, K. et al. Nanolithography based on real-time electrically controlled indentation with an atomic force microscope for nanocontact elaboration. *Nano Lett.* **3**, 1599–1602 (2003).

47. Takens, F. in *Dynamical Systems and Turbulence*, Vol. 898 of *Lecture Notes in Mathematics* (eds. Rand, D. A. & Young, L.-S.) 366–381 (Springer-Verlag, Heidelberg, 1981).

48. Grassberger, P. & Procacia, I. Characterization of strange attractor. *Phys. Rev. Lett.* **50**, 346 (1983).

49. Theiler, J. Spurious dimensions from correlation algorithms applied to limited time-series data. *Phys. Rev. A* **34**, 2427–2432 (1986).

50. Grassberger, P. & Procacia, I. Estimation of the Kolmogorov entropy from a chaotic signal. *Phys. Rev. A* **28**, 2591 (1983).

51. Eckmann, J.-P. & Ruelle, D. Ergodic theory of chaos and strange attractors. *Rev. Mod. Phys.* **57**, 617–656 (1985).

52. Lesne, A., Blanc, J.-L. & Pezard, L. Entropy estimation of very short symbolic sequences. *Phys. Rev. E* **79**, 046208 (2009).

53. Lempel, A. & Ziv, J. On the complexity of finite sequences. *IEEE Trans. Inform. Theory IT-22*, 75–81 (1976).

54. Lempel, A. & Ziv, J. A universal algorithm for sequential data compression. *IEEE Trans. Inform. Theory IT-23*, 337–343 (1977).

55. Kaspar, F. & Schuster, H. G. Easily calculable measure for the complexity of spatiotemporal patterns. *Phys. Rev. A* **36**, 842–848 (1987).

**Acknowledgements**

The authors thank Cyrile Deranlot and Stéphanie Girod for their assistance in film growth and sample preparation. This work was supported by the Agence Nationale de la Recherche under Contract No. ANR-17-CE24- 0008 (CHIPMuNCS), the Horizon2020 Research Framework Programme of the European Commission under Contract No. 751344 (CHAOSPIN), and the French RENATECH network. The Chaire Photonique is funded by the European Regional Development Fund (European Commission) (ERDF), Ministry of Higher Education and Research (FNADT), Moselle Department, Grand Est Region, Metz Metropole, AIRBUS-GDI Simulation, CentraleSupélec, and Fondation CentraleSupélec.
**Author contributions**

J.-V.K., S.P.-W., V.C. and M.S. designed the study. K.B. and V.C. fabricated the samples. T.D. designed and implemented the experimental setup. M.-W.Y., J.L. and T.D. performed the high-frequency electrical measurements. M.-W.Y., D.R. and J.-V.K. analysed and interpreted the data. M.-W.Y. performed the simulations and interpreted the results. M.-W.Y. and D.R. prepared the manuscript. All authors edited and commented on the manuscript.

**Competing interests**

The authors declare no competing interests.

**Additional information**

**Supplementary information** is available for this paper at https://doi.org/10.1038/s41467-020-14328-7.

**Correspondence** and requests for materials should be addressed to M.-W.Y. or J.-V.K.

**Peer review information** *Nature Communications* thanks Andrei Slavin and the other, anonymous, reviewer(s) for their contribution to the peer review of this work. Peer reviewer reports are available.

**Reprints and permission information** is available at http://www.nature.com/reprints

**Publisher's note** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_8.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2020