Abstract
An algorithm to compute the Lyapunov exponents of piecewise linear function-based multidirectional multiscroll chaotic oscillators is reported. Based on the regions in the piecewise linear functions, the suggested algorithm determines the individual expansion rate of Lyapunov exponents from -piecewise linear variational equations and their associated -Jacobian matrices whose entries remain constant during all computation cycles. Additionally, by considering OpAmp-based chaotic oscillators, we study the impact of two analog design procedures on the magnitude of Lyapunov exponents. We focus on analyzing variations of both frequency bandwidth and voltage/current dynamic range of the chaotic signals at electronic system level. As a function of the design parameters, a renormalization factor is proposed to estimate correctly the Lyapunov spectrum. Numerical simulation results in a double-scroll type chaotic oscillator and complex chaotic oscillators generating multidirectional multiscroll chaotic attractors on phase space confirm the usefulness of the reported algorithm.
1. Introduction
The deterministic, although unpredictable, behavior of the phenomenon known as chaos, which appears in nonlinear dynamical systems under certain conditions, is a key area of research in several fields of science, ranging from mathematics and biology to engineering [1–3]. Nevertheless, as the chaotic behavior is highly sensitive to small variations on initial conditions and parameter values, quantitative and qualitative characterizations need to be used to confirm it. Mainly, a nonlinear dynamical system is characterized by using fractal dimension, bifurcation diagrams, Poincaré maps, flat power spectrum, and Lyapunov exponents [4]. Lyapunov exponent (LE) depicts the rate of separation and unification of close trajectories in a dynamical system and provides the most characteristic description of a deterministic nonperiodic flow involved [5]; that is, a nonlinear dynamical system must have one positive LE, at least, for being classified as chaotic.
Chaotic systems are categorized into maps (discrete-time) and flows (continuous-time). Whereas the chaotic maps can be modeled by using only one-step recursive equation, chaotic flows are represented by a coupled set of first-order differential equations. Nevertheless, continuous-time chaotic systems can exhibit multiscroll in multidirections. In particular, we are interested in continuous-time chaotic oscillators with their nonlinearities being modeled by a piecewise linear (PWL) function [6]. These systems can generate one-directional (1D) -scrolls, two-directional (2D) grid scroll, and three-directional (3D) grid scroll chaotic attractors by adding breakpoints in the PWL function and increasing the number of PWL functions into the nonlinear system [7, 8]. Furthermore, the scroll number can be chosen arbitrarily, and the scrolls can be located anywhere in the phase space. This provides a systematic way to design multidirectional multiscroll chaotic oscillators as it was previously demonstrated in [9, 10]. Thereby, continuous-time PWL chaotic systems are regularly used as the core of engineering applications such as in liquid mixing [11], autonomous mobile robots [12], and synchronization methods with emphasis on encrypted communications [13, 14].
However, although most of them have been experimentally verified by using analog circuits constructed with both discrete electronics components [15–17] and integrated circuit technology [18, 19], a systematic approach does not exist to evaluate their chaotic regime considering design tradeoffs at electronic system level. For instance, as a function of the design requirements, an electronic circuit designer applies regularly two analog design techniques named as frequency-scaling (FS) and excursion level-scaling (ELS). As it has been reported in [9, 11–13, 15–18, 20] these techniques are applied to modify the frequency bandwidth and voltage/current dynamic range of the chaotic signals, respectively. Unfortunately, those practical modifications cause that the magnitude of LEs to increase or decrease, for example, a negative LE could be a positive value resulting in a false chaotic behavior. It would lead to an erroneous interpretation of the dynamics.
Recently, automatic analog circuit tools and evolutionary computation approaches have been introduced to either obtain new circuit topologies for continuous-time PWL chaotic oscillators or increase their complexity by optimizing the maximal Lyapunov exponent [9, 21, 22]. In both cases, they must evaluate a large design space to verify if potential solutions are in the chaotic regime. No matter what kind of design methodology is followed, that is, custom-made designs or automated approaches, an algorithm capable of determining the chaotic behavior of PWL function-based chaotic oscillators at electronic system level is vital to design further chaos-based applications.
To the best of our knowledge, Lyapunov spectrum, which means the collection of all LEs, has been predominantly researched on chaotic systems based on a continuous nonlinear function and by considering analytical relations among state-variables [5, 23–29]. Opposite of that are the algorithms focusing on computing Lyapunov spectrum of PWL function-based chaotic oscillators.
In this paper, we report an algorithm to compute the Lyapunov spectrum of continuous-time double-scroll type chaotic oscillators and complex chaotic oscillators generating multidirectional multiscroll chaotic attractors. The novelty of the proposed method is to determine the Lyapunov spectrum from the individual contributions of the -piecewise linear variational equations and their associated -Jacobian matrices associated with the regions of a given PWL function. It leads to having a constant entry for those matrices during all computation cycles, which reduces the complexity of the algorithm. Additionally, in a previous research [9, 20, 30], a scheme for generating PWL function-based multiscroll multidirectional chaotic attractors with an analog circuit designed by using operational amplifiers (OpAmp’s) was presented. Based on such work, we propose a renormalization factor to compute correctly the magnitude of Lyapunov vectors under variations of the frequency bandwidth and the voltage/current dynamic range of chaotic signals at electronic system level. Therefore, the suggested approximation introduces a systematic procedure that can be successfully applied to several well known continuous-time PWL function-based chaotic oscillators.
The paper is organized as follows. Section 2 contains the basic concepts of Lyapunov exponents. Section 3 introduces the proposed algorithm to determine the LEs from -PWL variational equations as well as its accuracy analysis. Besides, a description of its application is introduced when frequency bandwidth and dynamic range variations at electronic systems level are considered. Section 4 shows up the numerical simulations results for Chua’s circuit, 1D 6-scroll chaotic oscillator, and 2D grid scroll chaotic oscillator. Section 5 discusses the electronic design tradeoffs. Finally, the main conclusions are summed up in Section 6.
2. Computation of Lyapunov Exponents
As well known, the LEs characterize the expansion and contraction rate of the solutions of a dynamical system as a function of long-term sensitivity with respect to its initial conditions. Let us consider a continuous-time dynamical system depicted by whereandare vector fields of-order. To determine the LEs of (1) its variational equation is considered, which is defined by being the Jacobian matrix and the variational state-variable [25]. A solution of (2) under initial conditions can be written as , with matrix satisfying
Hereindicates the identity matrix. If we relate this solution to the evolution of the axes of an infinitesimal ellipsoid by means of for , where denotes an orthogonal basis; then the Lyapunov vectors are defined as the expansion rate of the axes length . Therefore, by solving (3) the LEs are computed with
3. Proposed Method Based on -PWL Variational Equations
A PWL functions-based autonomous chaotic system can be modeled by applying state-variables approach [13], which consists of a set of first-order differential equations describing the evolution of the system-variables as follows: where is the state matrix containing the linear part of such systems, is the input matrix including the PWL functions, is an input unit vector, and is the state-variables vector. Note that a PWL function can be represented by a set of affine equations. In each one of the regions, an affine equation is given by (6), which includes a slope and an offset . Because 1D-, 2D- and 3D-multiscroll chaotic oscillators depend on several PWL functions, we divide their contributions per region into two matrices as follows:
Here, matrix contains all slopes and matrix collects all offsets and and are the entries for those matrices, respectively, whose position is denoted by indexesandaccording to the effect caused in the state-variable system [6]. This means that index specifies which state-variable is being affected by PWL function while index indicates the state-variable that controls the PWL function. Finally, refers to the number of regions in PWL function. As a result, we transform the problem of solving a nonlinear system into a sequence of linear and purely algebraic equations, which can be computed with a first-order numerical solver as the Forward-Euler algorithm used herein. After that, by combining (5), (6), and (7) as follows: we obtain -PWL systems, where and are matrices of appropriate dimensions. In this framework, (8) gives a solution of (3) in the form of yielding-PWL variational systems and their associated -Jacobian matrices defined by . It is important to remark that we obtain -Jacobian matrices whose elements are constant for all computation cycles. Besides, based on the location of chaotic trajectory through the regions of PWL function, the algorithm integrates only one PWL system and one PWL variational system.
In this manner, to compute the LEs, our method proceeds to integrate (9) with initial condition to obtain the evolution of the-PWL variational systems . Here, is an orthogonal basis and are the Lyapunov vectors at time . Later, the Gram-Schmidt orthogonalization is applied to matrix to get . From the norm of the vectors , the LEs are finally accumulated with
3.1. Renormalization of the Lyapunov Spectrum
As it has been reported, FS and ELS design techniques modifie the frequency bandwidth and voltage/current dynamic range respectively [9, 20]. Those modifications produce the magnitude of LEs that grows according to the scaling applied, leading to inaccurate values. We address this issue by renormalizing of each certain number of orbital periods the Lyapunov spectrum. An orbital period refers to the fundamental frequency bandwidth of the chaotic oscillators, and that is determined by the passive and active electronic devices used to implement the chaotic oscillator physically. In this work, we renormalize each two orbital periods set with . Since eigenvalues are strongly related to frequency scaling procedures [20], the renormalization factor (RNF) is obtained from , with being the absolute minimum of all the eigenvalues of -Jacobian matrices in (9). Hence, (10) is rewritten as
The orthogonality among Lyapunov vectors is kept because the RNF affects the magnitude only but not the angle. To check the accuracy of (11), we consider the error in orthogonality [6]. This is calculated by where and are the matrices generated by the Gram-Schmidt method in and its transpose, respectively. When the absolute value of that difference tends towards to zero, the approximation of Lyapunov spectrum will be better. To guarantee small errors using (8), (9), and (11), the orthogonalization must be performed in short orbital periods of the dynamical solution of (1) [24]. Therefore, it highlights a straightforward approach to determine all LEs of PWL chaotic oscillators taking into account design considerations at electronic circuit level. As a consequence the proposed method, summarized in the flow diagram of Figure 1, can be implemented into the workflow of automated circuit synthesis tools [9].
Regarding the implementation, this algorithm has been developed in Maple 14.0. The numerical kernel, including the subroutines for the Forward-Euler and Gram-Schmidt methods (developed by the authors), has been coded in the native Maple programming language. Note that, from the flow diagram of Figure 1, other specific and general programming languages, for example, MATLAB, C, and Fortran, can be used to implement the algorithm.
4. Numerical Simulations
In this section, it is demonstrated that the proposed algorithm can be successfully applied to compute the LEs of double-scroll and complex chaotic oscillators generating -scrolls with multiple directions on phase space, at electronic system level. Numerical simulation results for Chua’s Circuit, 1D 6-scroll chaotic oscillator, and 2D 4-scroll chaotic oscillator are shown. Additionally, LEs are computed when the frequency bandwidth and voltage/current excursion of these chaotic oscillators are modified at electronic system level.
4.1. Chua’s Circuit
Chua’s circuit is the first autonomous electronic circuit where a chaotic waveform was observed experimentally, established numerically, and proved theoretically [18].
Chua’s circuit consists of five circuit elements: one linear resistor , one inductor , two capacitors and , and one nonlinear resistor as shown in Figure 2(a). This circuit produces chaotic oscillations as a function of its parameters. The Chua’s circuit is modeled using the state-variable approach as previously defined by (5) to apply our algorithm. Consequently, we obtain where the nonlinear resistor is approximated by a current-voltage (-) PWL function defined by
(a)
(b)
In Figure 2(b) it is observed that, for a voltage signal of the state-variable less than the hyperplane (HP1), the PWL function becomes a linear region with a negative slope . For absolute voltages between HP1 and HP2, it has two linear regions of negative slope .
It is possible to represent Chua’s circuit by means of the modified state-variable equation defined in (8) by incorporating these slopes and offsets in equation (7). Accordingly, we divide the original nonlinear system in 3-PWL systems and 3-PWL variational equations, which are linked with the three regions shown in Figure 2(b) as follows:
In each region, the algorithm proceeds to solve (15a), (15b), and (15c) by using the Forward-Euler method. By applying (11), we obtain the Lyapunov spectrum shown in Figure 3(a). Also, the orthogonalization error defined by (14) is plotted in Figure 3(b).
(a)
(b)
4.2. 1D Multiscroll Chaotic Oscillator
The algorithm can also be used to compute the Lyapunov spectrum of complex chaotic systems as shown herein. To increase the complexity of chaotic behavior, multidirectional multiscroll chaotic systems have been proposed [9, 13, 15]. For instance, a multiscroll chaotic oscillator generates more than two scrolls on phase space in contrast to Chua’s circuit. This is produced when the PWL function has multiple regions and, therefore, multiple equilibrium points. Besides, multidirectional means that -scrolls are propagated in one, two, or three directions, that is, a grid of scrolls. For this behavior to occurs, it is necessary to add more PWL functions to the chaotic system. In this context, we study the 1D- and 2D-multiscroll chaotic oscillators shown in Figure 4 [6, 9, 13]. Similar to Chua’s oscillator, we found the state-variables equation of those oscillators by applying Kirchhoff circuit laws as for Figures 4(a) and 4(b) respectively, in both cases , , are the state-variables; , , , , and are the system parameters. Besides, , are PWL functions described by with
(a)
(b)
Here, defines the slope of the PWL functions; delimits the hyperplane interfaces; sets the saturation levels in odd multiples; establishes the distance between the middle points of adjacent slopes in even multiples; and is the saturation voltage for the OpAmp’s. All those parameters are associated with the finite gain model of OpAmp’s [13, 20]. Besides, and are positive integers related to the number of scrolls. This means that (16) and (17) can generate 1D scroll and 2D scroll grid, respectively [6]. In particular, we analyze the Lyapunov spectrum of 1D 6-scroll chaotic oscillator generated by selecting . For this case, there are eleven regions as sketched in Figure 5(a). In this manner, we split (16) intoas a function of with respect to the hyperplanes given by respectively. Since just slopes contribute to obtain the modified state-variables equation defined in (8), we get five different PWL variational equations denoted as , , , , and in (20b), (20d), (20f), (20h), and (20j). On the other hand, represents a unique variational equation for the six saturated levels in Figure 5(a). By using (11), the algorithm integrates the -PWL systems and-PWL variational equations in (20a), (20b), (20c), (20d), (20e), (20f), (20g), (20h), (20i), (20j), and (20k) for computing the LEs and its error as shown in Figure 6.
(a)
(b)
(a)
(b)
4.3. 2D Multiscroll Chaotic Oscillator
The proposed algorithm also determines the LEs of multidirectional chaotic systems. Hence, we prove this point by computing the Lyapunov spectrum of 2D 4-scroll chaotic oscillator shown in Figure 4(b). As a function of the two PWL functions of Figure 5(b), we can derive 49-PWL systems and 49-PWL variational systems. For the sake of simplicity, those equations are not shown. Besides, they can be obtained analogous to 1D 6-scroll chaotic oscillator of Section 4.2 Numerical simulation results for the 2D case are sketched in Figure 7.
(a)
(b)
5. Discussion
In this section, the numerical aspects of our algorithm and how these are associated with the electronic design requirements are discussed. First, by designing Chua’s circuit of Figure 2 with the elements listed in Table 1, a frequency bandwidth of about 411e3 rad/s is obtained. Taking this frequency into account and considering a long-term evolution, the Lyapunov spectrum simulating 46 oscillations of the chaotic signals is computed and shown in Figure 3(a).
The LEs for Chua’s circuit are established nearly to 0.25 ms. Also, the orthogonalization can be performed every two orbital periods of the fundamental frequency bandwidthwhile the error remains limited to minimum values as shown in Figure 3(b). This reduces the number of operations per simulation cycle.
Otherwise, it is important to note that, if the Lyapunov spectrum is not renormalized, the maximum and zero Lyapunov exponents converge to inaccurate magnitudes and it will conduct an imprecise estimation of dynamics. This is more evident in the maximum and zero Lyapunov exponents, for example, without renormalization of Chua’s circuit they are and , respectively, where these variations are induced by the time-constants of the electronic circuits, causing rapid evolution of the chaotic trajectories. This issue is addressed by using a renormalization factor (RNF), which is dependent on the absolute minimum eigenvalue of Jacobian matrices. Evidently, this attenuation factor changes as a function of the frequency scaling done by a circuit designer. Furthermore, as mentioned in Section 2, the error of the algorithm is intrinsically associated with the orbital period used to orthogonalize the Lyapunov vectors by using Gram-Schmidt method. In Figures 3(b), 6(b), and 7(b), it is noted that the error also varies as the complexity of the PWL chaotic systems grows. In this case, complexity refers to the increasing number of equilibrium points in the attractor. For instance, the error for 2D 4-scroll chaotic oscillator with 47 equilibrium points is greater than Chua’s circuit with only three equilibrium points, as shown in Table 1. This behavior agrees with previous papers [5, 23–29].
Moreover, not only the frequency bandwidth of PWL chaotic systems can be scaled, that is, by reducing the hyperplanes in Figures 2(b) and 5, but also voltage excursion of chaotic signals can be modified. In Table 1 is noted that frequency scaling (FS) procedure is the most dominant when approximating LEs because it is strongly linked to the evolution speed of solution of the chaotic trajectories. To prove this fact, four different electronic designs are analyzed. First, the normalized case of the 1D 6-scroll chaotic oscillator is considered. By using the elements listed in Table 1, Normalized design, the frequency bandwidth was set to 1 rad/s. Therefore, the RNF used to compute the Lyapunov spectrum shown in Figure 6 is around 0.6862. As this value approximates to unity, it means that the magnitude of LEs is less affected by the changes performed at the circuit level. Then, by designing chaotic oscillators at low frequency, inaccurate values for LEs are avoided. However, as reported in [6], this electronic design is impractical since its frequency bandwidth is narrow and it needs power supplies of about ±45 V to handle the voltage dynamic range of chaotic signals.
In this manner, the 1D 6-scroll chaotic oscillator in Figure 4(a) is redesigned changing and with the values in Table 1 labeled as ELS design. As these passive elements are associated with the slope and offset of the PWL functions in Figure 5(a), the voltage excursion of the chaotic signals is modified with about ±250 mV while the frequency bandwidth remains normalized [9]. Hence, the dynamic range is useful in real applications but not the frequency bandwidth. By applying the proposed algorithm, the Lyapunov spectrum for this case is shown in Figure 8(a). Because only the slope of PWL functions is included in the modified state matrix in (8), a few elements of Jacobian matrices are altered. Since RNF depends on eigenvalues of those Jacobians, this suffers a slight modification but continues being close to unity.
(a)
(b)
Furthermore, the impact in the Lyapunov spectrum is studied when both the frequency bandwidth and voltage dynamic range of chaotic signals are varied. By keeping the voltage dynamic range of ELS design, the frequency bandwidth of 1D 6-scroll chaotic oscillator of Figure 4(a) is increased from 1 rad/s to 65e3 rad/s. This redesign, labeled as FS-ELS design, was obtained by modifying the capacitors with the value shown in Table 1. As the absolute minimum eigenvalue depends inversely on the time-constant defined by capacitors, the Lyapunov spectrum of Figure 8(b) was computed by using a RNF of 0.000004. This means that there exists a rapid growth of the dynamical solution of 1D 6-scroll chaotic oscillator. Associated with this, the chaotic system in (16) is multiplied by an equivalent constant to the magnitude of RNF affecting all elements of the modified state-matrices and its associated Jacobians in (20a), (20b), (20c), (20d), (20e), (20f), (20g), (20h), (20i), (20j), and (20k). Therefore, RNF increases as a function of the frequency bandwidth. It is important to remark that, at this oscillation frequency, the chaotic circuit is practical to be used in engineering applications, for example, analog communications for encrypted voice transmission [13].
Finally, the study is extended to multidirectional chaotic systems. In particular, the Lyapunov spectrum of 2D 4-scroll chaotic oscillator propagating four scrolls in -axis and -axis directions of Figure 4(b) is computed. By selecting the elements in Table 1, labeled as 2D-FS-ELS design, the frequency bandwidth is set to 105e3 rad/s and voltage dynamic range of ±250 mV. Analogous to 1D chaotic oscillator, frequency scaling was the dominant procedure. Accordingly, the LEs shown in Figure 7(a) were approximated using a RNF of 0.000002.
From the results summarized in Table 1 the proposed algorithm is vital to evaluate the chaotic behavior of PWL multidirectional multiscroll chaotic oscillators at electronic system level. The positive LE (LE1 in Table 1) of Chua’s circuit computed with this algorithm is quite similar in ~0.3% of relative error to that reported in previous paper [29]. This verifies our numerical results. Meanwhile the Lyapunov spectrums of chaotic oscillators of Figure 4, considering electronic design requirements, are introduced here.
6. Conclusions
A straightforward algorithm to compute Lyapunov exponents of PWL functions-based chaotic systems with a multidirectional orientation, specifically Chua’s circuit, a 1D 6-scroll chaotic oscillator, and a 2D 4-scroll chaotic oscillator, has been reported. The method was based on the solution of the-PWL variational systems generated from splitting the PWL functions in several linear regions and collecting them in a pair of matrices, called modified state-variable equation. Also, it has been demonstrated that the Jacobian matrices are constants for all computation cycles and only one variational system is solved in each integration step-size. In this manner, the computation of LEs for OpAmp-based PWL chaotic oscillators considering electronic design aspects was found out to be reliable as demonstrated by the numerical simulations of five different cases. Besides, small errors ensuring the accuracy of the algorithm were obtained. From an experimental standpoint, design tradeoffs between the frequency bandwidth and the renormalization factor used to avoid incorrect values of LEs have been observed. Taken into account other effects of electronic circuits when Lyapunov spectrum is estimated, such as slew-rate, input offset, and common-mode rejection ratio, could be a future direction of investigations.
Acknowledgments
This work has been partially supported by SEP PROMEP, CONACyT/MÉXICO, and VIEP-BUAP under Projects BUAP-PTC-359, 174330 and 2013-VIEP Grants, respectively.