A New Approach of Asymmetric Homoclinic and Heteroclinic Orbits Construction in Several Typical Systems Based on the Undetermined Padé Approximation Method

In dynamic systems, some nonlinearities generate special connection problems of non-Z 2 symmetric homoclinic and heteroclinic orbits. Such orbits are important for analyzing problems of global bifurcation and chaos. In this paper, a general analytical method, based on the undetermined Padé approximation method, is proposed to construct non-Z 2 symmetric homoclinic and heteroclinic orbits which are affected by nonlinearity factors. Geometric and symmetrical characteristics of non-Z 2 heteroclinic orbits are analyzed in detail. An undetermined frequency coefficient and a corresponding new analytic expression are introduced to improve the accuracy of the orbit trajectory. The proposed method shows high precision results for the Nagumo system (one single orbit); general types of non-Z 2 symmetric nonlinear quintic systems (orbit with one cusp); and Z 2 symmetric system with high-order nonlinear terms (orbit with two cusps). Finally, numerical simulations are used to verify the techniques and demonstrate the enhanced efficiency and precision of the proposed method.


Introduction
Heteroclinic orbits (HOs) play a central role in nonlinear dynamics and many scholars have undertaken the study of HOs and heteroclinic bifurcation [1].In mathematics, a HO (or heteroclinic connection) is a path in phase space which joins two different equilibrium points ( 0 and  1 ) of a dynamical system.A HO is contained in the stable manifold of  1 (   1 ( 1 )) and the unstable manifold of  0 (   0 ( 0 )) [2] as shown in Figure 1.There are many physical applications which correspond to the homoclinic and heteroclinic orbits [3][4][5], such as the pull-in phenomenon in MEMS [6] and the self-contact problem of elastic rods [7].
Existing research of HOs is focused mainly on classic loworder Z 2 symmetric system.The study of Z 2 symmetric HOs in an unperturbed system has resulted in the development of a number of analytical methods.Some examples include the generalized harmonic multiple scales method [8], the generalized hyperbolic perturbation method [9], and the elliptic Lindstedt-Poincaré (LP) method [10].HOs in classic loworder autonomous systems, which consider perturbations of a Z 2 symmetric system, are also well understood.For example, analysis of weak nonlinear systems employs techniques such as the elliptic perturbation method [11], the hyperbolic perturbation method [12], and the method developed by Izydorek and Janczewska [13].Similarly, in a strongly nonlinear system, the undetermined fundamental frequency method [14], the modification of the perturbationincremental method [15], the hyperbolic LP method [16], and the generalized harmonic function perturbation method [17] are all applicable analysis techniques.
In the aforementioned techniques, the HOs all meet the characteristics that the saddle points are symmetric about a center point.The more general asymmetric case remains to be further explored where the non-Z 2 symmetric HOs are neither symmetric about x-(or y-) axis nor symmetric about the origin.Special non-Z 2 symmetric HOs with some complicated properties (orbit with cusps) were discovered when 2 Mathematical Problems in Engineering concerning the number of limit cycles and their relative locations in polynomial vector fields.The existence of these orbits has also been proved recently [18][19][20].In fact, regardless of the number of cusps contained in the orbits, their properties can be discussed according to the symmetry along the -axis and -axis.In engineering applications, such as power systems, special non-Z 2 symmetric HOs have been found.Beyn [21] applied numerical methods to analyze the truncation error in bifurcation calculations and the non-Z 2 symmetric HOs, which rose the phenomenon of superconvergence, were discovered.Power system analysis by Pai and Laufenberg [22] revealed that local instabilities are due to the influence of non-Z 2 symmetric heteroclinic orbits, which means that the characteristics of the system will be a qualitative change.Petris ¸or [23] also found that non-Z 2 symmetric heteroclinic saddle points are present in the study of the dynamic characteristic of reversible magnetic vector fields.Since such orbits are more complicated than HOs with non-Z 2 symmetry in classic autonomous low-order systems, these orbits are still being solved by numerical methods [24][25][26].To the authors' knowledge, existing analytical methods are rarely used for the analysis of these orbits due to the complexity of their application.Therefore, a more effective analytical method is necessary.A new analytical method based on the undetermined Padé approximation method is proposed in this paper.The application of this new method results in a higher accuracy of the achieved HOs.This method not only is applicable to totally asymmetric HOs but also shows its superiority in symmetric and conventional asymmetric HOs.
This paper is organized as follows: in Section 2, the symmetry characteristics of non-Z 2 symmetric HOs are categorized.After that, the proposed method, the "undetermined Padé approximant method, " is introduced.In Section 3, analysis of HOs in three general systems (categorized by the complexity of the symmetry) demonstrates the validity and feasibility of the proposed method, which is verified against numerical simulations.Finally, the conclusions are presented in Section 4.

The Undetermined Padé Approximation Method
2.1.Symmetry Analysis of Non-Z 2 Symmetry HOs.Consider the following general non-Z 2 symmetric nonlinear dynamical system: Suppose that system (1) has a center (0, 0) and two saddle points  1 (  , 0) and  2 (  , 0) which are located on both sides of the center and  represents the perturbation parameter.The corresponding HOs of system (1) are shown in Figure 2. When  = 0, the unperturbed system is conservative and the orbit trajectories (the dotted line in Figure 2) have extreme values in the -axis direction corresponding to the saddle points: (2) Figure 1: Heteroclinic connection is a path joining two different equilibrium points.
The extreme values of the trajectories in the vertical direction are In other words, the maximum value of the orbit in the ẋ -axis direction lies in the center .
When a small perturbation is applied ( ̸ = 0), HOs still exist in this autonomous system.There will be a critical orbit close to the saddle points, as shown by the solid line in Figure 2. In this case, the maximum value of the orbit in the ẋ -axis direction does not lie on the center ; that is, The non-Z 2 symmetric HOs of simple autonomous system have been widely studied, as shown by the solid line in Figure 2(a).The unperturbed system has Z 2 symmetric HOs, shown by the dotted line in Figure 2(a), where the saddle points fulfill the criteria In other words, the distances between each saddle point and the center are equal.However, in some complex systems, the unperturbed system will also result in a non-Z 2 symmetric HO phenomenon where the saddle points satisfy the following condition: In other words, the distances between each saddle point and the center are not equal.As shown in Figure 2(b), the HOs of the unperturbed system (dotted lines) and the HOs of the autonomous system (solid lines) satisfy condition (6); that is, all HOs are non-Z 2 symmetric.The solid line trajectories shown in Figure 2(b) correspond to those totally asymmetric HOs discussed previously in Section 1.
A trajectory that does not possess Z 2 symmetry characteristics is classified as non-Z 2 symmetry.This means that the saddle points are asymmetric,   ̸ =   , or the maximum points of the trajectory are not on the ẋ -axis, or both cases happened simultaneously.There is no guarantee that analytical methods which are suitable for solving Z 2 symmetric HOs will return the required level of accuracy if applied to Conservative system Autonomous system Conservative system Autonomous system 1.0 0.5 0.0 a non-Z 2 symmetric HO.For example, when the saddle point   does not coincide with −  , the calculations will show that the obtained results cannot pass through both saddle points.If the convergence conditions are that the heteroclinic solution finally tends to one saddle point   (−  ), then the determined trajectory will shift to this saddle point   (−  ).Thus, complete orbits cannot go through another saddle point −  (  ).As shown in Figure 3, (1) represents that the analytical results calculated by employing the method discussed in [27] and (2) are numerical results.
In order to resolve this problem, the Padé approximation method [28,29] is used.The extended Padé approximation [27,[30][31][32] has proven its validity and reliability through its simple arithmetic process and accurate results in solving the non-Z 2 symmetric homoclinic orbit.The solving process of this superior method can be further enhanced and is presented in detail in the next section.

The Calculation Procedure.
In this section, the solving process of this analytical method is proposed.The differential equation, including the entire nonlinear part, is solved directly.In addition, the magnitude of the perturbation parameter  does not need to be ascertained in advance.It has universality and can directly solve any symmetric case of critical orbits.
A generic case, such as the solid line orbits shown in Figure 2(b), is taken for example.This trajectory, which is completely nonsymmetric, simultaneously satisfies conditions (4) and (6).Therefore, the solution achieved in this example is universal and versatile.
First, without loss of generality, system (1) is rewritten as the following general form (considering a one-degree-offreedom general nonlinear system): where (, /, ) is a nonlinear term including damping and  0 is the natural frequency of the system.It is obvious that system ( 7) is an autonomous system ( ̸ = 0) and we can assume that system (7) has HOs which are non-Z 2 symmetric.
Take a series solution of system (7) as with ( 0 ,  1 ) being the initial value of the trajectory point such that where the undetermined parameters  0 and  1 are arbitrary nonzero constants.The rest of the parameters of ( 8) can be expressed as a function of  0 and  1 as follows: Usually, the initial point is taken as one of the maximum points on the trajectory where it also fulfills its tangent perpendicular to the ẋ -axis; namely, If system ( 7) is a conservative system ( = 0), then the maximum point is on the ẋ -axis.In this case, condition (11) above can be simplified to According to the abovementioned initial conditions, we consider a new heteroclinic solution expression given by An undetermined constant  is introduced into expression (13).Compared with the traditional method, the undetermined constant  not only reflects the original frequency  0 but also contains the effect of the nonlinear term (, /, ).Expression (13) satisfied the following properties: The form of the expression ensures its convergence.Meanwhile, the convergence speed of the whole calculation process is greatly accelerated, reducing the amount of computation time required.This expression can resolve the orbits in both autonomous and conservative systems with the boundary conditions  0 and   /  .
When  tends to infinity, the heteroclinic solutions tend to   and   : For the upper orbit, the boundary conditions can be expressed as left boundary  0 →   and right boundary   /   →   .For the lower orbit, left boundary   /  →   ; right boundary  0 →   .The abovementioned relationships are illustrated in Figure 4.By equating the heteroclinic solution QPA  to series solution (8), A functional relationship between the coefficients   and   (heteroclinic solution QPA  ) and coefficient   (series solution) is achieved by comparing terms of ( 16) to the same order of .This gives the following: Equations (17) show that all unknown coefficients   and   of QPA  can be expressed as a function of  0 and  1 .By applying boundary conditions ( 14) and ( 15), the initial values  0 and  1 and the undetermined constant  can be calculated.When the determined values of  0 and  1 are substituted into (17), every coefficient of QPA  can be determined and even the higher order expressions are achievable.For instance, the specific expression of the fourthorder heteroclinic solution is given as follows: The higher the approximation order, the better the accuracy.Ordinarily, a calculation to second-order approximation (QPA 2 ) would meet the accuracy requirements (error range is Δ < 1%) for low-order Z 2 symmetric systems (only with cubic nonlinearity).Similarly, a calculation to fourthorder approximation (QPA 4 ) would meet very high accuracy, generally required for systems with non-Z 2 symmetric terms (i.e., with square or fourth nonlinearity) or systems with highorder nonlinear terms (i.e., with fifth nonlinearity and so on).

Some Examples
In this section, the analytical approach proposed is applied to several classical nonlinear systems.

The Nagumo System.
The HO in this system meets the characteristics that the two distances between each saddle point and the center are not equal (condition ( 6)), and the maximum value of the orbit in the ẋ -axis direction dose not lie in the center (condition ( 4)).
The Nagumo system [22] is defined as where  1 = √ 2(1/2 −  2 ) and 0 <  2 < 1/2.For all values of  2 , system (19) has one center ( 2 , 0) and a HO connecting the two saddle points ( 1 (0, 0) and  2 (1, 0)).The phase diagram is shown in Figure 5 with  2 = 2/5.The maximum point of the trajectory of autonomous system (19) is no longer at the center  and the saddle points are asymmetric ( 1 >  2 ).Using the present method, coefficient expressions (10) can be written as the following series of coefficients: Using condition (11), we find that The fourth-order heteroclinic solution can be achieved making use of (17).By applying boundary conditions ( 14) and ( 15), the initial value  0 = 0.5 and undetermined constant  = 0.49999 are calculated.The analytical heteroclinic solution (A.1) is also obtained.The superscript "+" indicates the solution with an initial value  1 > 0, which is the upper half of the trajectory.Similarly, the superscript "−" indicates the solution with an initial value  1 < 0 corresponding to the lower half of the trajectory.The corresponding phase diagram is shown in Figure 6 with  0 = 0.5 and  1 = 0.176776.The error of the maximum value of the orbit in the ẋ -axis direction between analytical results and numerical results is 0.00057%.

General Non-Z 2 Symmetric Nonlinear Quintic System.
The HOs (with one cusp) in this system meet the characteristics that the two distances between each saddle point and the center are not equal (condition ( 6)).When  = 0 the maximum value of the orbit in the ẋ -axis direction lies in the center (condition (3)).If  ̸ = 0, the maximum value of the orbit in the ẋ -axis direction dose not lie in the center (condition (4)).Very few analytical methods exist for solving the problem with a fifth-order nonlinearity.Even fewer methods are available to analyze the problem of complex  0.0 0.5 non-Z 2 symmetry systems.However, the undetermined Padé approximation can be employed to solve such problems with high efficiency and precision.Consider the following autonomous system containing non-Z 2 symmetric term: which has a rich physical meaning and is widely used in many fields.Equation (22) contains nonlinear terms from secondary to quintic, which is a typical non-Z 2 symmetric system.System (22) has five equilibrium points, so the corresponding stable solution and its phase diagram have a variety of possible cases and complex dynamic behavior.Therefore, system (22) has a good generality, which encompasses the systems discussed in the literature [9,16,17,30].Here, we study the cases in which all nonlinear terms of ( 22) exist.In this scenario, its system has a pair of non-Z 2 symmetric homoclinic orbits and a pair of non-Z 2 symmetric heteroclinic orbits.Those heteroclinic orbits have asymmetric saddle points and a maximum point which intersects the ẋaxis.
Firstly, we analyze unperturbed system (23) of system (22) where  = 0: The case when conservative system (23) has five equilibrium points is studied, and any singular point can be selected as the zero point.For the parameters used in Figures 7 and 8, system (23) has three centers and two saddle points, a pair of non-Z 2 symmetric homoclinic orbits, and a pair of non-Z 2 symmetric HOs (with one cusp).There are two kinds of orbits when the parameters   take different values: one is where   > |  | shown in Figure 7 (cusp towards right) and the other one is   < |  | shown in Figure 8 (cusp towards left).For brevity, the calculation process of Figures 7 and 8 is not presented here in this paper.By applying this method, the initial values  1 = ±0.904032and  = 1.20207 are obtained.The analytical heteroclinic solutions of system (23) are formulas (A.2) and (A.3).Similarly, analytical expressions of a pair of non-Z 2 symmetric homoclinic orbits can also be obtained by using this method.The analytical results for the upper HO in Figure 7 return the parameter values  0 = 0 and  1 = 0.904032 and the numerical results predict  0 = 0 and  1 = 0.900371.The error of the maximum value of the orbit in the ẋ -axis direction between the analytical results and the numerical results is 0.4%.Results for the lower orbit are similar to the upper orbit with an error of 0.4%.
In the perturbed system ( ̸ = 0), trajectories in autonomous system (22) meet characteristics (4) and (6).It should be noted here that the perturbation parameter  cannot be too large; otherwise, the phase diagram of system (22) will be significantly affected.When  1 = 0,  2 = −1, and  = 0.01, the value of heteroclinic bifurcation parameter is  = −0.155.By applying the analytical method developed in this paper, we obtain the analytic heteroclinic solutions (A.4) and (A.5) and the corresponding phase diagram of system (22) is shown in Figure 9.The analytical results for the upper orbit are  0 = −0.000588768and  1 = 0.908729 with the corresponding numerical results giving  0 = −0.000596024and  1 = 0.900226.The error of the maximum value of the orbit in the ẋ -axis direction between the analytical results and the numerical results is 0.9%.Similarly, the analytical results for the lower orbit are  0 = 0.000585934 and  1 = −0.904105with the numerical results returning  0 = 0.000590165 and  1 = −0.900514.The error of the maximum value of the orbit in the ẋ -axis direction between the analytical results and the numerical results is 0.4%.

General Z 2 Symmetric System with High-Order Nonlinear
Terms.The HOs (with two cusps) in this system meet the characteristics that the distances between each saddle point and the center are equal (condition (5)).If  = 0, the maximum value of the orbit in the ẋ -axis direction lies in the center (condition (3)).If  ̸ = 0, the maximum value of the orbit in the ẋ -axis direction does not lie in the center (condition (4)).
Consider the following triple-well Z 2 symmetric autonomous system: ẍ + ( 1  +  3  3 +  5  5 ) =  ( +  1  −  2  2 ) ẋ . ( When  = 0, this conservative system is Z 2 symmetric and has Z 2 symmetric HOs (with two cusps).As shown in Figure 10, when  1 = 1,  3 = −2, and  5 = 1, the saddle points of this conservative system are ±1.The application of this method is also available and precise and the analytical heteroclinic solutions are (A.6) and (A.7).The analytical results for the upper orbit are  0 = 0 and  1 = 0.586815 and the numerical results are  0 = 0 and  1 = 0.57735.The error of the maximum value of the orbit in the ẋ -axis direction between the analytical results and the numerical results is 1.6%.
When  ̸ = 0, the autonomous system has non-Z 2 symmetric HOs (with two cusps).Bifurcation curves, as shown in Figure 11, can be obtained by applying this method.When  1 = 0,  2 = −1, and  = 1, the value of the bifurcation parameter obtained by this method is  = −0.172755and

Figure 4 :
Figure 4: The boundary condition for (a) the upper orbit and (b) the lower orbit.

Figure 6 :
Figure 6: Comparison of HO calculated using the proposed analytical method (∘) and a numerical Runge-Kutta algorithm (solid).

Figure 9 :
Figure 9: Comparison of HOs for the present results and the Runge-Kutta procedure.

Figure 10 : 1 Figure 11 :
Figure 10: HOs for the present results and the Runge-Kutta procedure comparison.