DIRECT APPROACH TO DETECT THE HETEROCLINIC BIFURCATION OF THE PLANAR NONLINEAR SYSTEM

. In this paper, we present a novel way of directly detecting the heteroclinic bifurcation of nonlinear systems without iteration or Melnikov type integration. The method regards the phase and fundamental frequency in a hyperbolic function solution and bifurcation parameter as the unknown com- ponents. A global collocation point, obtained from the energy balance method, together with two special points on the orbit are used to determine these un- known components. The feasibility analysis is presented to have a clear insight into the method. As an example, in a third-order nonlinear system, an expres- sion for the orbit and the critical value of bifurcation are directly obtained, maintaining the precision but reducing the complication of bifurcation analy- sis. A second-order collocation point improves the accuracy of computation. For a broader application, the eﬀectiveness of this new approach is veriﬁed for systems with a large perturbation parameter and the homoclinic bifurcation problem evolving from the even order nonlinearity.


1.
Introduction. Homoclinic and heteroclinic (global) orbits are particular solutions for the governing equations connected, in the limit of infinite time, with a given unstable equilibrium point. Although they are distinct and isolated solutions, they play a fundamental role in the organization of the global dynamics of complex systems, as repeatedly verified in practice [1,20]. However, their detection is not easy, and the main difficulties lie in finding the critical value of bifurcation and analytic solutions for these manifolds. For a long period of time, the Melnikov method [19] was commonly applied to detect the critical value of global bifurcation. Following that, Belhaq and his coworkers [2,4,5] presented other methods to detect the critical value, which were based on the collision of a periodic orbit with a saddle point and letting the period of oscillator become infinite. Zhang et al. [22,23] applied the method of undetermined fundamental frequency to find the bifurcation values of homoclinic and heteroclinic strongly nonlinear systems. In general, most subsequent research shares a common approach with the Melnikov method for finding the critical value of bifurcation.
Moreover, the biggest challenge comes in obtaining analytic expressions for the global manifolds (the expressions for the homoclinic or heteroclinic orbits). Vakakis and Azeez [21] developed an iterative technique to approximate certain homoclinic orbits. Mikhlin [18] and Feng et al. [13] used the Padé and quasi-Padé approximation to construct these types of orbits. Other methods, such as the Homotopy perturbation [14], perturbation-incremental [6], hyperbolic perturbation [7,9], elliptic Lindstedt-Poincaré (LP) [3] and hyperbolic LP [8] methods and so on, have also been used to determine the homoclinic and heteroclinic solutions of both weakly and strongly nonlinear oscillators. Chen et al. [9] used the hyperbolic perturbation method, involving a new time scale in the hyperbolic function, to derive the expressions of homoclinic orbits, where the systems were characterized by quadratic, cubic, quartic, and strong nonlinearity.
The object of this paper is to find the critical value and the orbit expression for heteroclinic type orbits through a direct approach. By 'direct', we mean that the critical value and the expression for a heteroclinic orbit can be obtained in a single step, where there are no successive iterations or Melnikov type integration, which distinguishes it from other approaches. We start by introducing the phase, frequency and bifurcation parameter as unknown perturbation components to construct the heteroclinic orbit expression. Consequently, three special points (including a new global collocation point) on the orbit are calculated to find these unknown parameters. The new collocation point is located by globally extending the energy balance method (EBM) [15]. Then, using these three points, we obtain an analytic expression and the critical bifurcation value of the heteroclinic orbit. The numerical simulations verify the effectiveness of this method.
By using the kinetic and potential energy balance equation and finding a collocation critical point at π/4, the traditional EBM was introduced by He [15] and improved by other scientists [10,15,16,17,24]. They used this method to investigate periodic solutions of the general weakly nonlinear systems as well as strongly nonlinear systems (where the perturbation parameter in the equation does not have a small value). By combining the improved amplitude-frequency formulation and EBM, Davodi et al. [10] increased its accuracy, reliability and simplicity for solving nonlinear oscillating systems. Using the collocation method and the Galerkin-Petrov method, Khan and Mirzabeigy [16] improved the accuracy of the EBM for even larger amplitudes of oscillation. Younesian et al. [24] considered all odd-type forcing functions involved in the governing equation, where the natural frequency was obtained as a function of the initial conditions. To expand its engineering applications, Mayoof and Hawwa [17] used EBM in a physical model using the Duffing equation for the vibrational analysis of a curved single-walled carbon nanotube.
Differing from the general periodic analysis [10,15,16,17,24], and using a hyperbolic function expression, a global extension of EBM based on a heteroclinic collocation point is obtained. Furthermore, in the framework of our hyperbolic function expression, it is effective to increase the accuracy of the global EBM by using a higher-order collocation point.
The article is structured as follows. In section 2, some basic ideas of the method, such as the perturbation of the hyperbolic function with phase θ and He's EBM are introduced. Then, based on the hyperbolic form of the balance of kinetic and potential energy, a new global collocation point is computed. This forms an auxiliary equation in the following bifurcation analysis. A feasibility analysis is presented to give a clear insight into the new approach and answer the question why it works. In section 3, for the heteroclinic bifurcations of a third-order nonlinear system, all the undetermined parameters and the critical value of bifurcation are obtained in a single equation group. This maintains precision and reduces the complication of the bifurcation analysis. In section 4, we use the second-order collocation point to improve the accuracy of the direct approach and also to verify that the available range covers typical strongly nonlinear systems. Finally, the methodology in homoclinic bifurcation with even order nonlinearity is discussed to present the possible applications in broader areas.
2. Basic idea of the perturbed hyperbolic orbit and the global extension of EBM. Consider the following nonlinear oscillator governed by Eq. (1) For the different choices of c and components in f 1 (u), Eq. (1) may have homoclinic or heteroclinic orbit connections with the saddle points. For example, for c = 1, f 1 (u) = u 2 , the unperturbed system (ε = 0) processes a homoclinic orbit to a saddle point S(−1, 0) as shown in Figure 1(a); then, in c = 1, f 1 (u) = −u 3 , the unperturbed system processes a heteroclinic orbit through two saddle points S 1,2 (±1, 0) as shown in Figure 1(b). f 2 (µ, u,u) represents an arbitrary nonlinear function of its arguments, and µ refers to the bifurcation parameter. In this article, we mainly discuss the heteroclinic bifurcation, and homoclinic bifurcation can be detected similarly. With the phase portraits of the perturbed (ε = 0) and unperturbed (ε = 0) heteroclinic orbit as shown in Figure 2, we realize that the existence of perturbation mainly produces a phase hysteresis between u andu. It can be seen more clearly on the plot ofu shown in Figure 2 (c). As a result, using the basic hyperbolic functions solutions [7] u = ±a tanh ω 10 t, u = ±aω 10 sech 2 (ω 10 t),

LING-HAO ZHANG AND WEI WANG
we give the general heteroclinic solution with phase perturbation whereu represents the time derivative of the solution u, a is the amplitude, and ω 10 is the angular frequency; the phase θ represents the effect of perturbation in Eq. (1) assumed as a static or dynamic variable of time t. The feasibility analysis of that hypothesis is presented in section 2.3.
We regard a as a constant number that means the positions of saddle points do not change with perturbation. Then ω 10 , θ and the bifurcation parameter µ can be seen as the undetermined variables, to be decided by considering the special points on the orbit (see Figure 2).  [15]. Consider the nonlinear oscillator governed bÿ The Hamiltonian, therefore, can be written in the form We assume ω 10 as the unknown angular frequency, and substitute the trial function u = a cos ω 10 t, into Eq. (5). Assuming that solution as the exact solution of the system, according to He and other researchers, the Hamiltonian should be valid for all values of t. However, as Eq. (6) is only a type of approximation solution, the validation of Eq. (5) in terms of that trial function can be held only at some special points, the collocation point [15]. The whole energy always changes between kinetic energy and potential energy. Following the balance of kinetic and potential energy, that point is located at ω 10 t = π/4, which gives the frequency of oscillation as That result maintains a similar accuracy as compared to other perturbation methods.
To finish the global extension of EBM, begin with the heteroclinic bifurcation of Eq. (4). Its Hamiltonian comes to According to EBM, given a hyperbolic solution can be considered as the exact solution of the system, the Hamiltonian in Eq. (8) would be available for all values of t by appropriate choice of ω 10 , the unknown fundamental frequency. However, the unperturbed hyperbolic solution in Eq. (2) is only a type of approximation, which means Eq. (8) is only valid at special points. From the balance of kinetic and potential energy, we obtain According to Chen et al. [9], in a preliminary computation, ω 10 equals ω 0 √ 2 in the case of zero perturbation (ε = 0). That produces the global collocation point for heteroclinic bifurcation Based on the global collocation point, the Hamiltonian in Eq. (8) becomes a third supplementary equation. In the following sections, we will show how to use these special points to find the undetermined parameters and the critical value of bifurcation.

The feasibility analysis of the new approach.
Here the problem is about the heteroclinic bifurcation with Z 2 symmetry. So we consider the Duffing type nonlinearity in the conserve part and other cases can be discussed similarlÿ Multipleu on both sides of Eq. (11) and finish the integration where H is the Humiliation of the system. Consider a third order expression of f 2 (µ, u,u) so that the integration in terms of Eq. (2) is a 2 (a 2 γ 2,1 − 5µ)ω 10 tanh ω 10 t − 1 15 a 2 (a 2 γ 2,1 − 5µ)ω 10 sech 2 ω 10 t tanh ω 10 t + 1 5 a 4 γ 2,1 ω 10 sech 4 ω 10 t tanh ω 10 t.
Eq. (23) verifies that the dynamic phase θ may cover more complicated perturbation component in the governing equation so that Eq. (3) can be regarded as the general solution in other systems. In section 5, we use the new methodology to address the homoclinic bifurcation with even order nonlinearity in the conserve part.
From the discussions above, it is possible to concern a hyperbolic function solution to the heteroclinic bifurcation problem, with frequency ω 10 , phase θ and bifurcation parameter µ as the variables. Basically, in order to find those unknown variables, 3 independent equations are necessary. They come from substituting the orbit expression into the original system and Hamiltonian equation respectively, at 3 referent points (u(φ i ),u(φ i )) with phases φ i , i = 1, 2, 3. So the choice of those points become extremely important and decide the accuracy of heteroclinic analysis.
To maintain the accuracy of computation, those points should be chosen as the orbit passing through them very closely. As seen in Figure 2 (a) and (b), the system is Z 2 symmetry that means (u,u) with φ 1 = 0 leads to the first referent point. Then, from Figure 2 (c),u reaches its maximum value with a phase hysteresis between u andu that represents the second referent point φ 2 = −θ. Finally, the global EBM analysis in section 2.2, finds the third point φ 3 from the Hamiltonian equation (see Eq. (10)). So consider (u(φ 1 ),u(φ 1 )) and (u(φ 2 ),u(φ 2 )) in the original system, and (u(φ 3 ),u(φ 3 )) in the Hamiltonian and we will have the equation group with the unknown variables (µ, ω 10 , θ). The entire procedure can be seen in the flow chart Figure 3. The efficiency of the approach will be verified in section 3 and 4 in terms of the comparisons with the numerical simulation and Melnikov method.  3. Direct approach for the heteroclinic bifurcation of the third-order system. In order to clearly describe the methodology, we consider a third-order Duffing-van der Pol oscillation system where ε represents the perturbation and µ is the bifurcation parameter. For the heteroclinic bifurcation of the oscillator, we substitute the perturbed hyperbolic function solutions of Eq. (3) into Eq. (24). That produces aω 2 0 tanh ω 10 t − 2aω 2 10 sech 2 (ω 10 t + θ) tanh(ω 10 t + θ) = a 3 α tanh 3 ω 10 t + aεω 10 sech 2 (ω 10 t + θ)(µ − a 2 γ 2,1 tanh 2 ω 10 t). (25) Consider the special points on the orbit Combining Eqs. (27) and (29) together, in Table 1, we find solutions for the bifurcation problem corresponding to different initial values of Eq. (24). For example, substituting the parameter values G 1 and G 2 of Table 1 into Eq. (3), we have the following heteroclinic solutions and G 2 : In Figure 4 (a) and (b), we plot the analytic solutions of Eqs. (30) and (31) against the numerical simulation. In Table 1, we compare the results obtained by our direct approach with the numerical simulation and Melnikov integration under different initial values. It is noted that, in some cases, the parameter resulting from the direct approach may not closely coincide with the numerical method. That requires a higher-order approximation to improve the accuracy of computation.
Apply the result (µ, ω 10 , θ) obtained in section 3 to find the second-order collocation point φ 3,1 in Eq. (3), and then substitute it into the equation group (formed by Eqs. (27) and (29)) to calculate a more accurate group of parameter values. For example, consider the third group of initial values in Table 1, where the saddle points are S 1,2 (±0.894, 0). Substituting the first-order parameter values into Eq. (32), we have a new collocation point at In Figure 5, the analytic solution of Eq. (35) is compared with the result obtained by numerical simulation. In Table 2, we compare the results of G 3 and G 4 obtained by the first and second order direct approach with the numerical simulation. It shows that the direct approach, based on the second-order collocation points, is more accurate than the first-order approximation as compared to the numerical simulation.   In Figure 6, we compare the results of bifurcation parameter µ obtained by different approaches with the parameters of G 3 and the perturbation parameter ε changing from 1 to 10. It shows that the results coincide very well with the numerical simulation even with a large perturbation. That is an improvement as compared with the traditional Melnikov integration, which merely keeps the same value during the variation of parameter ε.

5.
Analysis of the even order homoclinic bifurcation problem. In this paper, we mainly discuss the heteroclinic bifurcation of the system with Z 2 symmetry. According to our research, this new approach, including the hyperbolic solution with phase perturbation, is also available for the homoclinic bifurcation analysis. Considering the following nonlinear system where, in the case of α > 0, Eq. (36) has a homoclinic orbit connected to the saddle point S(0, 0), as seen in Figure 1(a). Based on the hyperbolic function solution u = asech 2 ω 10 t, u = −2aω 10 tanh ω 10 tsech 2 ω 10 t, we give the following homoclinic solution with phase perturbation u = asech 2 ω 10 t, u = −2aω 10 tanh ω 10 tsechω 10 tsech(ω 10 t + θ), where, in the case ε = 0, the amplitude is a = 3ω 2 0 2α and the frequency is ω 10 = ω 0 /2. So the similar unknown variables (µ, ω 10 , θ) can be concerned for ε = 0 in terms of the direct approach through φ i , i = 1, 2, 3.
u & 1 2 Figure 7. Homoclinic orbit of the oscillator; 1-Direct method, 2-Numerical simulation method 6. Conclusions. In this paper, we present a direct approach to the global bifurcation analysis. It combines the perturbed hyperbolic function with the global extension of He's EBM. Thus the unknown parameters in the analytic orbit expression and the bifurcation parameter value give rise to the solutions of a group of equations. Taking the heteroclinic bifurcation case as an example, we consider the perturbed hyperbolic function solution and find the first and second order global collocation points by extending the EBM method. The principal merits of our approach can be summarized as follows: 1. No successive iteration, perturbation or Melnikov type integration are involved in the calculations, which enables this approach to be more easily interpreted by other researchers. 2. Analytical results agree with both the numerical simulation and the Melnikov method, and the accuracy of computation can be improved through the second-order collocation point analysis. In addition, the entire procedure can be programmed by using computational algebra such as Mathematica. 3. It should be a useful approach in performing bifurcation analysis in more complicated ordinary differential equations or the solitary wave research in the partial differential equations [11,12].
The strategy presented can be used in various fields. It shines a light on the global bifurcation analysis and the prediction of the chaotic movement of nonlinear systems, which will be the topic of further research.