Simulink implementation of a new optoelectronic integrated circuit: stability analysis and infinite-scroll attractor

In this paper a 6-D optoelectronic system consisting of an optical injected semiconductor laser driven by a resonant tunneling diode is reported. A stability analysis of the hybrid system is analytically and numerically performed and paramount role of the effective gain coefficient is stuck out in the framework of new stability control. As a result, this parameter allows improving the accuracy of the stability study by circumscribing locked and unlocked regions. Besides, a narrow area of stability is pointed up within the sea of unstable points from which a complex fractal attractor so-called infinite-scroll attractor is highligted. Thereby, Simulink shows generation effectiveness of infinite-scroll attractor erratically interpersed by laminar phases. Also dynamics of Lyapunov exponents has confirmed that it refers to a strange fractal attractor. Moreover chaos control is structurally carried out by direct current polarisation.


Introduction
Chaos is a rich nonlinear phenomenon characterized by interesting properties such as unpredictability, ergodicity, mixing property, complicated structure dynamics, deterministic dynamics, high sensitivity to initial condition (Abdullah et al. 2018) to name a few. Many natural and non-natural systems are commonly modeled by nonlinear differential equations exhibiting chaos. So, nonlinear process plays as a cornerstone in developing and understanding novel complex systems as well it has received significant attentions in various fields. Chaotic systems have been also attracted in several scientist fields and usually play a relevant role improving their performances. There are many reasons why the nonlinear dynamics have been intensively studied in recent years. For example in physics to mention merely a few, nonlinear dynamics offers a great opportunity to improve memristive systems (Wang et al. 2017;Ma et al. 2018), circuit-systems (Li and Yao 2016;Ma 2014;Sun et al. 2011) lasers (Mengue and Essimbi 2012a;Daly et al. 2013) RTD optoelectronic based-systems (Figueiredo 2008), etc. In recent years, several studies have been reported in the field of optoelectronic related to resonant tunneling diode-laser diodes (RTD-LD) systems with very interesting potential in improving laser features (Romeira 2008(Romeira , 2013Romeira et al. 2017;Ironside et al. 2019). However, the two-rate equations have been used to model the LD dynamics. Recently, Essebe et al. (Essebe et al. 2020) have taken into account the modified RTD-LD systems including a semiconductor laser with optical injection modeled by Mengue and Essimbi (Mengue and Essimbi 2011). This has lead to a 6-D optoelectronic system with furthermore complex dynamical behavior such as the generation of multiscroll chaotic attractors, coherent collapse, as well as hyperchaotic dynamics controlled by the effective gain coefficient and other parameters (Essebe et al. 2020).
As well-known, in nonlinear chaotic systems, a slight variation of parameters or initial conditions can give rise to different outputs dynamics due to instabilities for example (Abdullah, et al. 2018). In addition in lasers' dynamics, optical injection, optical feedback, electrical current modulation or other external excitation provides additional degree of freedom (Ohtsubo 2006;Mengue and Essimbi 2012b). Hence, despite interesting behavior which occurs in laser output due to additional degree of freedom, the system may be subject to undesired phenomena such as relaxation oscillations, coherence collapse, instabilities, bursting oscillation (Mengue and Essimbi 2012a), etc. that can deteriorate the laser outputs (Ohtsubo 2006). Stability analysis may provide suitable conditions improving the features of those systems. Understanding instabilities of nonlinear systems is important to access their performance or to control their dynamics. Stability analysis deals better with nonlinear systems (Daly et al. 2013;Ohtsubo 2006;Mengue and Essimbi 2012b). The stability analysis control includes the phase portrait analysis, Lyapunov stability theory using eigenvalues or Lyapunov exponents and Routh-Hurwitz criterion to name a few.
In this paper, the Routh-Hurwitz stability criterion, Simulink implementation and Lyapunov exponents' dynamics are applied to various control parameters to investigate the oscillatory behavior of this system. Using theoretical and numerical calculations, we investigate the effects of main parameters such as DC bias, optical injection strength and the effective gain coefficient (EGC) on dynamical behavior. The purpose of this paper is to study the stability analysis and infinite-scroll attractor behavior of this novel model. The control of the system is provided through bifurcation diagram to furthermore improve functionalities of RTD-LD oscillators. In recent years, Matlab-Simulink environment has taken growing consideration in the simulation and implementation of chaotic systems (Abdullah et al. 2018) and in designing new chaotic maps for secure communication/transmission (Hamsa 2018) with accuracy. So, in this paper, Matlab/ Simulink implementation will be built to understand well the formation of infinite scroll attractors.
The paper proceeds as follows: in Sect. 2 we present the optoelectronic model equations and the stability analysis according to Routh-Hurwitz criterion. Section 3 provides its infinite-scroll attractor and Simulink implementation. Section 4 discusses the results and the effect of main parameters of dynamical behavior and the last section is devoted to conclusion.

Optoelectronic model
Recently starting from the classical RTD-LD Liénard oscillator (Figueiredo 2008;Romeira 2008Romeira ,2013Romeira et al. 2017), we have derived the time dependent equations of a novel 6-D optoelectronic hybrid oscillator by replacing the 2-D LD by 3-D LD from Essebe et al. (2019). The model's circuit is presented in Fig. 1. From Ref. (Essebe et al. 2020), the normalized dimensionless differential equations of the system are given by Eqs. (1) Let us note that, the triplet (n, s, ) represents the optical output subsystem where n, s and are normalized carrier density, photon density and optical phase respectively. Equation (1) is obtained rescaling parameters as follows: is the normalized frequency of external AC perturbation. G 0 = G∕ 0 , with G the normalized differential gain, LC is the normalized natural angular-frequency of forced RTD circuit, refer to normalized frequency detuning and injection strength respectively. The parameter is the mean EGC (Essebe et al. 2020). The RTD consists of 2 AlAs barriers each 2 nm thick and 6 nm InGaAs quantum with measured values of valley-to-peak voltage difference of around 0.8 V and peak current densities up to 18 kAcm −2 . The device emitted in the 1550 nm optical communications window and can achieve up to 28 dB optical modulation in a 200 μm active length device area InGaAlAs RTD, with 500 nm as the width of the low doped layer on the collector side of the RTD (Figueiredo et al. 2001). The tunable semiconductor diode Laser with external cavity from 635 to 2450 nm, output power up to 1mW and high fiber coupling efficiency (Mengue and Essimbi 2011). (1) sin ( (t)) − ΔΩ

Analytical study of stability
The N-shape of the RTD is approximated by a current-voltage cubic law which presents a negative differential resistance region (Essebe et al. 2020), as depicted in Fig. 2. Here, we denote the first positive differential resistance PDR1, the negative differential resistance NDR and the second positive differential resistance PDR2. Also, the area (between the curve and V-axis) under the RTD current-voltage characteristic is S 1 and that above is S 2 . Later we shall discuss the dynamic of the system when S 1 ≈ S 2 (Fig. 2a), (Fig. 2c). The number of equilibria points will be analyzed as follows. Let ̇x =̇y =̇z =̇s =ṅ =̇= 0 , so we get the non-trivial steady states. This leads to numerical calculation of equilibria points solving the Eqs.
(2) giving the relationships between fixed points successively is the magnitude of the complex optical field),s 0 must always be a positive value i.e. y 0 > n 0 . Numerical calculation of Eqs.
(2) shows that this conditions is satisfied if and only if the frequency detuning obeys ΔΩ ≤ 2GHz i.e.ΔΩ ≤ 0.41952 . In addition, according to the boundary of sine function, the necessary condition for optical locking requires the following relation For constant biased laser, we found in Ref. (Mengue and Essimbi 2012b) that all the points situated in space (ΔΩ = ± k∕ √ s 0 ) inside two curves are in locking region. However for the described driven laser, the condition ΔΩ ≤ 0.41952 is added to insure that the photon number s 0 should be positive. Then possible locking regions can be found as presented in Fig. 3. From Eqs. (2), z 0 is always equal to zero due to the t-component from non-autonomous system. The steady state solutions x 0 and y 0 can be numerically calculated, and therefore we can obtain the corresponding n 0 , s 0 and 0 by solving numerically Eqs.
The Jacobian matrix of the system around the equilibrium set of point The corresponding characteristic equation of the system (3) can be written as The coefficients a 1 , a 2 , a 3 , a 4 and a 5 are given in Appendix. The Hurwitz matrix associated to Eq. (4) is given by the Eq. (5).
Equation (4) can then be solved in order to capture the type of fixed point stability. In the same way, from the Routh-Hurwitz criterion, the equilibria M 0 x 0 , y 0 , 0, n 0 , s 0 , 0 are stable if all roots of Eq. (4) have negative real parts; otherwise they are unstable fix points. The necessary condition for stability is satisfied if and only if all the principal minors of (5) are positive. Since the coefficient a 0 = 1 , this implies the following Routh-Hurwitz criterion stable locked solutions condition: In this stage, we can solve Eqs. (6) graphically using Maple software and the result is depicted in Fig. 4.

Dissipativity
The divergence of the system (1) is It is hard to directly determine this dissipativity because of its dependence to the system variables (x, n, s) . However, as a , b , B 0 , m , r , (with < 1) , r , and G 0 are positive, the variables (x, n, s) have always positive values and n is nearly equal to one. Hence, it clearly appears that the divergence ∇V = − (x, n, s) is less than zero for convenient values of (x, n, s) such as − x 0 = 1, n 0 = 0.99, s 0 = 6.4751 = −0.78516 wherein s 0 is determined by using Eq. (2).
In this respect, the system (1) could be dissipative (as we have shown in tracking section that ≈ 0.87 ) and the element volume of the system is contracted at time t by the flow in the volume V 0 e − t . This means that the asymptotic solution settles into an attractor in the 6-D phase space. We denote that the selected maps will not be the 6-D, 5-D or 4-D but the projection of the attractor of the system in 2-D or 3-D spaces generated by adjusting selective control parameters such the EGC, frequency detuning and r-parameter.

The Infinite-scroll attractor and Matlab/Simulink implementation
In this section, we use the fourth order Runge-Kutta (RK4) in the rate equation system (1) to generate the infinite scroll attractor and verify its effectiveness using MATLAB/ Simulink implementation.
Let us recall that there is not non-chaotic system that exhibits multiple scrolls. Multiscroll attractors can be found either the system is chaotic or hyperchaotic (Ma et al. 2018). In this section, we set the parameters as = 0.0264 , r = 6.65, = 0.68 , = 1∕i th ≈ 196, = 0.1512 , v dc = 1.145 v and v ac = 0.1441 and the system (1) being under relative strong external optical injection k = 0.020637 . We choose the normalized time varying and numerical implementations are implemented using the RK4 algorithm in Matlab. After an initial transient, the solution settles into irregular and non-periodic oscillation that never repeats exactly in time. The trajectories in phase diagram when the optical phase (t) is plotted against the carrier density n(t) in Fig. 5 show that, an attractor which the number of scroll that increases while the calculating time increasing. The system generates: two-scroll, five-scroll, and forty one-scroll attractors as the time is increasing, respectively in Fig. 5a-c.This stands for the infinity-attractor. We recall in addition that, according to refs. (Ma 2014) and (Ma et al. 2018), the number of scroll increases the complexity of the systems; it is useful that this number be stabilized to a fixed number. This fact is possible to be carried out by choosing > 0.68 when others parameters are fixed as indicated above or by decreasing the optical injection parameter k.
The MATLAB/Simulink implementation is used in addition in normalized system (1). A complex scheme is designed in a chip simply by choosing essentially small building blocks of Simulink such as multipliers, gains, summations, constants, functions, integrators, etc. as shown in Fig. 6. The behavior of (1) through Simulink scheme is investigated for same parameter values as in Fig. 5.
The first step is to verify the generation of infinite-scroll attractor. Simulink observations referring to the formation of infinite-scroll attractor, in the − n plane from XY Graph3 is depicted in Fig. 7a. Figure 7a (i) shows a three scroll attractor for t = 2,226, Fig. 7a (ii) presents a four scroll for t = 2,410, a five scroll attractor in Fig. 7a (iii) for t = 3,000 and a six scroll attractor in Fig. 7a (iv) for t = 3,360. It is shown that, the number of scrolls generated increases by the time of implementation. The Simulink observations confirm the RK4 numerical simulation results and the effectiveness of the infinite scroll attractor.
Second step consists to present optoelectronic attractors for other different parameters. It is interesting to discuss the optoelectronic phase portrait using the triplet (x, n, s) which is implied in the divergence in Sect. 3.2. Figure 7b-d depict the current-voltage, the light-voltage (optoelectronic) and the carried-photon projection of the chaotic attractors respectively.
Let us mention in addition that, the chaotic behavior of the infinite-scroll attractor is confirmed in respect to the frame of Lyapunov exponent dynamics. As depicted in Fig. 8, we have at least one positive Lyapunov exponent so; the system (1) the generated infinitescroll is chaotic.

Improved linear stability and dynamics of infinite-scroll
In this section, the stability properties of the optically injected semiconductor laser (LD) driven by a resonant tunneling diode (RTD) current is numerically investigated and  Page 11 of 16 388 dynamic implementation of these stability properties is performed.
In Fig. 3a it is shown that the locking regions extend when EGC increases from a certain value. EGC acts as a control stability parameter and enables to improve the accuracy of the stability study further. From Fig. 3b almost all possible values of 0 are negative. Their number should be larger as decreasing. Additionally, under weak injection strength and according to Eq. (2), 0 min,max = ± ∕2 − arctan( ∕ ) and under strong injection, 0 min,max = ± ∕2 similar to the Ref. (Daly et al. 2013). As = −1 − 1 , the stable phase condition at the boundaries of EGC will be such that, when 0 min,max ( < 0) = ± ∕2 and 0 min,max ( < 1) = − or 0. Note that, according to Eq. (6) the stable region corresponding to the frequency detuning in locking region shall be given from voltage-carrier steady states x 0 n n 0 as all states can be connected from Eq. 2. In Fig. 4 we illustrate the Routh condition for the first minor Fig. 4a, the second minor Fig. 4b, the third minor Fig. 4c, the forth minor Fig. 4d and the fifth minor Fig. 4e from the Hurwitz matrix (5). Figure 4f illustrates the superposition of aforementioned figures corresponding to all combined parts of Eq. (6). It is obvious to observe that, Fig. 4f does not display a wide intersecting region which satisfies Eq. (6). Nevertheless, it exists a shrunk region where the system has stable fixe points; this tiny region is circumscribed by a stability boundary denoted by the line A in Fig. 4f. Thereby, the study reveals a shrunk area of the system from which nonlinear dynamics can be undertaken as attraction basin and a sea of points in which the system drops to instability.
Let us in addition mentioned that, the local stability of all steady-state points for a given set of decisive parameters values v dc , r, can be evaluated from Eq. (4) directly. For the simplicity we set ΔΩ ≈ 0 according to the limit of locking region for low frequency detuning as depicted in Fig. 3.The pick of the control parameter for stability response through eigenvalues of the Jacobian matrix J M 0 for different equilibria points are illustrated in Table 1. Each equilibrium point has at least two complex-conjugate eigenvalues. Then Routh criterion applied to different control parameters indicates that all fixed points are almost unstable, and the system (1) may exhibit either chaotic or hyperchaotic dynamics. However for appropriate pick of control parameters, there exist stable fixed points. For instance when the RTD is biased in NDR, stable fixed points A 1 and A 2 the Table 1 have x 0 and n 0 values which effectively belong to the shrunk stable region too close to the line A on in Fig. 4f. So, it is obviously to point out that stable fixed points A 1 and A 2 also belong to stable phase locking regions from Fig. 3.
Moreover, in Fig. 5, when the time increases, optical phase values upgrade and the number of attractor scrolls infinitely grows. With this framework, the system depicts two-scroll attractor in Fig. 5a for (t = 1,000), five-scroll attractor in Fig. 5b for (t = 3000) and forty one-scroll attractor in Fig. 5c for (t = 20,000). We notice that scrolls are interspersed to laminar phases which firstly appear after the fifteen scroll, as the time evolves they erratically appear and after a relative longtime their appearance becomes flat before falling again. This aforementioned study proves that it refers to a so-called an infinite-scroll chaotic attractor therein the number of generated scrolls depend to time calculation and non-periodic reproduction of infinite-scroll as the time evolves in regard to laminar phases. So the system (1) generates an infinite-scroll attractor. Indeed the transformation M 0 x 0 , y 0 , z 0 , n 0 , s 0 , 0 → M 0 x 0 , y 0 , z 0 , n 0 , s 0 , 0 + 2k ,where k = 0, ± 1, ± 2, ± 3, … is an invariant since trigonometric functions involved in optical phase give rise to infinitely many equilibrium points that provide infinite-scroll attractor. Nowadays, infinite-scroll attractor is not mostly reported in nonlinear dynamical systems for the best knowledge of the authors.

Influence of the areas of RTD current-voltage curve and chaos control bifurcation
In this section, we investigate the dynamical behavior of the system for three different cases; when S 1 ≈ S 2 , S 1 < S 2 and S 1 > S 2 . For the comparison purposes, we investigate the phase portrait diagrams in (n − s) plane when others parameters are set as in Ref. (Essebe et al. 2020). Throughout all above sections, the case S 1 ≈ S 2 has been used. It is important to show by the current-voltage characteristic of the RTD; there is invariance by changing a-parameter to b-parameter. In Fig. 9, we give prediction of the driven laser dynamics firstly when S 1 ≈ S 2 in Fig. 9a, secondarily for S 1 < S 2 in Fig. 9b and finally when S 1 > S 2 in Fig. 9c. We use in addition the EGC as bifurcation parameter in order to magnify the influence of the areas S 1 and S 2 in the system (1). It is shown in Fig. 9a that the complexity of (1) decreases while is increasing. For example, the laser exhibits three-scroll chaotic attractor for = 0.6 in Fig. 9a(i), two-scroll attractor for = 0.75 in Fig. 9a(ii) and fails down unto one-scroll attractor for = 0.9 in Fig. 9a(iii). However in the case S 1 < S 2 we can observe relaxation dynamics for = 0.6 in Fig. 9b(i), a non-oscillatory dynamics for = 0.75 in Fig. 9b(ii) and crisis when = 0.9 in Fig. 9b(iii). When S 1 > S 2 , the system generates merely periodic limit cycles whatever the value of the EGC as depicted in Fig. 9c(i), 9c(ii), and 9c(iii). We find very important to highlight that as one of the area S 1 gradually increasing the system evolves from chaos to periodic dynamics as depicted from Fig. 9a and c. Nevertheless a more careful observation of Fig. 2 reveals that the overall shape of RTD current-voltage characteristics is not changed but everything suggests that a shift in the axes has occurred. Therefore, interestingly the Fig. 9 Influence of RTD current-voltage (I-V) characteristics with EGC . a-When S 1 ≈ S 2 : a(i) Threescroll chaotic attractor, a(ii) Two-scroll attractor and a(iii) One-scroll attractor. b-When S 1 < S 2 : b(i) relaxation dynamics, b(ii) non-oscillatory dynamics and b(iii) crisis. c-When S 1 > S 2 , the system generates periodic limit cycles whatever the value of . c(i), c(ii), and c(iii) stabilization resistor and voltage biased could be adjusted to get a desired dynamics such as chaos, hyperchaos or periodic oscillations resulting from each of the investigated cases namely S 1 ≈ S 2 , S 1 < S 2 and S 1 > S 2 according to the threshold characteristics of the semiconductor laser.
Since chaos has a strong dependence on the initial conditions, its time-future oscillations being unforeseen, chaos control becomes a very useful tool in nonlinear dynamics. In this section we also apply an external control using the polarization control technique to stabilize electrical and optical outputs so as to control chaos. We investigate to the bifurcation diagram of y-state variable with the variation of DC polarisation v dc and fixed AC excitation v ac . From Fig. 10 we can notice that periodic and aperiodic oscillations occur in electriccal current Fig. 10a and in optical output Fig. 10b if the system is biased in the NDR of the RTD. When the DC biasing is beyong the NDR region (PDR1 or PDR2 region), the dynamical behavior becomes periodic. However, the system could exhibite either slight chaotic or periodic windows in NDR region. This sticks out the relevance of DC polarisation on the chaos control in optoelectronic systems.

Conclusion
We have investigated stability analysis of a novel OEIC in the framework of weak and strong optical injection. An analytical study has been performed from Hurwitz matrix and by using Routh-Hurwitz criterion to highlight some stable fixed points within the shrunk stable region inside the sea of unstable points. Besides, dissipativity of the system has been analytically examplified to firm divergence of Lyapunov exponents and nonlinear dynamics characteristic of the system. Numerical implementation has allowed to show that the laser locking regions are drastically affected by slight variations of EGC-parameter and DC bias. Moreover, from numerical calculations, the system displays stable and unstable equilibria points for negative or low frequency detuning values. Later on, generation and complexity of an infinite-scroll chaotic attractor has been discussed and stick out by using Lyapunov bifurcation diagrams of electrical driving current and optical output, respectively exponents' dynamics and Simulink implementation. At last, RTD current-voltage characteristics is studied under slight variations of control parameters namely the areas under and above the curve, further bifurcation diagrams and phase portraits are used to control DC polarization of the system. The momentous role played by EGC ( a new decisive parameter for this kind of OEIC) so as to improve the accuracy of stability analysis and the chaos control outlined prove that this study can be relevant for experimental implementations of optoelectronic devices.