Abstract
In the present study, the chaotic response of the nonlinear magnetostrictive actuator (GMA) vibration system is investigated. The mathematical model of the nonlinear GMA vibration system is established according to J-A hysteresis nonlinear model, quadratic domain rotation model, Newton’s third law, and principle of GMA structural dynamics by analyzing the working principle of GMA. Then, the Melnikov function method is applied to the threshold condition of the chaotic response of the system to obtain the sense of Smale horseshoe transformation. Furthermore, the mathematical model is solved to investigate the system response to the excitation force and frequency. Accordingly, the corresponding displacement waveform, phase plane trajectory, Poincaré map, and amplitude spectrum are obtained. The experimental simulation is verified using Adams software. The obtained results show that the vibration equation of the nonlinear GMA vibration system has nonlinear and complex motion characteristics with different motion patterns. It is found that the vibration characteristics of the system can be controlled through adjusting the excitation force and frequency.
1. Introduction
Giant magnetostrictive material (GMM) is a functional material that can generate magnetostrictive deformation at the magnetic field. Compared with other smart materials, GMM has a high Curie temperature, high magnetic-mechanical coupling coefficient, fast response, and large magnetostrictive strain [1, 2]. Moreover, the giant magnetostrictive actuator (GMA) is a new generation of actuators that utilizes the magnetostrictive effect for driving the GMM as a driving element to convert electromagnetic energy into mechanical energy. Considering its superior characteristics, GMA is widely applied to drive electrohydraulic servovalves, control the mechanical vibration, operate microelectromechanical systems, and perform the energy harvesting [3–8].
However, since the GMM rod has the problem of hysteresis nonlinear phenomenon, there is a nonlinear correlation between the applied magnetic field and the output strain of the GMA, which seriously destroys the nonlinear stability of the GMA [9–11]. In order to guide the design of the GMA structure, overall performance estimation, control, and applying the GMA in different applications, it is necessary to carry out research on the nonlinear stability of the GMA. To this end, researchers have proposed diverse schemes, including the robust control [12], adaptive control algorithm [13], optimal control [14], and H∞ control [15], to reduce the influence of nonlinear factors on output stability.
Chaotic motion is a special form of motion in nonlinear systems. When the GMA is applied in the working process of a nonlinear chaotic system, the long-term prediction and control of GMA systems becomes almost impossible, which significantly hinders the application of the GMA in many fields. Therefore, studying the chaotic characteristics of nonlinear vibration systems is of significant importance in the GMA application.
Haghighi and Markazi [16] used the Melnikov method to study the chaotic behavior of micromechanical resonators under effects of double-sided electrostatic forces and designed a robust adaptive fuzzy control algorithm to suppress the chaotic motion. Luo et al. [17] studied the chaotic behavior and adaptive control of magnetic-field electromechanical transducers and showed these chaotic characteristics in phase diagrams with different fractional values. Moreover, Kouomou and Woafo [18] studied the dynamic characteristics of a nonlinear electrostatic transducer with two outputs and found a chaotic motion in the system. Zeng and Li [19] carried out experiments and found that at low frequencies, excessive excitation current or low system stiffness results in chaos. Zhang et al. [20] proved that electrostatic actuators based on the micro-electromechanical systems (MEMS) have a rich bifurcation and chaotic behavior under the parameter excitation. Zhao et al. [21] analyzed the performed research studies on rare earth giant magnetostrictive actuators based on experiments and theories and found that there is unpredictable chaos in their working process. However, reviewing the literature indicates that effects of various parameters on the chaotic motion characteristics of the GMA system has neither been systematically investigated through theoretical, numerical simulation, and experiments nor has it been studied to control the vibration characteristics of the system by controlling each parameter. In fact, only few research studies have been carried out so far to investigate the chaotic characteristics of the GMA with nonlinear vibration systems. While the chaotic responses of other nonlinear vibration systems are comprehensively investigated in numerous research studies [22–28].
The present study intends to initially establish the mathematical model of the GMA nonlinear vibration system based on the analysis of the working principle of the GMA. Then, the Melnikov method will be applied to study the threshold of the chaos in the system. Moreover, the chaotic motion characteristics of the system will be analyzed through multiple points of view, including the bifurcation diagram, Maximum Lyapunov exponential graph, displacement waveform diagram, phase plane locus map, Poincaré map, and the amplitude spectrum. Finally, the numerical simulation will be performed using the Adams software. It is expected that the obtained results can provide theoretical basis and technical support for the structural stability design of the GMA.
2. Structure and Working Principle of the GMA
Figure 1 schematically illustrates the structure of the GMA. The working principle of the GMA can be described as follows: the output displacement of the GMA is controlled by adjusting the input current of the excitation coil. The bias magnetic field generated by the permanent magnet eliminates the GMA frequency doubling phenomenon. Moreover, the outer casing, disc spring, upper end cover, and the cross-shaped output rod constitute a pretensioning device so that the GMM rod obtains a higher expansion and contraction amount. The outer sleeve of the water-cooling chamber, the phase change material, and the coil bobbin constitute a phase change temperature with a control device. Furthermore, the coil bobbin, outlet pipe, water-cooling cavity, and inlet pipe form a cooling water circulation device. The device for controlling the temperature during the phase change and the cooling water circulation device are combined to bring the internal heat of the GMA to the external environment, where the thermal error output can be effectively suppressed. Furthermore, the upper and lower magnetic rings, the magnetic guide sleeve, the permanent magnet, and the cross-shaped output rod are magnetic conductive materials and form a closed magnetic circuit system with the GMM rod, which can reduce the magnetic leakage and reduce the interference of the magnetic field in the driver on the external equipment.

3. Modelling the GMA Nonlinear Vibration System
During the movement of the GMA nonlinear vibration system, the GMM rod has no displacement at one end, while the other end has displacement x(t), velocity , and acceleration . Figure 2 shows the equivalent mechanical system model of the GMA.

Under the joint action of the excitation magnetic field and the prepressure disc spring, the nonlinear dynamic differential equation of the GMA is established according to the working principle of GMA and the principle of mechanical dynamics.
Output displacement of the GMM rod can be expressed as the follows:where ζ and lM denote the GMM rod strain and the effective length of the GMM rod, respectively.
The GMM rod strain can be calculated from the following expression [29]:where σ, EyH, and λ are the stress of the GMM rod, elastic modulus of the GMM rod, and the magnetostrictive strain of the GMM rod, respectively. Moreover, cD and ρ are the internal damping coefficient and the mass density of the GMM rod, respectively.
Based on the quadratic domain rotation model [30], the correlation between the magnetostrictive strain λ of the GMM rod and the magnetization M can be mathematically expressed as follows:where λs and Ms are the saturation magnetostriction coefficient of the isotropic material and the saturation magnetization, respectively. It should be indicated that the total magnetization M generated by the excitation magnetic field is composed of the reversible magnetization Mrev and the irreversible magnetization Mirr.
According to Newton’s third law, the output force F of the GMM rod is equal but in opposite directions. This can be mathematically expressed in the form below:where Fd, F, and σ0 are the output force of the load, output force of the GMM rod, and prestressed stress, respectively. Moreover, AM denotes the equivalent cross-sectional area of the GMM rod.
In order to perform the force analysis of the GMM rod, the force acting on the GMM rod can be described as follows:where mL, cL, and kL are the equivalent mass of the load, equivalent damping coefficient, and the equivalent stiffness coefficient, respectively. Moreover, b is the third stiffness term coefficient of the disc spring.
Under temperature and quasi-static conditions, the Jiles–Atherton hysteresis nonlinear model can be expressed as follows [31, 32]:
The parameters from equations (6)–(10) are Man = nonhysteresis magnetization; He = effective magnetic field inside the GMM; H0 = Amplitude of the excitation magnetic field; = the equivalent parameter of internal magnetic domain interaction of the GMM rod; a = hysteresis magnetization shape parameter; c = reversible coefficient; χr = the influence factor of the internal magnetic field distribution of the material; ω = the frequency of excitation magnetic field.
When the hysteresis of the GMM rod is not considered, the total magnetization M is equal to the hysteresis magnetization Man:
Tyler Expansion of equation (11) results in the following expression:
Ignoring high-order terms and substituting equation (9), equation (12) can be rewritten as follows:
So,
When equations (1), (2), and (5) are implemented into equation (4), the nonlinear dynamic equation of the GMA can be written as follows:wherewhere mM, cM, and kM are the equivalent mass, equivalent damping coefficient, and equivalent stiffness coefficient of the GMM rod, respectively.
Substituting equations (3), (10), and (14) into equation (15), results in the following expression:wherewhere F0 and Ω denote the excitation force of the GMA system and the excitation frequency of the GMA system, respectively.
4. Analyzing Chaotic Characteristics of the System
In order to study the nonlinear vibration problem of the system, it is necessary to normalize equation (17). To this end, a new variable is defined as the following:
Applying the defined variable into equation (17) yields the following expression:where
The chaotic response analysis of the GMA nonlinear vibration system is studied by using the small parameter ε in the front crown of r and f in (20):
Defining function,
Turning equation (22) into a standard form:
When ε = 0, equation (25) is an undisturbed Hamilton system, which is described as the following:
Equation (25) shows that the equilibrium points are (0, 0), (1, 0), and (−1, 0). It is observed that (0, 0) is the center, and (1, 0) and (−1, 0) are saddle points. Therefore, the system has heteroclinic orbit I (0, −1) and (−1, 0), as shown in Figure 3. Since system (25) is a Hamilton system, the Hamilton function is defined as follows:

Since the heteroclinic orbit I crosses the equilibrium points (1, 0) and (−1, 0), h = 1/4 at this time, the heteroclinic orbit I can be expressed as follows:
Available from equation (27),
Only will be studied below:
Equation (29) is the analytical expression of the isomorphic orbit I. The Melnikov function is described as follows:
Substituting equation (29) into equation (30), the following equation is obtained:where csch is the Hyperbolic Cosecant function and is the critical point at which chaos occurs:
Since, so
Substituting r and f into equation (33) and using the Melnikov method indicate that can be used to obtain the threshold condition for the chaotic threshold in the sense that the perturbed system under consideration may produce Smale horseshoe transformation:
For a sufficiently small ε, if M (τ0) has a simple zero at τ0 ∈ (0, T), the system can occur with heteroclinic orbital bifurcation and meaning that chaotic behavior can occur.
5. Numerical Simulation of System Chaotic Response
Table 1 shows the set values of parameters in equation (17). The response of the GMA nonlinear vibration system is calculated with the excitation force and the excitation frequency change bifurcation diagram to determine the range of values of each parameter when the GMA nonlinear vibration system generates chaos. Moreover, the 4th-order Runge–Kutta method is utilized to draw the displacement waveforms, phase plane trajectories, Poincaré maps, and amplitude spectra of the excitation force and excitation frequency in the GMA nonlinear vibration system. The chaotic response of the system with the excitation and excitation frequency are studied separately.
5.1. Influence of the Excitation Frequency on the Chaotic Response of the System
In order to determine the value of the excitation frequency parameter, according to the experience, first set F0 = 1 kN, excitation frequency Ω = 0∼1 Hz, step size ΔΩ = 0.01, and the initial condition as [0 0 0.05]. Figure 4 shows that the fourth-order Runge–Kutta method is used to find the bifurcation diagram of the response of the GMA nonlinear vibration system with the excitation frequency.

The bifurcation diagram (Figure 4) and the maximum Lyapunov exponent diagram (Figure 5) describe the system path for the excitation frequency Ω changing in the following order: chaotic motion ⟶ period 4 ⟶ chaotic motion ⟶ degenerate into period 3 ⟶ multiple paroxysmal chaos ⟶ chaotic motion ⟶ period 4 ⟶ chaotic motion ⟶ multiple paroxysmal chaos ⟶ degenerate into period 2 ⟶ paroxysmal chaos ⟶ chaotic motion ⟶ degenerate into period 1 motion.

Figure 4 shows that when the excitation frequency Ω is within the range of 0.05 to 0.21, the system is chaotic and the corresponding maximum Lyapunov exponent in Figure 5 is positive. Figure 6 illustrates the distribution of a sample Lyapunov exponent, when the excitation frequency Ω is set to 0.2. On the contrary, when the excitation frequency Ω is within the range of 0.21 and 0.32, the system moves in a period of 4 and the maximum Lyapunov exponent is negative, as shown in Figure 5. Furthermore, when the excitation frequency Ω is between 0.32 and 0.35, transient chaotic motions appear in the system and the maximum Lyapunov index is positive, as shown in Figure 5. When the excitation frequency Ω is between 0.35 and 0.42, the system degenerates into period 3 motion, while the corresponding maximum Lyapunov index is negative, as shown in Figure 5. Since then, the system has multiple paroxysmal chaoses. When the excitation frequency Ω is within the range of 0.42 to 0.48, the system is chaotic and the corresponding maximum Lyapunov index is positive, as shown in Figure 5. Moreover, when the excitation frequency Ω is between 0.48 and 0.52, the system degenerates into a period 4 motion, while the corresponding maximum Lyapunov exponent is negative, as shown in Figure 5. When the excitation frequency Ω is within the range of 0.52 to 0.55, chaotic motions appear in the system and the maximum Lyapunov index is positive, as shown in Figure 5. Meanwhile, the system generates multiple paroxysmal chaoses, and then it is degraded into periodic 2 motions. When the excitation frequency Ω is between 0.55 and 0.60, the system is in the period 2 state and the maximum Lyapunov index is negative. Accordingly, Figure 7 illustrates the distribution of a sample Lyapunov exponent, when the excitation frequency Ω is set to 0.58. It is observed that the phase plane trajectory map is a closed curve. Moreover, displacement waveforms are regular and the cycle is stable. Moreover, it is found that the Poincaré map has two attractive points. The amplitude spectrum has a single frequency component, so the system has a typical two-cycle motion. When the excitation frequency Ω is within the range of 0.60 to 0.72, the system is chaotic and the corresponding maximum Lyapunov index is positive. Figure 8 illustrates the distribution of a sample Lyapunov exponent, when the excitation frequency Ω is set to 0.62. Finally, the excitation frequency Ω is about 0.72, the system undergoes paroxysmal chaoses once again, and then the system eventually degenerates into a steady state motion with a period of 1 along with the inverse bifurcation and echoes the largest Lyapunov index in Figure 5 at the same time.

(a)

(b)

(c)

(d)

(a)

(b)

(c)

(d)

(a)

(b)

(c)

(d)
Figures 6 and 8 shows the displacement waveforms, phase plane trajectories, Poincaré maps, and amplitude spectra of the chaotic motion of the system when the excitation frequencies are Ω = 0.2 and Ω = 0.62, respectively. At this point, it is observed that the displacement waveform diagrams in Figures 6(a) and 8(a) are irregular and the periodic vibration is unstable. The phase plane trajectory diagrams in Figures 6(b) and 8(b) are filled with the phase space and cannot be closed for a long time. It should be indicated that the phase trajectory motion is reciprocating motion and the period is infinitely long. The Poincaré maps in Figures 6(c) and 8(c) are neither a finite point set nor a closed curve. There are many frequency components in the amplitude spectrum response in Figure 6(d), including the frequencies of 0.28 Hz, 0.56 Hz, 0.63 Hz, and 1.2 Hz, which include a large proportion. It should be noted that the GMA nonlinear vibration system can change the line spectrum component of the system response, and the amplitude spectrum is a continuous curve and has wide frequency characteristics in the ranges of 0∼0.33 Hz, 0.41∼0.59 Hz, 0.62∼0.92 Hz, and 1.1∼1.57 Hz. The amplitude spectrum response in Figure 8(d) has multiple frequency components, such as frequencies occupying a large proportion at 0.16 Hz, 0.32 Hz, 0.51 Hz, and 1.14 Hz. It should be noted that the GMA nonlinear vibration system can change the line spectrum component of the system response, and the amplitude spectrum is a continuous curve and has wide frequency characteristics in the ranges of 0∼0.22 Hz, 0.24∼0.37 Hz, 0.42∼0.89 Hz, and 0.99∼1.69 Hz. Meanwhile, the maximum Lyapunov exponent corresponding to Figure 5 is positive, so it is concluded that the system is in a chaotic motion.
5.2. Influence of the Incentive Force on the Chaotic Response of the System
When the excitation frequency is Ω = 0.4 Hz, the excitation force is F0 = 0∼1 kN, the step size is ΔF0 = 1, and the initial condition is [0 0 0.06]. Figure 9 shows that the fourth-order Runge–Kutta method is used to find the bifurcation diagram of the response of the GMA nonlinear vibration system with the excitation force.

The bifurcation diagram (Figure 9) and the maximum Lyapunov exponent diagram (Figure 10) describe the system path to chaos, when the excitation force F0 changes in the following order: periodic motion ⟶ jumping motion ⟶ quasi-periodic motion ⟶ multiple paroxysmal chaos ⟶ chaotic motion ⟶ degenerate into period 3 ⟶ multiple paroxysmal chaos ⟶ chaotic motion ⟶ periodic motion.

When the exciting force F0 in Figure 9 is within the range of 0.06 to 0.30, the system obtains a motion with a steady state period of 1, and the maximum Lyapunov index corresponding to Figure 10 is negative. Figure 11 illustrates that when the exciting force F0 is set to 0.1. The phase plane trajectory of the system shows that the phase plane trajectory is a closed curve, while the corresponding displacement waveform is regular and the period is stable. Moreover, it is observed that the Poincaré map is an isolated point and the amplitude spectrum has a single frequency component. Therefore, the system has a typical single-cycle motion. It should be indicated that when the excitation force F0 exceeds 0.25 kN, the system jumps. On the contrary, when the excitation force F0 is between 0.30 and 0.35, the system jumps from the steady state periodic motion to the quasi-periodic motion, and the maximum Lyapunov index corresponding to Figure 10 is positive. Figure 12 illustrates that when the excitation force F0 is set to 0.3. The phase plane trajectory of the system indicates that it is filled with the phase space, while it is closed. The Poincaré graph is a finite point set forming a semiclosed curve. Moreover, the displacement waveform graph is the superposition of two periodic signals and the amplitude spectrum has some isolated and incommensurable frequency peaks. When the excitation force F0 is within the range of 0.35 to 0.5, the system has a chaotic motion and the maximum Lyapunov index corresponding to Figure 10 is positive. Then, the system generates multiple paroxysmal chaos, which is accompanied by the division and degeneration to the period 3 movement. When the excitation force F0 is between 0.5 and 0.53, the system is in a short period 3 state and the maximum Lyapunov index corresponding to Figure 10 is negative. Furthermore, when the excitation force F0 is between 0.53 and 0.76, a chaotic motion appears in the system and the maximum Lyapunov index corresponding to Figure 10 is positive. Figure 13 illustrates that the excitation factor F0is set to 0.6. Finally, when the excitation force F0 is between 0.76 and 1.00, the system returns to the steady state period of 1 and echoes the largest Lyapunov index in Figure 10 at the same time.

(a)

(b)

(c)

(d)

(a)

(b)

(c)

(d)

(a)

(b)

(c)

(d)
Figure 13 shows the displacement waveform, phase plane trajectory, Poincaré map, and amplitude spectrum of the chaotic motion of the system when the excitation force is F0 = 0.6. At this point, the displacement waveform diagram in Figure 13(a) shows an irregular reciprocating cycle. The phase plane trajectory in Figure 13(b) has a chaotic multiturn phase set curve. Moreover, the Poincaré map in Figure 13(c) is linearly distributed from infinite points. The amplitude spectrum response in Figure 13(d) has multiple frequency components, such as the frequencies of 0.27 Hz, 0.55 Hz, and 1.2 Hz, which has a large proportion. It should be indicated that the GMA nonlinear vibration system can change the line spectrum component of the system response, and the amplitude spectrum has a continuous curve and has wide frequency characteristics in the range of 0∼0.36 Hz, 0.4∼0.82 Hz, and 1.0∼1.63 Hz. Meanwhile, the maximum Lyapunov exponent corresponding to Figure 10 is positive, so it is concluded that the system is in a chaotic motion.
6. Adams Simulation Experiment
In this section, the GMA system is simulated with the Adams software. Moreover, Figure 14 shows that a single-degree-of-freedom mass-linear spring system test model is established under Adams/View.

It should be indicated that the mass is set to 1 kg, the mass and the base are constrained with a straight line, and the freedom is retained in the vertical direction. Moreover, the excitation force is set with the Adams Force, acting at the center of the mass and the editing force function is 600 cos (0.4 t). The upper and lower ends of the nonlinear spring act on the centroid of the mass and the base. Moreover, the initial offset of the nonlinear spring is zero. Figure 15 shows that the Dynamic simulation of the GMA system is conducted; the mass displacement versus time is plotted and it is compared with the 4th-order Runge–Kutta method.

Figure 15 shows that the displacement-time history curve of the mass obtained by the Adams simulation test has no steady state. Moreover, it is observed that the vibration is chaotic and long-term unpredictable. Figure 16 shows that the data obtained from the Adams simulation is extracted, subjected to the Fourier transform, and converted to the frequency domain to obtain its spectrum.

Figure 16 shows that the frequencies are concentrated at 0.16 Hz, 0.32 Hz, 0.47 Hz, and 1.15 Hz, which are the main frequency components. It should be indicated that these frequencies are the results of the system natural frequency and subharmonic response. Moreover, it is observed that the spectrum curves in the range of 0∼0.18 Hz, 0.26∼0.38 Hz, 0.42∼0.58 Hz, 0.59∼0.67 Hz, and 1.10∼1.30 Hz are almost continuous with obvious broadband characteristics. This shows that there are many other frequency components in the response of the simulation.
By comparing the spectrogram of Figures 6(d), 8(d), and 13(d) and the spectrogram of Figure 16, it is observed that the frequency components of the Adams simulation without Runge–Kutta numerical solution are complex and the displacement-time history curves of the two results illustrated in Figure 15 are not completely consistent. This is because Adams’s simulation of the nonlinearity, in other words the nonlinear part, is linearly equivalent near its equilibrium position and the stiffness does not change continuously at any time. However, the results of the Adams simulation can still illustrate the broadband response characteristics of the single frequency input multifrequency output of the GMA nonlinear vibration system.
7. Conclusion
Main achievements and conclusions of the present study can be briefly expressed as follows:(1)The chaotic response of the GMA nonlinear vibration system is comprehensively investigated, and the mathematical model of the GMA nonlinear vibration system is established. Then, the Melnikov function method is applied to obtain the threshold conditions for the chaotic response of the system in the sense of Smale horseshoe transformation. Numerical simulations are applied to obtain the vibration response of the system. To this end, ADAMS software is employed to verify the simulation results.(2)It is concluded that changing the excitation frequency and the excitation force parameters can induce the chaotic motion of the system or even avoid the chaotic motion. Moreover, it can control the vibration characteristics of the system.(3)In the present study, the theoretical basis of the chaotic characteristics of the GMA nonlinear vibration system is established, which is beneficial to the further popularization and application of GMM materials and is helpful to guide the avoidance of harmful chaotic vibrations in the GMA applications.
Data Availability
No data were used to support this study.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.