1 Introduction

The study on bifurcation behaviors of nonlinear dynamical systems is a research field that has rapidly developed recently, attracting the attention in a variety of research areas^[1-2]. Lyapunov stability and the Floquet theory can provide the stability criterion of single solution near the equilibrium point, while the bifurcation behaviors will not be studied comprehensively and discussed in detail. The singularity theory is regarded as an efficient method to study the bifurcation characteristics of the reduced equations, which can use the unformed and specific methods to deal with a variety of complex bifurcation problems^[3-5]. The classification of bifurcation problems with co-dimension no more than 7 was provided by Keyfitz^[6]. Golubitsky and Guillemin^[7] and Martinet^[8] studied the unfolding of smooth mapping germs based on strong equivalence, providing the versal unfolding theorem of the various forms. Futer et al.^[9] applied techniques from the singularity theory to the multi-parameter bifurcation problem with symmetry on the bifurcation parameter and the state variables. The finite determinacy theorems, the stability problems, and the structural stability of the universal unfolding (UU) are discussed in Ref. [10].

The singularity theory can be divided into three kinds of problems as follows: the recognition problem, the unfolding problem, the co-dimension problem, which were summarized by Golubistky and Schaeffer^[4-5]. Uppal et al.^[11] used the singularity theory to analyze the dynamical behavior of the continuous stirred tank reactor and classify it for variable reactor residence time. Jin and Matsuzaki^[12] used the singularity theory to study the bifurcation behaviors of double pendulum system and obtained the dynamical characteristics of various parameter spaces. The singularity theory was used in the autonomous system of restrained pipe conveying fluid to analyze bifurcations and chaotic motions^[13]. The traditional perturbation theory can only obtain a sort of concrete response curves for the research on periodic solutions of the nonlinear vibration system^[14-15]. Chen and Langford^[16] proposed a new method to study the periodic bifurcation solution to the nonlinear vibration equation with parametric excitation, providing six typical bifurcation response patterns in different parametric regions. The singularity theory can be applied in various nonlinear systems, e.g., the restricted system and the hysteretic system^{[1, 17]}. The UU is a core part in the singularity theory. Researchers in various areas paid widespread attention to studies on the UU^[18-25]. In this paper, we study engineering unfolding (EU) to compare with Chen's previous work^[21] based on the arch structure model.

The research on the UU is complex, and many parameters have no clear physical meanings when we find all the unfolding parameters. The concerned parameters in engineering can be regarded as actual unfolding parameters to study the concerned bifurcation behaviors (the jumping phenomenon, the hysteresis phenomenon, etc.) in engineering, in order to avoid the dangerous parameter region^[26]. The EU parameters consider main parameters in actual engineering systems, e.g., eccentricity and clearance in rotor systems. The EU possesses no universality compared with UU, but EU is enough to discuss the concerned bifurcation behaviors of actual nonlinear dynamical systems^[27-30]. Hou and Chen^[27] used the singular method to study the bifurcation behaviors of a rotor system of aero-engine with the constant maneuver load. Qin and Chen^[28] generalized one state variable singularity to two-state variable singularity and analyzed the bifurcation systems with two parameters. The two-state variable singularity method is applied to actual rotor systems^{[27, 29-30]}. Based on the multi-parameter stability theory^[31], the main EU parameters of nonlinear dynamical systems can be recognized. Meanwhile, the unfolding parameters can be ranked according to what affects the dynamical characteristics of the systems.

The motivation of this paper is to propose the method to select main EU parameters based on the effects of parametric variation of dynamical systems on the eigenvalues of the Frechet matrix. In Section 2, the basic theory of the eigenvalue analysis and the effects of system parameters are introduced. The arch structure model is established in Section 3. The method to select bifurcation is applied to the arch model in Section 4. Section 5 compares EU with UU based on the transient sets and the corresponding bifurcation characteristics. Finally, the conclusions are drawn in Section 6.

2 Basic theory

In this section, we will discuss the eigenvalue analysis in Subsection 2.1, and then the effects of parametric variation on the eigenvalues will be discussed in Subsection 2.2.

2.1 Eigenvalue analysis

Consider the eigenvalue problem,

(1)

where G is an m× m real matrix, and u is the right eigenvector. Therefore, we can get the eigenvalue λ,

(2)

where I is the unit matrix.

The multiplicity of an eigenvalue as a root of characteristic equation is called the algebraic multiplicity k. The maximum number that the linearly independent eigenvectors correspond to the eigenvalues is the geometric multiplicity k_g. As usual, k_g ≤ k. If the algebraic multiplicity is k, and the geometric multiplicity is 1, the eigenvalue is derogatory. If the algebraic multiplicity is 1, the eigenvalue is simple. If the algebraic and geometric multiplicities are both k, the eigenvalue is semi-simple. The simple and semi-simple cases are non-derogatory.

2.1.1 Derogatory eigenvalue case

First, consider a derogatory eigenvalue λ, the number of the linearly independent eigenvectors k, giving Eq. (3),

(3)

where u₀, u₁, …, u_k−1 is the Jordan chain with the length k, u₀ is the proper vector, and u₁, …, u_k-1 is the correlation vector. The formula (3) can be written as a matrix form,

(4)

where

(5)

(6)

Similarly, we can consider the left eigenvector and the Jordan chain problem,

(7)

where v is the left eigenvector of eigenvalue λ, and v^T is the transposition of v.

In the non-derogatory eigenvalue case, the left and right eigenvectors satisfy the orthogonal relation,

(8)

2.1.2 Simple eigenvalue case

If λ is a simple eigenvalue, then the eigenvector u satisfies

(9)

Similarly, the left and right eigenvectors satisfy the orthogonal relationship as follows:

(10)

2.1.3 Semi-simple eigenvalue case

If the eigenvalue is semi-simple and has l-multiplicity, the linearly independent eigenvectors u₁, …, u_l satisfy the following equation:

(11)

which can be written as a matrix form, i.e.,

(12)

where

(13)

(14)

The left and right eigenvectors satisfy the orthogonal relationship as follows:

(15)

where δ_ij is the Kronecker delta symbol.

The Jordan normal form can be calculated for all the eigenvalues,

(16)

(17)

2.2 Effects of dynamical system parameters

In the general autonomous dynamical system, it can be written as a one-order state equation,

(18)

where X is the state vector, G is the Frechet derivative, and F(X) is the nonlinear vector function of X. In the actual problem, the matrix G usually contains parameters. Consider that G is smooth and depends on the parametric vector p=(p₁, …, p_n).

Suppose that λ (p₀) is the simple eigenvalue of G(p₀), calculating the derivation of λ and u to the parameter vector p at p₀. Calculate the derivation of Eq. (9) at both ends. Then, we can get

(19)

where u₀ is the right eigenvector of the eigenvalue λ₀ =λ (p₀). Equation (19) can be transformed as

(20)

Taking the pre-multiplication v₀^T on both sides of Eq. (20) and v₀^TG₀ =v₀^Tλ₀, we can obtain the following equation:

(21)

Therefore, the case of eigenvalue λ varies with parameters of the system. Based on Eq. (20), we can get . G₀ -λ₀I is a singular matrix, and we cannot calculate the inverse matrix directly. Therefore, we import an addition term cu₀. The orthogonal condition satisfies

(22)

Differentiating the system parameters of the formula (22) yields

(23)

Pre-multiply the conjugate vector v₀ and add the formula (20). Then, we can get

(24)

where H₀ is nonsingular,

(25)

Therefore, we can obtain the following equation:

(26)

3 Formulation of arch structure model

In this section, the arch structure model is established. Consider an arch structure in Fig. 1(a). Assume that acting forces and vibration are both in the plane range of the shown figure, and the shaft length is invariable. Taking its micro unit, the forces are shown in Fig. 1(b). The dynamical equation of the y-direction can be established by the Newton's second law as follows:

(27)

where

(28)

(29)

In Fig. 1, q_x and q_y are distributed forces. T_s is the acting force to the arch structure imposed by other objects. K_s is the elastic supporting stiffness offered by other objects around the arch structure. q(x) is the horizontal distribution force. y(x, t) is the deflection curve when the beam of arch structure is on vibration. S is the section area of arch structure. y₀ is the originally static deflection amplitude of the arch beam. u(t) is the deflection amplitude of the arch beam imposed by the dynamical load. u₀₁ is the static deflection amplitude imposed by the static load. ϕ(x) is the vibration mode function matched by the finite element results of the arch structure beam. y₀ is the deflection amplitude of equilibrium point when the beam of arch structure is vibrating. Substitute Eqs. (28) and (29) into Eq. (27), consider the damping force of the beam and the relationship of deformation geometry, multiply ϕ(x) on both sides, and integrate ϕ(x) in the domain (0, l). Then, we also make the assumption that

where q_0x, q_0y, and T_0s are the static load amplitudes, and q_1x, q_1y, and T_1s are the dynamical load amplitudes. Therefore, we can get the static balance equation,

(30)

where W is the mode value calculated by the finite element method in the direction of spring deformation, and EJ is the flexural stiffness of the arch. Therefore, the deflection of vibration balance position of arch structure can be calculated. Meanwhile, we can get the nonlinear vibration equation^[21],

(31)

where is the amplitude caused by the live load, and the parameters D_1i(i=1, …, 7) are determined by the parameters of the arch structure, initial state data, amplitudes of the dead load, the live load, and the vibrational shape function of the arch structure^[21] (which can be obtained by the finite element method), which are expressed in detail in Appendix A.

4 Selection method of bifurcation parameters

Based on the multi-parameter stability theory^[31], we propose a method to find main bifurcation and unfolding parameters of non-autonomous systems. On the basis of the mathematical aspect, the issue to seek parameters is translated into the effect of parameters on the solution structure of the system. When the variation of solution structure is reflected by the eigenvalues of the Frechet derivative matrix, the main parameters can be found by studying the effect of parameter variation on critical eigenvalues of the Frechet derivative matrix.

Applying the averaging method, the first-order approximate solution, in the case of 1/2 sub-harmonic resonance, is obtained as^[21]

(32)

where A is the amplitude of response.

Consider the case of the 1:2 internal resonance . Substituting into Eq. (32) yields

(33)

We can get the Frechet derivative matrix as the form in Eq. (33).

It is evident that (0, 0) is an equilibrium of Eq. (33), and the Frechet derivative matrix is

(34)

We follow the parameters in Ref. [21], i.e., , D₁₄ = γ, and . The eigenvalues of the matrix G are 1 and , and the corresponding right and left eigenvectors are

(35)

According to Eq. (21), the relationship of eigenvalue λ =1 varies with parameters near ,

(36)

Similarly, we can talk about the case ,

(37)

Based on Eqs. (36) and (37), the parameter D₁₄ has the greatest effect on eigenvalues 1 and in comparison with other parameters (D₁₀, ω). Therefore, D₁₄ is chosen as the bifurcation parameter.

5 Bifurcation analysis based on EU and UU

Bifurcation behaviors of EU and UU will be studied in this section, respectively. The geometrical shapes of bifurcation will be listed based on the EU and UU. Then, the comparison of transient sets and bifurcation characters will be discussed.

5.1 Geometrical forms of bifurcation based on EU

We know that the EU considers the actual engineering parameters. Therefore, the EU of a bifurcation equation contains the main dynamical bifurcation behaviors when the original system is subjected to a small perturbation.

In Eq. (3), consider A=x, eliminate the phase angle θ, and x is the state variable. As mentioned in Section 4, λ is the bifurcation parameter (λ =D₁₄). Then, we obtain the following bifurcation equation:

(38)

Let , and c₁ = 64D₁₀² ω² = 8/5. Then, we can get

(39)

where E₁ and E₂ are unfolding parameters, and they are the combination of actual physical parameters.

D₁₃ and D₁₆ are actual physical parameters. Therefore, they are nonnegative. In the physical parameter space, we let ε₁=D₁₃, and ε₂ =D₁₆ for the convenience of expression. Two parameters satisfy

(40)

In the unfolding parameter space spanned by E₁ and E₂, the corresponding parameters satisfy

(41)

It is clear that there is a one-to-one correlation between the physical parameters ε₁ and ε₂ in the unfolding parameter space. Therefore, the unfolding parameter is equivalent to the physical parameter. The dynamical characteristics of arch structure can be reflected by the unfolding parameters. We will use unfolding parameters to analyze singularity of the system.

We know that the formula of the transient set T of a universal unfolding G : R×R×R^k → R is regarded as follows^[4]:

(ⅰ) Bifurcation point set

(ⅱ) Hysteresis set

(ⅲ) Double limit point set

(ⅳ) Transition set

Based on the engineering significance of the state variable x in the original problem, it should be demanded that x≥0. The transition sets and the corresponding bifurcation diagrams are shown in Figs. 2(a) and 2(b), respectively.

5.2 Bifurcation behaviors of UU

The UU of a bifurcation equation can reveal all possible bifurcation behaviors corresponding to the EU. The study on the UU should consider the recognition conditions of strong equivalence of germs^[4]. It is evident that Eq. (38) is strongly equivalent to the following equation (see Ref.[21] in detail):

(42)

Substitute and into Eq. (42), and Eq. (42) can be expressed as

(43)

Similar to Subsection 5.1, the state variable and bifurcation parameter are chosen as x and λ, respectively^[21]. Take the germ equation, i.e.,

(44)

The research process of UU of the germ g can be summarized briefly as follows:

(ⅰ) Prove the restricted tangent space form of the germ g.

(ⅱ) Prove the co-dimension of the germ g.

For further details, the proof process of UU is in Ref. [21]. The co-dimension of the germ g is 5, and the UU of g is given by

(45)

where the parameters α_i (i=1, …, 5) are unfolding parameters.

It is very difficult to discuss the bifurcation behaviors of Eq. (45) in detail, as the UU is a high co-dimension problem with the co-dimension 5. Therefore, we consider three of five parameters as constants when we study the transient sets and bifurcation diagrams. On the basis of transient set formula in Subsection 5.1, the transient sets of the unfolding parameters in the parameter plane and the corresponding bifurcation diagrams in parametric domains are plotted when c=10 and β =3/2.

Here, we choose two groups of transient sets and bifurcation diagrams of ten in Ref.[21] to compare with the EU in Subsection 5.1.

5.3 Comparison of transient sets and bifurcation characteristics

The UU reflects all the unfolding parameters of dynamical systems, and the co-dimension is verified. Therefore, the UU can reflect all the bifurcation behaviors of the system. The unfolding parameters of UU have no clear physical significance, and the parameters contain only mathematical meaning.

The unfolding parameters of EU are equal to the actual physical parameters, which are verified in Subsection 5.1. Therefore, the EU has the actual physical significance.

The transient sets of EU and UU are compared based on Subsections 5.1 and 5.2. In Fig. 2, the transient set is a plane curve related to unfolding parameters E₁ and E₂ (E₁ varies with E₂). Therefore, one single transient set of EU can reveal all the bifurcation characteristics of the concerned parameters in actual engineering. The geometric forms contain the normal shape (see② and ⑥ in Fig. 2), the hysteretic shape (see③ and ④ in Fig. 2), and the fragmentary shape (see⑤ in Fig. 2). The normal shape can be regarded as the stable solution case.

The transient sets can describe the bifurcation variation cases in different regions. For example, the bifurcation is the normal shape in the region ②. The hysteretic shape occurs in the region ③ when the parameter runs out the boundary, and the unstable case occurs. From the region ③ to the region ④, the more complex bifurcation behaviors occur. The bifurcation behavior converts back into the normal shape from the region ⑤ to the region ⑥. The results are very clear in Fig. 2, and the EU can show the transitions of dynamical behaviors, which can be used to optimize the parameters in the nonlinear dynamical model. The parameters in the regions ③ and ④ should be avoided, and it is very important in the theory of nonlinear dynamics.

In comparison with the EU, the results of UU are very complex. The co-dimension of the germ g is 5, the germ equation contains five unfolding parameters in Eq. (45). There are ten groups of transient sets, and two of them are listed in Figs. 3 and 4. The geometric forms of bifurcation contain two normal shapes (⑤ and ⑥ in Fig. 3, ③ and ④ in Fig. 4), the hysteretic shape (①-④ in Fig. 3, ①, ②, and ⑤ in Fig. 4), and the fragmentary shape (⑧ and ⑨ in Fig. 3).

Similar to the EU, the UU can show the transitions of dynamical behaviors. For example, the normal shape turns into the hysteretic shape from ④ to ⑤ in Fig. 3 and ② to ④ in Fig. 4. There are ten groups of transient sets, and we only list two here. Therefore, it is more difficult to show all the bifurcation behaviors of UU than EU.

The concerned bifurcation behaviors are reflected by the analysis of EU, and the dangerous parametric-region should be avoided in engineering design. Therefore, the UU can be replaced by the EU in actual engineering systems when the bifurcation behaviors are analyzed. The EU provides great theoretical guidance for optimum design of parameters in actual engineering systems.

6 Conclusion

A method for seeking main bifurcation parameters of generally nonlinear dynamical systems has been proposed based on the discussion of solution structure of the Frechet matrix. The EU and UU of the arch structure model have been studied via the analysis of bifurcation behaviors. The equivalence of unfolding parameters and physical parameters has been verified. The numerical results have verified that the EU can replace the UU via the comparison of transient sets and bifurcation behaviors, and the EU can be generalized to the more complex engineering system. The EU can reflect the concerned bifurcation behaviors in the engineering system, e.g., the normal shape and the hysteretic shape. In the future, studies on this research will be carried out by the present authors in three aspects as follows: the first is to prove the equivalence of physical parameters and unfolding parameters based on more comprehensive cases (unfolding parameters contain more complex combinations of physical parameters), the second is to study the high-order approximation methods of dynamical equation, and the third is to apply the method of selecting bifurcation parameters to more complex dynamical systems (rotor-bearing systems).

Acknowledgements The authors appreciate the comments of the editors and reviewers and the support of the China Scholarship Council.