Double-zero degeneracy and heteroclinic cycles in a perturbation of the Lorenz system

Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 Contents lists available at ScienceDirect Communications in Nonlinear Science and Numerical Simulation journal homepage: www.elsevier.com/locate/cnsns Research paper Double-zero degeneracy and heteroclinic cycles in a perturbation of the Lorenz system A. Algaba a, M.C. Domínguez-Moreno a, M. Merino a, A.J. Rodríguez-Luis b,∗ aDepartamento de Ciencias Integradas, Centro de Estudios Avanzados en Física, Matemática y Computación, Universidad de Huelva, 21071 Huelva, Spain bDepartamento de Matemática Aplicada II, E.T.S. Ingeniería, Universidad de Sevilla, Camino de los Descubrimientos s/n, 41092 Sevilla, Spain article info Article history: Received 17 September 2021 Received in revised form 8 March 2022 Accepted 24 March 2022 Available online 30 March 2022 Keywords: Lorenz system Normal form Double-zero bifurcation Global connections abstract In this paper we consider a 3D three-parameter unfolding close to the normal form of the triple-zero bifurcation exhibited by the Lorenz system. First we study analytically the double-zero degeneracy (a double-zero eigenvalue with geometric multiplicity two) and two Hopf bifurcations. We focus on the more complex case in which the doublezero degeneracy organizes several codimension-one singularities, namely transcritical, pitchfork, Hopf and heteroclinic bifurcations. The analysis of the normal form of a Hopf-transcritical bifurcation allows to obtain the expressions for the corresponding bifurcation curves. A degenerate double-zero bifurcation is also considered. The theoretical information obtained is very helpful to start a numerical study of the 3D system. Thus, the presence of degenerate heteroclinic and homoclinic orbits, T-point heteroclinic loops and chaotic attractors is detected. We find numerical evidence that, at least, four curves of codimension-two global bifurcations are related to the triple-zero degeneracy in the system analyzed. ©2022 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/). 1. Introduction It is familiar that the dynamics exhibited by a three-dimensional autonomous system can be extremely complicated. The best known example and one of the most studied is the Lorenz system [1,2]. Although it was derived from a simple model of convection in the atmosphere, it also appears in many other fields (see, for instance, [3–11]). Due to their physical meaning, in most cases, its three parameters take positive values. However, in some situations, it also makes sense for them to be negative (see, for example, [6,7]). Therefore, it is interesting to consider that the three parameters can take any real value. Despite the many works dedicated to knowing and explaining its extraordinarily complex dynamical behavior (for instance, we cite some of them devoted to chaotic behavior [12–18], heteroclinic loops called T-points [19–21], manifolds study [22–24], invariant algebraic surfaces [25,26] and resonances and periodic orbits [27,28]; see also references therein) there are still many aspects to be understood. Although numerical tools are essential to be able to carry out its study, the analysis of the local bifurcations of the equilibria provides valuable initial information. In this way, the full analysis of the pitchfork, Hopf and Bogdanov–Takens bifurcation exhibited by the origin of Lorenz system has already been carried out (see, for instance, [29–32]). ∗Corresponding author. E-mail address: [email protected] (A.J. Rodríguez-Luis). https://doi.org/10.1016/j.cnsns.2022.106482 1007-5704/©2022 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons. org/licenses/by-nc-nd/4.0/). A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 However, the study of the Hopf-pitchfork and the triple-zero degeneracies presents an additional difficulty since the origin is not an isolated equilibrium. In this case, it would be necessary to adapt the standard techniques of normal forms of the reduced system on the center manifold. This complication has encouraged us to start the study of the triplezero bifurcation for systems invariant to the change (x,y,z)→(−x,−y,z) (symmetry exhibited by the Lorenz system), considering the corresponding normal form [33]. The triple-zero bifurcation is of codimension three when the eigenvalue has geometric multiplicity two. The study of this case, both in systems without symmetry (such as Rössler’s) and in Z2-symmetric systems (such as Chua’s) illustrates the richness and complexity of the dynamics that appears as a consequence of this degeneracy (see [34–36] and references therein). However, if the triple-zero eigenvalue has geometric multiplicity one, then the bifurcation is of codimension five [33], which reduces to codimension three if the above symmetry condition is assumed. Given the enormous difficulty of this problem, as a first step, we propose an unfolding with three parameters (see Eq. (3)). This system is of interest to us not only for small values of the parameters (on which we can do a local analysis) but also for large values (which are reached by numerical continuation). In fact, this system represents a wide family of Lorenz-like systems studied in the literature [37–42] (some of them well known, such as the Shimizu–Morioka system), although not only for small values of the parameters. Concretely, Kokubu and Roussarie demonstrate in this system the existence of non-degenerate heteroclinic cycles, which give rise to complex dynamical behaviors, in a region of the parameter space where there are small, unity order, and large parameters [42, Theorem 2.3]. We are interested in finding heteroclinic cycles in other regions of the parameter space and, in particular, to see its possible relationship with the triple-zero degeneracy. A starting point to study the existence of heteroclinic cycles is the analysis of codimension-two local bifurcations which are close to the triple-zero. In this way, we find in the unfolding three local bifurcations of codimension two, Bogdanov– Takens (double-zero eigenvalue with geometric multiplicity one), double-zero (eigenvalue with geometric multiplicity two) and Hopf-zero. The Bogdanov–Takens bifurcation has already been analyzed [43]. The first part of this work is devoted to studying the double-zero bifurcation. The global connections found, whose existence is guaranteed by the local analysis (Bogdanov–Takens and double-zero bifurcations), can be numerically continued, with which it is possible to study their degeneracies and their relationship with the triple-zero degeneracy (see Fig. 10). The corresponding results appear in the second part of this paper. The final objective that we intend to achieve in this long process is the study of the triple-zero bifurcation in Lorenzlike systems obtained by perturbations of the Lorenz system, in such a way that an isolated equilibrium undergoes the triple-zero degeneracy. Thus, the Lorenz system will appear as the limit of a certain family of Lorenz-like systems (we trust that this procedure will allow to obtain interesting information about its dynamics). The paper is organized as follows. In Section 2we introduce the system under study and determine the local bifurcations it can exhibit, namely pitchfork, transcritical, Hopf, Bogdanov–Takens, Hopf-pitchfork, double-zero and triplezero. The study of the double-zero degeneracy (double-zero eigenvalue with geometric multiplicity two) is carried out in Section 3. The presence of a nonlinear degeneration in the double-zero bifurcation is also discussed (some useful lemmas appear in Appendix A). Moreover, Hopf bifurcations of the nontrivial equilibria are analyzed and expressions for the loci in the parameter space where the first Lyapunov coefficient vanishes are provided (all the details appear in Appendix B). In Section 4we conduct a numerical study starting near the double-zero bifurcation, taking advantage of the analytical results previously obtained. In this way we detect a very complex bifurcation scenario where, among other things, degenerate homoclinic and heteroclinic connections, T-point heteroclinic loops and chaotic attractors are present. We also show numerical evidence that at least four curves of codimension-two global bifurcations are related to the triple-zero degeneracy present in the unfolding analyzed. Finally, a section with conclusions is included. 2. A normal form for the Lorenz system The Lorenz system is given by (see [1,2]) ⎧ ⎨ ⎩ ˙ x=σ(y−x), ˙ y=ρx−y−xz, ˙ z= −bz +xy, (1) where σ,ρand bare real parameters. These equations are invariant under the change (x,y,z)→(−x,−y,z); therefore, the z-axis is invariant. The equilibria are the origin (0,0,0) and, when b(ρ−1) >0, a pair of symmetric nontrivial equilibria (±√b(ρ−1),±√b(ρ−1), ρ −1). The origin can exhibit the following linear degeneracies: •a pitchfork bifurcation when ρ=1, σ= 0,−1, b= 0 [32, Sect. 3]; •a Hopf bifurcation for σ= −1, ρ > 1, b= 0 [30]; •a Bogdanov–Takens bifurcation when σ= −1, ρ=1, b= 0 [31]; •a Hopf-pitchfork bifurcation for σ= −1, b=0, ρ > 1; •a triple-zero bifurcation if σ= −1, ρ=1, b=0. Unfortunately the above Hopf-pitchfork and triple-zero bifurcations cannot be analyzed by the usual techniques because the origin is a non-isolated equilibrium when b=0. To circumvent this obstacle in the case of the triple-zero 2 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 bifurcation, in this paper we consider a three-parameter unfolding, that is close to the normal form of the triple-zero bifurcation exhibited by the Lorenz system [33], given by ⎧ ⎨ ⎩ ˙ x=y, ˙ y=ε1x+ε2y+Axz +Byz, ˙ z=ε3z+Cx2+Dz2, (2) where ε1, ε2, ε3≈0 and A,B,C,Dare real parameters. System (2) exhibits a triple-zero bifurcation when ε1=ε2=ε3= 0. These equations are also invariant under the change (x,y,z)→(−x,−y,z). We remark that several systems studied in the literature appear as particular cases of (2) for certain parameter choices. Thus, if A= −1, C=1, B=D=0 we obtain the Shimizu–Morioka system [37,38] (for ε1=1, ε2= −λ,ε3= −α) and a low-order model of magnetoconvection [39] (for ε1= −λ,ε2=k,ε3= −1). If A= −a k,C=h,B=D=0, by means of the change x=˜ x,y=a( ˜ y−˜ x), z=˜ z, a Lorenz-like system is obtained [40,41] (for ε1=a b,ε2= −a, ε3= −c). Finally, if A= −1, C=1, B= −γ,D=δ,ε1=˜ A,ε2= −1, ε3= −β, the resulting system [42, Eq. (2.7)] has non-degenerate heteroclinic cycles that connect the equilibria located on the z-axis, in certain regions where ε1is sufficiently large, ε2= −1 and ε3,B<0 and D>0 are close to zero (see [42, Theorem 2.3]). If A C = 0, since system (2) is invariant to the changes (x,y,z,t, ε1, ε2, ε3,A,B,C,D)−→ (−x,y,−z,−t, ε1,−ε2,−ε3,−A,B,C,D), (x,y,z,t, ε1, ε2, ε3,A,B,C,D)−→ (−x,y,z,−t, ε1,−ε2,−ε3,A,−B,−C,−D), we can take without loss of generality A= −1,C=1 ⎧ ⎨ ⎩ ˙ x=y, ˙ y=ε1x+ε2y−xz +Byz, ˙ z=ε3z+x2+Dz2. (3) System (3) can have up to four equilibria, namely E1=(0,0,0), E2=(0,0,−ε3/D)if D= 0 and E3,4= (±√−ε1(ε3+Dε1),0, ε1)if ε1(ε3+Dε1)<0. Note that E1and E2are placed on the z-axis, that is an invariant set. The Jacobian matrix of system (3) at the trivial equilibrium E1is (010 ε1ε20 0 0 ε3),(4) and then, its characteristic polynomial is given by P(λ)=λ3+p1λ2+p2λ+p3, where p1= −(ε2+ε3),p2=ε2ε3−ε1,p3=ε1ε3. Therefore, the origin exhibits the following bifurcations: •a pitchfork bifurcation when ε1=0, ε2= 0, ε3= 0. The equilibria E3,4appear when ε1(ε3+Dε1)<0. •A transcritical bifurcation (involving also E2) when ε3=0, ε1= 0, ε2= 0, D= 0. •A Hopf bifurcation when ε2=0, ε1<0, ε3= 0. •A Bogdanov–Takens bifurcation (a double-zero eigenvalue with geometric multiplicity one) when ε1=0, ε2=0, ε3= 0. It is of homoclinic type when ε3<0 and of heteroclinic type if ε3>0. •A Hopf-zero bifurcation when ε2=0, ε3=0, ε1<0. •A double-zero bifurcation (a double-zero eigenvalue with geometric multiplicity two) when ε1=0, ε3=0, ε2= 0. •A triple-zero bifurcation (a triple-zero eigenvalue with geometric multiplicity two) when ε1=ε2=ε3=0. As a first step in the study of system (3), the Bogdanov–Takens bifurcation has been analyzed in [43] whereas the double-zero bifurcation is considered in this paper. The study of the Hopf-zero and the triple-zero bifurcations will be carried out in the future. Before starting with the study of the double-zero bifurcation, we are going to see that all the bifurcations that the equilibrium E2exhibits can be easily deduced from those of E1. Indeed, to study the equilibrium E2=(0,0,−ε3/D), D= 0, we translate it to the origin by means of the change x=¯ x,y=¯ y,z=¯ z−ε3 D, that transforms system (3) into ⎧ ⎨ ⎩ ˙ ¯ x=¯ y, ˙ ¯ y=(ε1+1 Dε3)¯ x+(ε2−B Dε3)¯ y−¯ x¯ z+B¯ y¯ z, ˙ ¯ z= −ε3¯ z+¯ x2+D¯ z2, (5) with D= 0. 3 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 Since system (3) is symmetric to the change (x,y,z,t, ε1, ε2, ε3,B,D)−→ (x,y,z−ε3 D,t, ε1+ε3 D, ε2−Bε3 D,−ε3,B,D),(6) it is immediate to obtain the stability and bifurcations of E2from the stability and bifurcations of E1. 3. A double-zero bifurcation of the origin The double-zero bifurcation exhibited by the equilibrium E1occurs when ε1=0, ε3=0, ε2= 0. Consequently, to study this codimension-two bifurcation, we fix the value of ε2=∆−1= 0 (we use ∆to facilitate the writing of the expressions that will appear later) and assume that ε1and ε3are close to zero. The results obtained in this section are summarized in the following theorem and in Fig. 1. Theorem 3.1. In system (3) the equilibrium E1exhibits a double-zero bifurcation when ε1=ε3=0,ε2=∆−1= 0. In a neighborhood of this codimension-two degeneracy, there are only periodic solutions (arising from a Hopf bifurcation) when D>0and ∆<0. In this case, the local codimension-one bifurcations involved are (see Fig. 1): (1) A transcritical bifurcation T, of equilibria E1and E2, when ε3=0. (2) Two pitchfork bifurcations P1, when ε1=0, and P2, for ε3= −Dε1+O(ε2 1). They involve E1or E2and E3,4. (3) A Hopf bifurcation hof the equilibria E3,4, when ε3= −2Dε1+O(ε2 1), if B = 2∆. It is supercritical when B >2∆and subcritical if B <2∆. (4) A heteroclinic loop between E1and E2, when ε1= − 1 2Dε3+∆(B−2∆)(D−∆) 2D2(2D−3∆)ε2 3+O(ε3 3),(7) if B = 2∆. It is attractive when B >2∆and repulsive if B <2∆. This loop is formed by two heteroclinic connections, one is structurally stable (since it goes from Ejto Eion the invariant z-axis) and the other one (from Eito Ejand placed outside this axis) is more relevant as it is structurally unstable: Heij indicates that the connection between Eiand Ejis outside the z-axis. Concretely, He21 exists when ε3>0and He12 if ε3<0. To study the double-zero bifurcation exhibited by the equilibrium E1first we apply the change of variables x=Y+Z,y=ε2Z,z=X,(8) that transforms system (3) into ⎧ ⎨ ⎩ ˙ X=ε3X+(Y+Z)2+DX2, ˙ Y= −∆ε1Y−∆ε1Z+∆XY +(∆−B)XZ, ˙ Z=∆ε1Y+(∆ε1+1 ∆)Z−∆XY −(∆−B)XZ. (9) Remark that the hyperbolic manifold of equilibrium E1(which is semi-hyperbolic at the double-zero) is attractive if ∆<0 and repulsive if ∆>0 (since ˙ Z≈(1/∆)Z). Moreover, note that system (9) is symmetric to the change (X,Y,Z)→(X,−Y,−Z), so the center manifold inherits this symmetry, that is, it will be of the form Z=YΨ(X,Y2), where Ψis a smooth function. Therefore, the reduced system on the center manifold is invariant to the transformation (X,Y)→(X,−Y), so Y=0 is an invariant curve of the reduced system. Consequently, it is enough to analyze the reduced system for Y≥0. To determine the bifurcations that arise from this singularity we begin by studying the reduced system up to second order on the center manifold Z=∆2XY +··· {˙ X=ε3X+DX2+Y2+··· , ˙ Y= −∆ε1Y+∆(1 −∆2ε1)XY +··· .(10) The change of variables X=−1 Dz−ε3 2D,Y=1 √|D|r,D= 0,(11) transforms (10) into {˙ r=µ1r+arz, ˙ z=µ2−sgn(D)r2−z2,(12) where µ1=−∆ 2D(2Dε1+ε3)+O(|ε1, ε3|2), µ2=ε2 3 4,a=−∆ D+O(|ε1, ε3|).(13) Therefore, sgn(a)= −sgn(∆)sgn(D), i.e., a>0 if ∆D<0 and a<0 if ∆D>0. 4 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 Fig. 1. Bifurcation set of system (3) in a neighborhood of the double-zero DZ exhibited by the equilibrium E1, for D>0 and ∆<0 when: (a) B>2∆; (b) B<2∆. The curves, according to Theorem 3.1, correspond to the following bifurcations: T, transcritical; P1and P2, pitchfork; h, Hopf (supercritical when B>2∆and subcritical if B<2∆); He12 and He21, heteroclinic loop (attractive if B>2∆and repulsive when B<2∆). To study system (12) we can take advantage of the analysis of the Hopf-saddle–node degeneracy performed in [44] (we remark that (12) would correspond to a Hopf-transcritical bifurcation since µ2=ε2 3/4 is always positive), whose Eq. (7.4.9) is {˙ r=µ1r+arz, ˙ z=µ2+br2−z2,(14) with b= ±1. In this system there are six possible configurations (see [44, Sect.7.4] and [45, Sect.20.7]), labeled as I (b=1, a>0), IIa (b=1, −1<a<0), IIb (b=1, a≤ −1), III (b= −1, a>0), IVa (b= −1, −1<a<0) and IVb (b= −1, a<−1). Between them, the cases of interest are IIa–IIb and III, due to the presence of periodic solutions. On the one hand, in case III, a Hopf bifurcation exists when µ2>0 and the periodic orbit disappears in a heteroclinic connection. On the other hand, in cases IIa–IIb, also a Hopf bifurcation exists when µ2<0. Comparing (12) and (14) we see that b= −sgn(D), that is, b= +1 if D<0 and b= −1 if D>0. Consequently, we deduce that system (12) is in case III when D>0 and ∆<0. Conversely, system (12) is in cases IIa–IIb, when D<0 and ∆<0. However, as µ2=ε2 3/4>0, there is no Hopf bifurcation (because, in cases IIa–IIb, it only exists when µ2<0). Note that the trivial cases I and IV (IVa–IVb), where only equilibria appear, occur when D<0 and ∆>0 and when D>0 and ∆>0, respectively. Therefore, when D<0 and ∆<0 (cases IIa–IIb), the following local bifurcations are present (see [44, Fig. 7.4.4]): 5 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 •A transcritical bifurcation when µ2=0, that is, when ε3=0. •A pitchfork bifurcation for µ2 1−a2µ2=0, that is, for ε1=0 and for ε3= −Dε1+O(ε2 1). As stated above, in this case there is no Hopf bifurcation and, consequently, no heteroclinic connection. In the rest of this section we will analyze system (12), when D>0 and ∆<0 (case III), by following the steps given in [44, Sect. 7.4] for the Hopf-saddle–node bifurcation (recall that our case would correspond to a Hopf-transcritical bifurcation). The equilibria of system (12) are (0,±√µ2)and ⎛ ⎝±√µ2−µ2 1 a2,−µ1 a⎞ ⎠when a2µ2> µ2 1, that in terms of the original parameters are given by (0,±ε3 2)and (±√−Dε1(ε3+Dε1) (∆2ε1−1)2+O(|ε1, ε3|3),−2Dε1−ε3 2(1 −∆2ε1)), ε1(ε3+Dε1)<0. Note that E1corresponds to the equilibrium (0,−ε3/2) and E2to (0, ε3/2). From the study carried out in [44, Sect.7.4] we can deduce that the following local bifurcations are present (see Fig. 1): •A transcritical bifurcation Tof the equilibria placed on the z-axis when µ2=0, that is, when ε3=0. •A pitchfork bifurcation for µ2 1−a2µ2=0 that is, P1for ε1=0 and P2for ε3= −Dε1+O(ε2 1). •A Hopf bifurcation hof the equilibria placed outside the z-axis when µ1=0, µ2>0, that is, when ε3= −2Dε1. However, the second-order terms are not enough to unfold neither the Hopf bifurcation nor the heteroclinic connection because, for the values of the parameter where both bifurcations appear, the system is integrable, i.e., it has an analytic first integral. In fact, for a=2, it has a polynomial first integral (only case considered in [44]). Hence, it has a continuum of periodic orbits that ends in a degenerate heteroclinic connection (it corresponds to a Hopf bifurcation in which all the Lyapunov coefficients are zero, i.e., a center singular point) [44, Fig. 7.4.9]. Thus, to unfold and calculate the curves of Hopf bifurcation and heteroclinic connections (involving the two equilibria located on the z-axis), as well as their stabilities, we need the reduced system up to third order. As the center manifold up to third order is Z=Y(∆2X+∆3Y2+(D−B+2∆)∆3X2), we obtain the reduced system up to order three {˙ X=ε3X+DX2+Y2+∆2XY2, ˙ Y= −∆ε1Y+∆(1 −∆2ε1)XY −[∆−B+∆2(D−B+2∆)]∆2X2Y−∆4ε1Y3.(15) The change X=−1 Dz−ε3 2D,Y=1 √D−∆2ε3 r,(16) transforms (15) into {˙ r= ˜µ1r+˜ arz +˜ cr3+˜ drz2, ˙ z= ˜µ2−r2−z2+˜ er2z,(17) where ˜µ1= −∆ε1−∆ 2Dε3+∆3 2Dε1ε3+∆2(B−∆) 4D2ε2 3+O(|ε1, ε3|3), ˜µ2=ε2 3 4+O(|ε1, ε3|3),˜ a=−∆ D+O(|ε1, ε3|),˜ c=O(|ε1, ε3|),(18) ˜ d=∆2(B−∆) D2+O(|ε1, ε3|),˜ e=2∆2 D+O(|ε1, ε3|). Applying Lemma A.1 of Appendix A, we obtain, modulo translation and scaling, ⎧ ⎪ ⎨ ⎪ ⎩ ds dτ= ˆµ1s+ˆ asw, dw dτ= ˆµ2−s2−w2+ˆ fw3, (19) where ˆµ1= ˜µ1+O(|ε1, ε3|3),ˆµ2= ˜µ2+O(|ε1, ε3|3),ˆ a=˜ a+O(|ε1, ε3|), ˆ f=2∆(B−2∆) 3D+O(|ε1, ε3|).(20) 6 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 Thus, in system (19), we can analyze both the Hopf bifurcation and the heteroclinic connection. On the one hand, the Hopf bifurcation exhibits a first-order degeneracy when B=2∆. To prove it, we obtain from system (19) the third-order reduced system (˙ s ˙w)=(0−1 1 0 )(s w)+(ˆ a √2ˆ as1 sw −ˆ a 2√2ˆ as1 s2−1 √2ˆ as1 w2+˜ f √2ˆ as1 w3),(21) where s1=√(ˆµ2ˆ a3− ˆµ2 1ˆ a−ˆ fˆµ3 1)/ˆ a3. By using the recursive algorithm developed in Ref. [46], we compute the normal form of the Hopf bifurcation up to third order ˙ r=a1r3+··· ,˙ θ=1+b1r2+··· ,(22) where the first Lyapunov coefficient is given by a1=3ˆ f/(8√2ˆ as1).Therefore, the Hopf bifurcation is supercritical when a1<0, i.e., if ˆ f<0, which occurs when B>2∆, whereas it is subcritical if a1>0, i.e., when ˆ f>0 or, equivalently, when B<2∆. Consequently, a degenerate Hopf bifurcation appears when a1=0, i.e., when B=2∆. On the other hand, when ˆ f= 0, that is, if B= 2∆, Melnikov’s method guarantees the existence of a heteroclinic connection whose curve is approximated by [47–49] (in [44] only the Hamiltonian case, ˆ a=2, is considered) µ1= − 3ˆ a2ˆ f 2(3ˆ a+2) µ2+O(µ2 2), that in terms of our parameters corresponds to ε1= − 1 2Dε3+∆(B−2∆)(D−∆) 2D2(2D−3∆)ε2 3+O(ε3 3).(23) Concretely, He21 exists when ε3>0 and He12 if ε3<0. Note that very recently, using the nonlinear time transformation method (NTT) [50–53], the expression of the heteroclinic curve in system (19) has been computed at any order for all ˆ a∈R[49, Eq.(3.20)]. With this we have completed the proof of the results stated in Theorem 3.1. Recall that all the curves mentioned in this theorem are drawn in Fig. 1. 3.1. Degenerate case B =2∆ Next we show the steps that should be followed to study the degenerate case B=2∆. In this situation, the Hopf bifurcation and the heteroclinic connection are not determined by system (19), and therefore we need the reduced system up to fifth order on the center manifold. Following the same steps as in the calculation up to order three, we obtain ⎧ ⎪ ⎨ ⎪ ⎩ ds dτ= ˆµ1s+ˆ asw+O(|s, w|4), dw dτ= ˆµ2−s2−w2+O(|s, w|4), (24) where ˆµ1,ˆµ2and ˆ aare given in (20). Applying Lemmas A.2 and A.3 of Appendix A, we obtain, modulo translation and scaling, ⎧ ⎪ ⎨ ⎪ ⎩ d¯ s dτ= ¯µ1¯ s+¯ a¯ s¯w, d¯w dτ= ¯µ2−¯ s2− ¯w2+¯ h¯w5, (25) where ¯µ1= ˆµ1+O(|ε1, ε3|3),¯µ2= ˆµ2+O(|ε1, ε3|3),¯ a=ˆ a+O(|ε1, ε3|),(26) and the expression for ¯ his not included because it is too long. When ¯ h= 0, the application of the NTT method provides the following expression for the heteroclinic curve [54] ¯µ1= − 15¯ a3¯ h 2(3¯ a+2)(5¯ a+2) ¯µ2+O(¯µ2 2), that in terms of the original parameters corresponds to ε1= − 1 2Dε3−15∆2¯ h 8D(2D−3∆)(2D−5∆)ε2 3+O(ε3 3). 7 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 3.2. Hopf bifurcations of E3,4and E2in the whole (ε1, ε2, ε3)-space The study made so far for the Hopf bifurcation of the equilibria E3,4is valid when ε1and ε3are small, close to the double-zero bifurcation. In Appendix B.1 we have completed the analysis of this Hopf bifurcation for any value of ε1and ε3. Specifically, we provide the value of the first Lyapunov coefficient in Eq. (B.4) and give an expression for the Hopf surface in the (ε1, ε2, ε3)-parameter space in Eq. (B.5). These results will be useful in the numerical study of system (3) carried out in Section 4. On the other hand, in Appendix B.2 we consider the Hopf bifurcation of the equilibrium E2which is related to a Bogdanov–Takens bifurcation of E2, as we will see in Section 4. In Eq. (B.6) the value of the first Lyapunov coefficient is given. 4. Numerical study We have seen in the previous section that the presence of a double-zero bifurcation in system (3) guarantees the existence of a heteroclinic loop which initially emerges as a planar phenomenon (in the same way that a Bogdanov–Takens bifurcation ensures the presence of homoclinic or heteroclinic connections). However, as soon as we move away from the point DZ, the heteroclinic connection will start to develop a three-dimensional structure. Thus, the numerical continuation of the heteroclinic curve away from the double-zero will make possible to determine if any point of degeneracy appears. This situation, in most cases, will imply a greater richness in the dynamics of the system. The objective of this section is twofold. On the one hand, in Section 4.1 we are going to study the bifurcations that appear from the degenerate double-zero bifurcation (B=2∆, where ∆=1/ε2), completing the local information provided in Theorem 3.1 and shown in Fig. 1. On the other hand, in Section 4.2, we will find regions of the parameter space in which the system (3) exhibits a very complex dynamics. 4.1. Degenerate DZ when B = −1and D =1 Taking advantage of the theoretical information obtained in Section 3, we are going to carry out a numerical study of system (3) by means of the continuation software AUTO [55]. Our goal is to obtain bifurcation sets in the vicinity of the double-zero bifurcation DZ exhibited by the equilibrium E1=(0,0,0) when (ε1, ε3)=(0,0) and ε2= 0. Specifically, we will take slices ε2=constant in the (ε1, ε2, ε3)-parameter space, fixing B= −1 and D=1. We know from Section 3that the double-zero bifurcation DZ is degenerate when ε2=2/B= −2. Then we will consider values on both sides of this singularity, for instance, ε2= −1.5 and ε2= −2.5, that correspond, respectively, to the cases B>2∆and B<2∆. 4.1.1. ε2= −1.5(case B >2∆) Initially, we are going to obtain the bifurcation set for ε2= −1.5. Thus, as it is usual, in order to continue numerically the bifurcation curves in a parameter plane we previously draw some bifurcation diagrams. Then, we also fix ε3= −0.5 and study the evolution of the nontrivial equilibrium E3=(√ε1(0.5−ε1),0, ε1)versus ε1. It exhibits a Hopf bifurcation hwhen ε1≈0.2297234. The bifurcation diagram corresponding to the asymmetric stable periodic orbit emerged from his drawn in Fig. 2(a). This periodic orbit disappears, when ε1≈0.2304689, in a heteroclinic cycle He12 between the saddle equilibria E1=(0,0,0) and E2=(0,0,0.5). Recall that this loop is formed by two heteroclinic connections, one is structurally stable (since it goes from E2to E1on the invariant z-axis) and the other one is more relevant as it is structurally unstable, of codimension one (the connection from E1to E2is placed outside this axis). For this reason this heteroclinic cycle is labeled with the superscript 12, that is, Heij indicates that the connection between Eiand Ejis outside the z-axis. The projection of He12 onto the (x,z)-plane appears in Fig. 2(c). Remark that throughout this work, with the aim of simplifying the notation, we will label the heteroclinic cycle (in fact, due to the symmetry, a pair of heteroclinic cycles exists) and the heteroclinic bifurcation with the same symbol, although they are two different objects. In addition, when necessary, we will use superscripts to indicate the equilibria that are involved in a certain bifurcation, or in a degeneration of it. If we fix ε3= −0.9, we obtain the bifurcation diagram shown in Fig. 2(b). The periodic orbit emerged from his of saddle type, later it becomes stable in the saddle–node bifurcation sn (ε1≈0.397233) and, finally, it disappears in a heteroclinic cycle He12 (ε1≈0.397237) similar to the previous case. Note that the heteroclinic loop is attractive for both values of ε3. Now we can numerically compute for ε2= −1.5 the loci where the bifurcations detected in Fig. 2 occur in the (ε1, ε3)-plane. Thus, in the partial bifurcation set drawn in Fig. 3(a), we can observe the curves h(Hopf bifurcation of the equilibria E3,4), He12 (heteroclinic cycle between the equilibria E1and E2) and sn (saddle–node bifurcation of periodic orbits). Moreover, three straight lines intersect at the double-zero point DZ, namely P1(pitchfork bifurcation of the origin, ε1=0), P2(pitchfork bifurcation of E2,ε1= −ε3) and T(transcritical bifurcation between E1and E2,ε3=0). Note that in these pitchfork bifurcations the equilibria E3,4emerge. Going into more detail, near the origin in the fourth quadrant the equilibria E3,4arise from the curve P1as attractive nodes (their three real eigenvalues are negative). By increasing the value of ε1they become attractive foci until they 8 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 Fig. 2. For ε2= −1.5,B= −1,D=1, bifurcation diagram of the asymmetric periodic orbit emerged from the Hopf bifurcation h, exhibited by the equilibrium E3, when: (a) ε3= −0.5; (b) ε3= −0.9. (c) Projection onto the (x,z)-plane of the heteroclinic cycle He12 of panel (a) that exists when ε1≈0.2304689. undergo a Hopf bifurcation h. On the other hand, in this quadrant the equilibrium E1is always a saddle–node (λ1,3< 0, λ2>0). This implies that, at any point in this region there exists, apart from the heteroclinic connection between E1and E2, another structurally stable (codimension-zero) heteroclinic connection between the equilibria E1and E3(see Fig. 3(b) for (ε1, ε3)=(0.15,−0.5)). When the Hopf curve his crossed, E3becomes a repulsive saddle-focus and then the unstable manifold of E1connects with the attractive periodic orbit arising from this Hopf bifurcation (see Fig. 3(c) for (ε1, ε3)=(0.2304,−0.5)). If ε1continues to increase, the periodic orbit grows and approaches the z-axis until the heteroclinic connection between E1and E2is formed (see Fig. 2(c)) for (ε1, ε3)≈(0.2304689,−0.5). According to Eq. (B.5) of Appendix B.1, for a fixed value ε2=constant, the Hopf bifurcation curve hthat emerges from the point DZ is given by the equation ε3= 3ε2 1−ε2 2+2ε1(1 −ε2)+√(ε1−ε2)4+4ε1(ε2 1+ε1−ε2 2) 2(−ε1+ε2). In the fourth quadrant this curve is unbounded and the following inequalities hold −2(ε1+1) < ε3<−2ε1. The first one corresponds to the condition p2>0 (a necessary condition for the existence of the Hopf bifurcation). The second inequality is due to the fact that the line ε3= −2ε1is an asymptote of the Hopf curve when ε1→ +∞ and ε2=constant. Moreover, in accordance with the bifurcation diagrams of Fig. 2(a)–(b), a degenerate point Dh appears on h, when (ε1, ε3)≈(0.373821,−0.842653), since the first Lyapunov coefficient a1of the normal form vanishes (see Fig. 3(a)). On the other hand, the curve of heteroclinic connections He12 is also unbounded. It emerges from the point DZ to the right of the Hopf curve hsince B>2/ε2(this agrees with the bifurcation diagram of Fig. 2(a)). However, when (ε1, ε3)≈(0.38281,−0.86445), both curves intersect and then, they change their relative position. Moreover, the curve He12 exhibits a degeneracy DHe12 at (ε1, ε3)=(7/16,−1). To analyze it in the fourth quadrant of the parameter plane, we denote the eigenvalues of the Jacobian matrix at the origin E1as λ1, λ3<0< λ2, where λ3=ε3, λ2,1= ε2±√ε2 2+4ε1 2, 9 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 Fig. 9. For ε2= −1,B= −0.1,D=0.01: (a)–(d) partial bifurcation set in the second quadrant of the (ε1, ε3)-plane. Projection onto the (x,z)- plane of the T-point heteroclinic loop between E2and the equilibria E3,4(note that because of the symmetry a pair of the corresponding orbits exists): (e) principal T-point when (ε1, ε3)≈(−8.3738877,0.08418346); (f) secondary T-point when (ε1, ε3)≈(−8.4159326,0.08503685). (g) For (ε1, ε3)=(−8.2,0.084), geometric Lorenz attractor with initial conditions (x0,y0,z0)=(0.2,0,−8.2). On the curve H2 1there exists a degenerate point DH2 1, which occurs at (ε1, ε3)≈(−8.533961,0.085536), where the condition δE2=1 (resonant eigenvalues) is fulfilled (see Fig. 9(c)). Moreover, as we can see in Fig. 9(b)–(c), the curve H2 2ends at the point TP, which occurs at (ε1, ε3)≈(−8.37388775, 0.08418346), where there exists a (principal) heteroclinic T-point loop between E2and E3,4[19,59]. Its projection on the (x,z)-plane is drawn in Fig. 9(e). As can be seen in Fig. 9(b)–(d), three curves of global connections emerge from the point TP:He34 (heteroclinic connections between E3and E4), H3(homoclinic connections to E3and E4) and a spiral shaped curve H2 2of homoclinic connections to E2. Remark that the curves H3and H2 3end at a secondary T-point heteroclinic loop TPsbetween E2and E3,4, which occurs at (ε1, ε3)≈(−8.4159326,0.08503685), whose projection on the (x,z)-plane is drawn in Fig. 9(f). In this scenario of global bifurcations the presence of chaotic motions is guaranteed [16,17]. For instance, when (ε1, ε3)=(−8.2,0.084) (a point placed in the range of Fig. 9(c)), the attractor obtained with initial conditions (x0,y0,z0)=(0.2,0,−8.2) is represented in Fig. 9(g). Note that its shape is different from that of the other attractors drawn in Fig. 8(c)–(f). 16 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 Fig. 10. Partial bifurcation set in a neighborhood of the triple-zero bifurcation TZ, for B= −0.1 and D=0.01: (a) in the (ε1, ε2, ε3)-space; (b) projection onto the (ε2, ε1)-plane; (c) projection onto the (ε2, ε3)-plane. Seven curves of codimension-two bifurcations appear: three are local (BT, HZ,DZ) and four are global (DHe1,DHe2,DH1 1and DH2 1). 4.2.1. Triple-zero bifurcation Finally, to have a global idea of the bifurcation scenario and the main organizing centers in this region of the parameter space of system (3), we have drawn in Fig. 10 a partial bifurcation set in a neighborhood of the triple-zero degeneracy TZ, when B= −0.1 and D=0.01. The curves in the (ε1, ε2, ε3)-space appear in Fig. 10(a), their projection onto the (ε2, ε1)-plane in Fig. 10(b) and onto the (ε2, ε3)-plane in Fig. 10(c). Specifically, we draw the loci where the following codimension-two bifurcations occur: •BT, Bogdanov–Takens bifurcation of equilibrium E1. •HZ, Hopf-zero bifurcation of E1. •DZ, double-zero bifurcation (a double-zero eigenvalue with geometric multiplicity two) of E1. •DHe1and DHe2, degenerate heteroclinic connections because the equilibria E1and E2, respectively, change from saddle–node to saddle-focus. •DH1 1and DH2 1, degenerate homoclinic connections of E1and E2due to the existence of resonant eigenvalues δE1=1 and δE2=1, respectively. Finally we emphasize that similar configurations of curves of codimension-two global bifurcations around a triple-zero degeneracy TZ have been found in previous works, both in the study of some Z2-symmetric control systems particularized in the Chua’s equation [35] as well as in the analysis of non-symmetric electronic circuits [60]. 5. Conclusions In order to continue advancing in the analysis of the Lorenz system (and of so many other quadratic systems called Lorenz-like systems) it is necessary to combine analytical and numerical tools. Given the impossibility of studying by standard methods the triple-zero bifurcation that the equilibrium at the origin of the Lorenz system exhibits when σ= −1, ρ=1, b=0, in this work we carry out a partial study of system (3), an unfolding of the normal form of 17 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 the triple-zero bifurcation. We remark that several systems studied in the literature appear as particular cases for certain parameter choices [37–42]. The theoretical analysis of the double-zero bifurcation (a double-zero eigenvalue with geometric multiplicity two) has allowed to determine the different cases that may appear. For this we have established a connection between the reduced system on the corresponding two-dimensional center manifold and the planar normal form of the Hopf-zero bifurcation [44, Sect. 7.4]. For its greater interest, we have focused on the case that leads to a more complex dynamical behavior (existence of Hopf bifurcation and heteroclinic connection), and we have provided the expressions of the bifurcation curves organized by the double-zero degeneracy. Taking as a starting point the theoretical results obtained, we have carried out a numerical study that has allowed to find bifurcations of codimension one (Hopf, transcritical, pitchfork, saddle–node of periodic orbits, heteroclinic, homoclinic), two (double-zero, Bogdanov–Takens, degenerate Hopf, degenerate heteroclinic, degenerate homoclinic, T- point) and three (degenerate double-zero, triple-zero) in the region considered of the parameter space. On the other hand, we have also found zones of existence of chaotic attractors. Various aspects of the bifurcation set of system (3) are worth studying in the future. For example, to complete the study of degeneracies that appear on the bifurcation curves and to determine the new curves that emerge from those degenerate points. The double-zero bifurcation has allowed to find a heteroclinic connection. The degeneracies of this heteroclinic cycle are connected with the triple-zero bifurcation TZ (see Fig. 10). This codimension-three linear degeneracy deserves further study whose results can be applied to the Lorenz system. Declaration of competing interest The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. Acknowledgments We thank the reviewers for their careful reading of the manuscript and their very constructive remarks which have helped a lot to improve the presentation of the results. This work has been partially supported by the Ministerio de Economía y Competitividad, Spain (project MTM2017-87915-C2-1-P, co-financed with FEDER funds), by the Ministerio de Ciencia, Innovación y Universidades, Spain (project PGC2018-096265-B-I00, co-financed with FEDER funds) and by the Consejería de Economía, Innovación, Ciencia y Empleo de la Junta de Andalucía, Spain (FQM-276, TIC-0130, UHU-1260150 and P20_01160). Appendix A. Some useful changes of variables The first lemma is obtained directly from the results shown in [44, Sect. 7.4]. We state it here for clarity and completeness. Lemma A.1. System {˙ r=arz +cr3+drz2, ˙ z= −r2−z2+er2z+fz3,(A.1) for a = 0, by means of the change of variables s=r, w =z+hr2+iz2,dτ=(1 +jz)−1dt,(A.2) where h=c a,i=d+ae +2(a+1)c 3a,j=ae −2d+2(a+1)c 3a,(A.3) is transformed into ⎧ ⎪ ⎨ ⎪ ⎩ ds dτ=asw+O(|s, w|4), dw dτ= −s2−w2+ˆ fw3+O(|s, w|4), (A.4) where ˆ f=f−j. Lemma A.2. System {˙ r=arz +cr3z+drz3, ˙ z= −r2−z2+er4+fr2z2+mz4,(A.5) 18 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 for a = 0,−1/2, is orbitally equivalent to ⎧ ⎪ ⎨ ⎪ ⎩ ds dτ=asw+O(|s, w|5), dw dτ= −s2−w2+O(|s, w|5), (A.6) by means of the change of variables s=r(1 +gr2+hz2), w =z(1 +irz +jz2),dτ=(1 +kr2)−1dt,(A.7) where g= −a2e+(a+1)c+af −2am −d 2a(2a+1) ,h=d−am 2,j=m, i=a(2a+5)m+(1 −2a)d+ae −2af −c 2(2a+1) ,k=ae −am −c+d 2a,(A.8) Lemma A.3. System {˙ r=arz +cr5+dr3z2+erz4, ˙ z= −r2−z2+fr4z+gr2z3+hz5,(A.9) for a = 0, is orbitally equivalent to ⎧ ⎪ ⎨ ⎪ ⎩ ds dτ=asw+O(|s, w|6), dw dτ= −s2−w2+ˆ hw5+O(|s, w|6), (A.10) where ˆ h=h+2(1 −4a2)c−2(a+1)d+2e−3ag +a(1 −2a)f 5a. Appendix B. Hopf bifurcations of the nontrivial equilibria In this appendix we analyze the Hopf bifurcations of the nontrivial equilibria E3,4and E2of system (3). The results obtained, valid for any value of ε1,ε2and ε3, will be useful in the numerical study carried out in Section 4. B.1. Hopf bifurcation of the equilibria E3,4 Due to the symmetry, we will only consider in our study one of these equilibria, namely E3=(√−ε1(ε3+Dε1),0, ε1), that exists when ε1(ε3+Dε1)<0. Firstly, we translate the nontrivial equilibria E3to the origin by means of the change X=x−√−ε1(ε3+Dε1),Y=y,Z=z−ε1, that transforms system (3) into ⎧ ⎨ ⎩ ˙ X=Y, ˙ Y=ξ2Y−ξ1Z−XZ +BYZ, ˙ Z=2ξ1X+ξ3Z+X2+DZ2, (B.1) with ξ1=√−ε1(ε3+Dε1), ξ2=ε2+Bε1, ξ3=ε3+2Dε1.(B.2) Therefore, the characteristic polynomial of the linearization matrix of system (B.1) at the origin is given by p=λ3+p1λ2+p2λ+p3, where p1= −(ξ2+ξ3),p2=ξ2ξ3,p3=2ξ2 1. Thus, a Hopf bifurcation of E3occurs when p1p2=p3,p2>0, p1= 0, that is, when ξ2 1= −1 2ξ2ξ3(ξ2+ξ3), ξ2ξ3>0, ξ2+ξ3<0, ξ2<0, ξ3<0. 19 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 Considering system (B.1) at the critical values where the Hopf bifurcation occurs, we perform the following linear change to put the matrix of the linear part into the corresponding canonical form (X Y Z)=⎛ ⎝ 0 1 1 ω00−K2 1 K2K3K4⎞ ⎠(u v w), where ω0=√ξ2ξ3,K1=√−(ξ2+ξ3),K2=√2ξ2 K1 ,K3=√2ω0 K1 ,K4=√2ω0K1 ξ2 . Thus, with this change, system (B.1) becomes (˙ u ˙v ˙w)=⎛ ⎝ 0−ω00 ω00 0 0 0 −K2 1⎞ ⎠(u v w)(B.3) +1 C⎛ ⎝ A1u2+A2v2+A3w2+A4uv+A5uw+A6vw A7u2+A8v2+A9w2+A10uv+A11uw+A12vw −A7u2−A8v2−A9w2−A10uv−A11uw−A12vw⎞ ⎠, where C=ξ2K2+ξ3K2+K3ω0−K4ω0= 0, A1=K2(K3ω0−K4ω0−ξ2K2−ξ3K2), A2=K32E2+K32E3−K32+K3K4+E2+E3, A3= −K3K4E2+2K42E2−K3K4E3+2K42E3−K3K4+K42+E2+E3, A4=2K2K3E2+2K2K3E3−K32ω0+K3K4ω0−K2K3+K2K4, A5= −K2K3E2+3K2K4E2−K2K3E3+3K2K4E3−K3K4ω0 +K42ω0−K2K3+K2K4, A6= −K32E2+3K3K4E2−K32E3+3K3K4E3−K32+K42+2E2+2E3, A7=2ω0K22, A8=ω0K32+K3K2+ω0, A9=K4ξ2K2+K4ξ3K2+K42ω0+K4K2+ω0, A10 =K2(3K3ω0+K2), A11 =K2(ξ2K2+ξ3K2+3K4ω0+K2), A12 =K3ξ2K2+K3ξ3K2+2K4ω0K3+K3K2+K4K2+2ω0. Now, considering the second-order approximation to the center manifold w=α1u2+α2uv+α3v2+··· , we obtain the reduced system up to third order on the center manifold. And, by using the recursive algorithm developed in [46], we obtain the normal form for the Hopf bifurcation to third-order (22), where the first Lyapunov coefficient is given by a1=√ξ2ξ3N1 4ξ2ξ3(ξ2+ξ3)(ξ22+6ξ2ξ3+ξ32)(ξ22+3ξ2ξ3+ξ32),(B.4) with N1=2ξ25+31 ξ24ξ3+39 ξ23ξ32+ξ22ξ33−ξ2ξ34+4ξ24+56 ξ23ξ3 +25 ξ22ξ32−12 ξ2ξ33−ξ34+2ξ22ξ3−24 ξ2ξ32−2ξ33. Since ξ2ξ3>0, we obtain that a1=0 if, and only if, N1=0. In Fig. B.1 we have drawn the curve N1=0 in the (ξ2, ξ3)-plane near the origin. It is a bounded curve connecting the points (−2,0) and (0,−2). In terms of the original parameters of the system, the Hopf surface of the equilibria E3,4is given by ε3=2ε1−2(B+2D)ε1ε2−B(B+4D)ε2 1−ε2 2 2(Bε1+ε2) +√(Bε1+ε2)4+4ε1(ε1−2D(Bε1+ε2)ε1−(Bε1+ε2)2) 2(Bε1+ε2),(B.5) in the zones where the above inequalities are fulfilled. 20 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 Fig. B.1. Curve N1=0 in the (ξ2, ξ3)-plane in the vicinity of the origin. It corresponds to the degeneracy a1=0 in the Hopf bifurcation of the equilibria E3,4. Fig. B.2. For ε2= −1, B= −0.1, curve on the Hopf surface where M1=0 in the (D, ε1, ε3)-space. It corresponds to the degeneracy a1=0 in the Hopf bifurcation of equilibrium E2. B.2. Hopf bifurcation of the equilibrium E2 Now we consider the equilibrium E2=(0,0,−ε3/D). According to Section 2, it exhibits a Hopf bifurcation when Bε3=Dε2, if ε1+ε3/D<0 and ε3= 0. Repeating the previous calculations we obtain that the first Lyapunov coefficient, in terms of the original parameters of the system, is given by a1=B2M1 8√−Bε1+ε2 BDε2(4B2ε1+4Bε2−D2ε2 2),(B.6) with M1=8B2ε1+8Bε2−D2ε22−2Dε2. If B= 0, we obtain that a1=0 if, and only if, M1=0. For ε2= −1 and B= −0.1<0, in Fig. B.2 we have drawn, in the (D, ε1, ε3)-space, the curve on the Hopf surface where M1=0. This curve exists between the points (0,−10,0) (where a triple-zero bifurcation in a continuum of equilibria occurs) and (2,−10,20) (that corresponds to a degenerate Bogdanov–Takens bifurcation) [43]. Also, if we set a value D∈(0,2), there is a unique degeneracy a1=0 on the Hopf curve of E2when ε1<−10, and there are no degenerations if D<0. References [1] Lorenz EN. Deterministic non-periodic flows. J Atmospheric Sci 1963;20:130–41. http://dx.doi.org/10.1175/1520-0469(1963)020<0130:DNF>2. 0.CO;2. [2] Sparrow C. The Lorenz equation: bifurcations, chaos and strange attractors. New York: Springer; 1982, http://dx.doi.org/10.1007/978-1-4612- 5767-7. 21 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 [3] Haken H. Analogy between higher instabilities in fluids and lasers. Phys Lett A 1975;53:77–8. http://dx.doi.org/10.1016/0375-9601(75)90353-9. [4] Knobloch E. Chaos in the segmented disc dynamo. Phys Lett A 1981;82:439–40. http://dx.doi.org/10.1016/0375-9601(81)90274-7. [5] Gorman M, Widmann PJ, Robbins KA. Nonlinear dynamics of a convection loop: a quantitative comparison of experiment with theory. Physica D 1986;19:255–67. http://dx.doi.org/10.1016/0167-2789(86)90022-9. [6] Elgin JN, Molina-Garza JB. Traveling wave solutions of the Maxwell–Bloch equations. Phys Rev A 1987;35:3986–8. http://dx.doi.org/10.1103/ PhysRevA.35.3986. [7] Knobloch E, Proctor MRE, Weiss NO. Heteroclinic bifurcations in a simple model of double-diffusive convection. J Fluid Mech 1992;239:273–92. http://dx.doi.org/10.1017/S0022112092004403. [8] Cuomo KM, Oppenheim AV. Circuit implementation of synchronized chaos with applications to communications. Phys Rev Lett 1993;71:65–8. http://dx.doi.org/10.1103/PhysRevLett.71.65. [9] Poland D. Cooperative catalysis and chemical chaos: a chemical model for the Lorenz equations. Physica D 1993;65:86–99. http://dx.doi.org/ 10.1016/0167-2789(93)90006-M. [10] Hemati N. Strange attractors in brushless DC motors. IEEE T Circuits-I 1994;41:40–5. http://dx.doi.org/10.1109/81.260218. [11] Alexeev I. Lorenz system in the thermodynamic modelling of leukaemia malignancy. Med Hypotheses 2017;102:150–5. http://dx.doi.org/10. 1016/j.mehy.2017.03.027. [12] Tucker W. The Lorenz attractor exists. C R Acad Sci 1999;328:1197–202. http://dx.doi.org/10.1016/S0764-4442(99)80439-X. [13] Osinga HM, Krauskopf B. Visualizing the structure of chaos in the Lorenz system. Comput Graph 2002;26:815–23. http://dx.doi.org/10.1016/ S0097-8493(02)00136-X. [14] Barrio R, Serrano S. Bounds for the chaotic region in the Lorenz model. Physica D 2009;238:1615–24. http://dx.doi.org/10.1016/j.physd.2009. 04.019. [15] Yajima T, Nagahama H. Tangent bundle viewpoint of the Lorenz system and its chaotic behavior. Phys Lett A 2010;374:1315–9. http: //dx.doi.org/10.1016/j.physleta.2010.01.025. [16] Barrio R, Blesa F, Serrano S. Global organization of spiral structures in biparameter space of dissipative systems with Shilnikov saddle-foci. Phys Rev E 2011;84:035201. http://dx.doi.org/10.1103/PhysRevE.84.035201. [17] Barrio R, Shilnikov AL, Shilnikov LP. Kneadings, symbolic dynamics and painting Lorenz chaos. Int J Bifurcation Chaos 2012;22:1230016. http://dx.doi.org/10.1142/S0218127412300169. [18] Leonov GA, Kuznetsov NV, Korzhemanova NA, Kusakin DV. Lyapunov dimension formula for the global attractor of the Lorenz system. Commun Nonlinear Sci Numer Simul 2016;41:84–103. http://dx.doi.org/10.1016/j.cnsns.2016.04.032. [19] Glendinning P, Sparrow C. T-points: a codimension two heteroclinic bifurcation. J Stat Phys 1986;43:479–88. http://dx.doi.org/10.1007/ BF01020649. [20] Algaba A, Fernández-Sánchez F, Merino M, Rodríguez-Luis AJ. Analysis of the T-point-Hopf bifurcation in the Lorenz system. Commun Nonlinear Sci Numer Simul 2015;22:676–91. http://dx.doi.org/10.1016/j.cnsns.2014.09.025. [21] Creaser JL, Krauskopf B, Osinga HM. α-Flips and T-points in the Lorenz system. Nonlinearity 2015;28:R39–65. http://dx.doi.org/10.1088/0951- 7715/28/3/R39. [22] Doedel EJ, Krauskopf B, Osinga HM. Global bifurcations of the Lorenz manifold. Nonlinearity 2006;19:2947–72. http://dx.doi.org/10.1088/0951- 7715/19/12/013. [23] Doedel EJ, Krauskopf B, Osinga HM. Global invariant manifolds in the transition to preturbulence in the Lorenz system. Indag Math 2011;22:222–40. http://dx.doi.org/10.1016/j.indag.2011.10.007. [24] Doedel EJ, Krauskopf B, Osinga HM. Global organization of phase space in the transition to chaos in the Lorenz system. Nonlinearity 2015;28:R113–39. http://dx.doi.org/10.1088/0951-7715/28/11/R113. [25] Llibre J, Zhang X. Invariant algebraic surfaces of the Lorenz system. J Math Phys 2002;43:1622–45. http://dx.doi.org/10.1063/1.1435078. [26] Llibre J, Messias M, da Silva PR. Global dynamics of the Lorenz system with invariant algebraic surfaces. Int J Bifurcation Chaos 2010;20:3137–55. http://dx.doi.org/10.1142/S0218127410027593. [27] Algaba A, Gamero E, Merino M, Rodríguez-Luis AJ. Resonances of periodic orbits in the Lorenz system. Nonlinear Dyn 2016;84:2111–36. http://dx.doi.org/10.1007/s11071-016-2632-5. [28] Algaba A, Merino M, Rodríguez-Luis AJ. Superluminal periodic orbits in the Lorenz system. Commun Nonlinear Sci Numer Simul 2016;39:220–32. http://dx.doi.org/10.1016/j.cnsns.2016.03.004. [29] Algaba A, Fernández-Sánchez F, Merino M, Rodríguez-Luis AJ. Centers on center manifolds in the Lorenz, Chen and Lü systems. Commun Nonlinear Sci Numer Simul 2014;19:772–5. http://dx.doi.org/10.1016/j.cnsns.2013.08.003. [30] Algaba A, Domínguez-Moreno MC, Merino M, Rodríguez-Luis AJ. Study of the Hopf bifurcation in the Lorenz, Chen and Lü systems. Nonlinear Dyn 2015;79:885–902. http://dx.doi.org/10.1007/s11071-014-1709-2. [31] Algaba A, Domínguez-Moreno MC, Merino M, Rodríguez-Luis AJ. Takens–Bogdanov bifurcations of equilibria and periodic orbits in the Lorenz system. Commun Nonlinear Sci Numer Simul 2016;30:328–43. http://dx.doi.org/10.1016/j.cnsns.2015.06.034. [32] Algaba A, Domínguez-Moreno MC, Merino M, Rodríguez-Luis AJ. A review on some bifurcations in the Lorenz system. In: Carmona V, et al., editors. Nonlinear systems, vol. 1. Mathematical theory and computational methods. Cham: Springer; 2018, p. 3–36. http://dx.doi.org/10.1007/ 978-3-319-66766-9_1. [33] Gamero E. On the normal form of the triple-zero degeneracy with geometric multiplicity two. Dynam Cont Dis Ser A 2001;8:531–50. [34] Freire E, Gamero E, Rodríguez-Luis AJ, Algaba A. A note on the triple-zero linear degeneracy: Normal forms, dynamical and bifurcation behaviors of an unfolding. Int J Bifurcation Chaos 2002;12:2799–820. http://dx.doi.org/10.1142/S0218127402006175. [35] Algaba A, Merino M, Freire E, Gamero E, Rodríguez-Luis AJ. Some results on Chua’s equation near a triple-zero linear degeneracy. Int J Bifurcation Chaos 2003;13:583–608. http://dx.doi.org/10.1142/S0218127403006741. [36] Gamero E, Freire E, Rodríguez-Luis AJ, Ponce E, Algaba A. Hypernormal form calculation for triple-zero degeneracies. Bull Belg Math Soc Simon Stevin 1999;6:357–68. http://dx.doi.org/10.36045/bbms/1103065855. [37] Shimizu T, Morioka N. On the bifurcation of a symmetric limit cycle to an asymmetric one in a simple model. Phys Lett A 1980;76:201–4. http://dx.doi.org/10.1016/0375-9601(80)90466-1. [38] Shil’nikov AL. On bifurcations of the Lorenz attractor in the Shimizu-Morioka model. Physica D 1993;62:338–46. http://dx.doi.org/10.1016/0167- 2789(93)90292-9. [39] Rucklidge AM. Chaos in a low-order model of magnetoconvection. Physica D 1993;62:323–37. http://dx.doi.org/10.1016/0167-2789(93)90291-8. [40] Liu C, Liu T, Liu L, Liu K. A new chaotic attractor. Chaos Solitons Fractals 2004;22:1031–8. http://dx.doi.org/10.1016/j.chaos.2004.02.060. [41] Mello LF, Messias M, Braga DC. Bifurcation analysis of a new Lorenz-like chaotic system. Chaos Solitons Fractals 2008;37:1224–55. http: //dx.doi.org/10.1016/j.chaos.2007.11.008. [42] Kokubu H, Roussarie R. Existence of a singulary degenerate heteroclinic cycle in the Lorenz system and its dynamical consequences: Part I. J Dyn Differ Equ 2004;16:513–57. http://dx.doi.org/10.1007/s10884-004-4290-4. 22 A. Algaba, M.C. Domínguez-Moreno, M. Merino et al. Communications in Nonlinear Science and Numerical Simulation 111 (2022) 106482 [43] Algaba A, Domínguez-Moreno MC, Merino M, Rodríguez-Luis AJ. A degenerate Takens–Bogdanov bifurcation in a normal form of Lorenz’s equations. In: Lacarbonara W, et al., editors. Advances in Nonlinear Dynamics. Proceedings of the Second International Nonlinear Dynamics Conference (NODYCON 2021), vol. 1. Cham: Springer; 2022, p. 699–709. http://dx.doi.org/10.1007/978-3-030-81162-4_60. [44] Guckenheimer J, Holmes PJ. Nonlinear oscillations, dynamical systems, and bifurcations of vector fields. New York: Springer; 1983, http: //dx.doi.org/10.1007/978-1-4612-1140-2. [45] Wiggins S. Introduction to applied dynamical systems and chaos. New York: Springer; 2003, http://dx.doi.org/10.1007/b97481. [46] Gamero E, Freire E, Ponce E. Normal forms for planar systems with nilpotent linear part. In: Seydel R, et al., editors. Bifurcation and chaos: analysis, algorithms, applications. International Series of Numerical Mathematics, vol. 97, Basel: Birkhäuser; 1991, p. 123–7. http: //dx.doi.org/10.1007/978-3-0348-7004-7_14. [47] Keener JP. Infinite period bifurcation and global bifurcation branches. SIAM J Appl Math 1981;41:127–44. http://dx.doi.org/10.1137/0141010. [48] Chow SN, Li C, Wang D. Normal forms and bifurcation of planar vector fields. Cambridge: Cambridge University Press; 1994, http://dx.doi.org/ 10.1017/CBO9780511665639. [49] Qin BW, Chung KW, Algaba A, Rodríguez-Luis AJ. High-order approximation of heteroclinic bifurcations in truncated 2D-normal forms for the generic cases of Hopf-zero and non-resonant double Hopf singularities. SIAM J Appl Dynam Syst 2021;20:403–37. http://dx.doi.org/10.1137/ 0141010. [50] Algaba A, Chung KW, Qin BW, Rodríguez-Luis AJ. A nonlinear time transformation method to compute all the coefficients for the homoclinic bifurcation in the quadratic Takens–Bogdanov normal form. Nonlinear Dyn 2019;97:979–90. http://dx.doi.org/10.1007/s11071-019-05025-2. [51] Algaba A, Chung KW, Qin BW, Rodríguez-Luis AJ. Computation of all the coefficients for the global connections in the Z2-symmetric Takens–Bogdanov normal forms. Commun Nonlinear Sci Numer Simul 2020;81:105012. http://dx.doi.org/10.1016/j.cnsns.2019.105012. [52] Qin BW, Chung KW, Algaba A, Rodríguez-Luis AJ. High-order analysis of global bifurcations in a codimension-three Takens–Bogdanov singularity in reversible systems. Int J Bifurcation Chaos 2020;30:2050017. http://dx.doi.org/10.1142/S0218127420500170. [53] Qin BW, Chung KW, Algaba A, Rodríguez-Luis AJ. Analytical approximation of cuspidal loops using a nonlinear time transformation method. Appl Math Comput 2020;373:125042. http://dx.doi.org/10.1016/j.amc.2020.125042. [54] Qin BW, Chung KW. Personal communication (2019). [55] Doedel EJ, Champneys AR, Dercole F, Fairgrieve T, Kuznetsov Y, Oldeman BE, Paffenroth R, Sandstede B, Wang X, Zhang C. Auto07-P: Continuation and bifurcation software for ordinary differential equations (with HomCont). Technical report, Concordia University; 2012. [56] Krupa M, Melbourne I. Asymptotic stability of heteroclinic cycles in systems with symmetry. Ergod Theor Dyn Syst 1995;15:121–47. http://dx.doi.org/10.1017/S0143385700008270. [57] Krupa M, Melbourne I. Asymptotic stability of heteroclinic cycles in systems with symmetry, II. P Roy Soc Edinb A 2004;134:1177–97. http://dx.doi.org/10.1017/S0308210500003693. [58] Geng F, Xu Y. Bifurcations of heteroclinic loop accompanied by pitchfork bifurcation. Nonlinear Dyn 2012;70:1645–55. http://dx.doi.org/10. 1007/s11071-012-0563-3. [59] Fernández-Sánchez F, Freire E, Rodríguez-Luis AJ. T-points in a Z2-symmetric electronic oscillator. (I) Analysis. Nonlinear Dynam 2002;28:53–69. http://dx.doi.org/10.1023/A:1014917324652. [60] Algaba A, Merino M, García C, Reyes M. Degenerate global bifurcations in a simple circuit. Int J Pure Appl Math 2009;57:265–78. 23