Analysis of Free Pendulum Vibration Absorber Using Flexible Multi-Body Dynamics

Structures which are commonly used in our infrastructures are becoming lighter with progress inmaterial science.These structures due to their light weight and low stiffness have shown potential problem of wind-induced vibrations, a direct outcome of which is fatigue failure. In particular, if the structure is long and flexible, failure by fatigue will be inevitable if not designed properly. The main objective of this paper is to perform theoretical analysis for a novel free pendulum device as a passive vibration absorber. In this paper, the beam-tipmass-free pendulum structure is treated as a flexiblemultibody dynamic system and theANCF formulation is used to demonstrate the coupled nonlinear dynamics of a large deflection of a beam with an appendage consisting of a mass-ball system. It is also aimed at showing the complete energy transfer between twomodes occurringwhen the beam frequency is twice the ball frequency, which is known as autoparametric vibration absorption. Results are discussed and compared with findings of MSC ADAMS. This novel free pendulum device is practical and feasible passive vibration absorber in the mitigation of large amplitude wind-induced vibrations in traffic signal structures.


Introduction
Many mechanical systems can be modeled as a beam with a lumped mass, such as a wing of an airplane with a mounted engine, a robot arm carrying a welding tool, or a traffic light.Understanding the dynamics of those systems having flexible and slender beams is of great importance in vibration analyses to prevent catastrophic failures of the structures.Therefore, there is an extensive amount of experimental and numerical work on the responses of beams in the nonlinear dynamics and vibration field.
There is widespread interest in pendulum modeling and the use of the pendulum as a vibration absorber.This interest ranges from the dynamics of Josephson's Junction in solid state physics [1] to the rolling motion of ships [2,3] and the rocking motion of buildings and structures under earthquakes [4].
The first studies in multibody systems were on the dynamics of the rigid bodies which were related to gyrodynamics, the mechanism theory, and biomechanics.A good review of this topic is given by Schiehlen [27].One of the first formalisms is given by Hooker and Margulies [28] in which they analyzed the satellites interconnected with spherical joints.Another formulation was published in 1967 by Roberson and Wittenburg [29].Wittenburg [30] wrote the first textbook on multibody dynamics in which he explained rigid body kinematics and dynamics as well as general multibody systems.In 1988, Nikravesh [31] provided information about the computer-aided analysis of multibody systems in his textbook.Haug [32] provided basic methods of the computer-aided kinematics and dynamics for spatial and planar systems.Many more authors provided textbooks in the field of kinematic and dynamic simulations of multibody systems such as Roberson and Schwertassek [33], Huston [34], and García de Jalón and Bayo [35].
Until now, we discussed papers and textbooks that were related to the multibody systems consisting of rigid bodies.However, in many applications, bodies undergo large deformations, which necessitate the modeling of the flexible bodies.Flexible multibody systems have attracted many researchers and several flexible multibody formulations have been established such as the floating frame of reference method, incremental finite element corotational method, and the large rotation vector method.Agrawal and Shabana [36] proposed the component mode synthesis method in which each elastic component is identified by three sets of modes: rigid body, reference modes, and normal modes.Rigid body modes are used to describe the rigid body translation and large rotations of a body reference, reference modes are used to define a unique displacement field, and the normal modes are used to define the deformation relative to the body reference.An alternative formulation was proposed by Yoo and Haug [37] in which a lumped mass finite element structural analysis formulation is used to generate deformation modes.In the floating reference frame formulation, a mixed set of absolute and local deformation coordinates are used to define the configuration of the deformable body [38][39][40].This method became the most widely used approach due to its straightforward nature.However, the mass matrix, centrifugal, and Coriolis forces appear to be highly nonlinear.The incremental finite element approach uses rotation angles as nodal variables, which lead to linearized kinematic equations.Therefore, models obtained by using incremental finite elements cannot describe the exact rigid body displacements [41].In order to solve this problem in the incremental finite element approach, a different approach called the large rotation vector formulation has been proposed.In this method, finite rotations are employed instead of infinitesimal rotations, which results in an exact modeling of the rigid body displacements [42].
Most of the methods explained above suffer from highly nonlinear terms inside the mass matrix, centrifugal, and Coriolis forces.Therefore, a new approach called the absolute nodal coordinate formulation (ANCF) was proposed for the solution of large deformation problems [40,[43][44][45][46][47][48][49][50].In this formulation, instead of the angle of rotations, absolute slopes are used as nodal variables.
In this paper, the beam-tip mass-ball structure is treated as a flexible multibody dynamic system and the ANCF formulation is used to demonstrate the coupled nonlinear dynamics of a large deflection of a beam with an appendage consisting of a mass-free pendulum system.This novel free pendulum device is practical and feasible passive vibration absorber in the mitigation of large amplitude wind-induced vibrations in traffic signal structures.

Formulation of Equations of Motions for
Flexible Multibody Dynamics 2.1.Displacement Field.In this paper, a planar beam element is used to model flexible beam under investigation.Referring to Figure 1, the global position vector r of an arbitrary point  on the element is defined in terms of the nodal coordinates and the element shape function as [44] where S is the global shape function and e is the vector of element nodal coordinates defined as The elements of the vector of nodal coordinates are defined as [44] where  is the beam element length and  is the axial coordinate that defines the position of an arbitrary point on where (5)

Mass Matrix.
Kinetic energy of the finite element can be written as Substituting (1) into (6) yields where  is the element volume,  is the mass density of the beam element material, and M  is the mass matrix of the element.Using (7), the mass matrix of the element can be calculated as 2.3.Generalized Elastic Forces.In order to develop the equations of motion of the beam element, generalized elastic forces, Q  , corresponding to the beam element have to be defined.If the strain energy of the element is , then the vector of elastic forces is defined as In this paper, continuum mechanics approach is used to derive the strain energy of the element in the absolute nodal coordinate formulation.This approach uses arc length to define the element deformation instead of using the local coordinate system.The strain energy of the element can be written as where   is the strain energy due to the longitudinal deformation and   is the strain energy due to the bending.The longitudinal strain,   , and the bending curvature, , can be defined as [51] Using (1), ( 9), (10), and (11), one can find the vector of elastic forces as where stiffness matrices K  and K  are defined as [51] K  =  ∫  0   S  , More detailed derivation for the elastic forces is well explained in [51].

Generalized Gravity Forces.
Let  be the distributed gravity force applied on an arbitrary point on the element.Then, the virtual work done due to this external force can be defined as where r is the virtual change in the position vector of the point of application of the force.Using (1), one can write the virtual work term as where Q  = S   is the vector of generalized forces associated with the element nodal coordinates.Using (15), virtual work due to the distributed gravity force of the beam element for the planar case can be written as Therefore, using ( 16), the vector of generalized distributed gravity forces can be written as be the set of generalized coordinates of the flexible body, where  is the number of coordinates.If  ℎ is the number of the constraints, where  ℎ ≤ , then the vector of constraint equations can be written in the form of and the Jacobian of the constraint equations can be defined as Using the vector of Lagrange Multipliers, , one can write generalized constraint forces, Q  , as [40] where C q (q, ) is given by (19) and  is defined as Input data

Kinematic constraints
Boundary conditions

Rigid body data
Use Newton-Raphson method to solve the constraint equation ( 18) for the dependent coordinates Using Newmark method, find the displacement and velocity vectors of the independent coordinates in terms of the unknown acceleration vector Apply Wehage coordinate partitioning technique to find the vector of velocities associated with the dependent coordinates in terms of the unknown acceleration vector Combine the equations of motion of the system given in equation ( 22) with the displacement and velocity vectors found in previous steps and solve them using the Newton-Raphson method for the unknown accelerations and Lagrange Multipliers Using Newmark method, find the displacement and velocity vectors of the independent coordinates Using Wehage coordinate partitioning technique, find the vector of velocities associated with the dependent coordinates Initial conditions Is the simulation time reached?

Stop
No Yes

Equations of Motion.
Using the principles of virtual work in dynamics and the expression of the kinetic and strain energies given in ( 6) and (10), the system equations of motion in augmented form can be written as [40] [ where M is the constant symmetric mass matrix,  is the vector of Lagrange Multipliers, Q  is the vector of body elastic forces, Q  is the vector of the externally applied body forces such as gravity, magnetic, and other forces, and the vector Q  is given by [40] where subscript  denotes partial differentiation with respect to time.

Modeling Beam-Tip Mass-Free Pendulum System
3.1.Coordinate Systems.Referring to Figure 2, the following coordinate systems can be defined for the beam-tip mass-free pendulum system: (i) -: inertial coordinate system (ii)  1 - 1 : body coordinate system of the tip mass, where origin is rigidly attached to the center of mass of the tip mass (iii)  2 - 2 : ball coordinate system, where origin attached to the end of the beam and the coordinate system rotates with the motion of the ball (iv)  3 - 3 : body coordinate system of the ball, where origin is rigidly attached to the center of mass of the ball.

Generalized Coordinates.
The beam-tip mass-free pendulum system consists of three bodies, among which the beam is assumed to be flexible, and the tip mass and ball are assumed to be rigid.The ANCF beam is modeled using three finite elements.Referring to Figure 3, the ANCF beam element has four nodes; each node has four degrees of freedom.Therefore, the total degree of freedom of the beam is 16.The vector of system generalized coordinates can be defined as

Kinematic Constraints.
The connection between the free end of the beam and the tip mass is modeled using a fixed joint.Referring to Figure 4, the following constraint equations for the fixed joint between the two bodies can be written:  F( ball ) Figure 16: FFT of the system for the forcing frequency of 4.13 Hz.
where R  is the absolute displacement vector of the center of mass of the tip appendage, r  is the absolute displacement vector of the beam free end, A  is the transformation matrix between the coordinates - and  1 - 1 , and u  is the displacement vector of point A relative to point B in  1 - 1 coordinate system.Constraint equations between the ball and the tip mass can be defined such that the velocity of the contact point C on the ball has to be equal to the velocity of the contact point C on the tip mass.Therefore, referring to Figure 5, one can write the following equations: or where ṙ  is the absolute velocity vector of the beam free end A, Ṙ is the absolute velocity vector of the center of mass of the ball, u 1 is the displacement vector of point C relative to point A in  2 - 2 coordinate system, u 2 is the displacement vector of point C relative to point D in  3 - 3 coordinate system, Ȧ1 is the time derivative of the transformation matrix A 1 between the coordinates - and  2 - 2 , and Ȧ2 is the time derivative of the transformation matrix A 2 between the coordinates - and  3 - 3 .

Identifying Dependent and Independent Coordinates.
The system has five dependent coordinates and seventeen independent coordinates.For the numerical analysis, the vectors of independent and dependent coordinates are selected as 3.5.Impact Force.Referring to Figure 6, if | − | ≥ 65 ∘ , the system will have impact between the tip mass and the ball.The average contact force,   av , between the bodies during the collision can be defined as where Δ is the contact time and Δ is the change in momentum of one of the colliding bodies which can be calculated using the conservation of momentum law and the concept of coefficient of restitution.
3.6.System Parameters.Numerical integration parameters, rigid body parameters, and flexible body parameters are given in Table 1.

Numerical Solution
The equations of motion of a multibody system consisting of interconnected rigid and deformable bodies are a combined set of ordinary differential and algebraic equations.These kinds of equation sets are called differential algebraic equations (DAEs) in literature.The solutions to DAEs are not as straightforward as ordinary differential equations.Specialized numerical techniques have been developed for the solution of DAEs.In this paper, the direct integration approach based on the Wehage coordinate partitioning technique [52] and the Newmark [53] and the Newton-Raphson methods is used for the solution of DAEs of the multibody system under investigation.The computational algorithm for the dynamic analysis of the multibody system is given in Figure 7.

ANCF Results
. System parameters given in Table 1 are used for the numerical solution, and detailed system   dynamics, including frequency response curves, time history curves, FFT curves, phase plane curves, and energy curves, are plotted for various base excitation frequencies.For each of the numerical analyses, the frequency of the beam-tip mass system is set to 4.13 Hz, and in order to maintain the condition of autoparametric interaction the frequency of the ball is tuned to one half of the beam-tip mass system frequency (i.e., it is set to 2.065 Hz).Before having detailed discussions on the unlocked ball cases, one can refer to Figure 8, which shows the system dynamics when the ball is locked.For the passive absorber case, the ball is locked inside the housing track and the system is excited at its natural frequency.Figure 9 shows the phase plane curve of the system, where a periodic response can be observed.
Figure 10 shows the frequency response ( in meters) curves of the beam and the ball ( in degrees) for the forcing amplitude of 1 mm peak-to-peak.Points A-B on the beam response curve and points C-D on the ball response curve are important as they define the complete energy exchange region.From Figure 10, it is evident that the amplitude of the beam is decreased and the energy is transferred to the ball when the primary resonance case is reached.Moreover, the second peak observed in the beam response is due to the beating phenomena and will be explained later in this paper.In this paper, the detailed system dynamics for the following important frequency values are given: 3.70 Hz is the frequency before the autoparametric region, 4.13 Hz is the frequency at which complete energy transfer occurs, and 4.25 Hz is the frequency where beating phenomena is observed.Figures 11 and 12 show system (ball is free to move) dynamics for the forcing frequency of 3.70 Hz.The time history curves of the beam and the ball, shown in Figure 11, are examined, in which the beam has oscillation, whereas the ball does not have any significance.The phase plane curves shown in Figure 12 are obtained for the full simulation time.Therefore, curves start from the static equilibrium position.After passing the transient response region, they go to steady state, where darker circular patterns on the curves are obtained.Beam  versus ẏ and ball  versus φ curves prove that the responses are periodic.The noninteractive motion between the beam and the ball can be verified by examining  versus  and  versus  curves in Figure 12.From these curves, it is clear that both the beam and the ball frequency ratios are one.Investigating the FFT curves in Figure 11 reveals that the dominant frequency information is due to the base excitation.In this figure, two peaks are observed, in which one of them is the base excitation frequency, which is dominant, and the other is the beam-mass frequency.Since the ball is not oscillating, its frequency information cannot be observed.
Among the numerical analysis performed for the system, the most important results are given in Figures 13-19 for the forcing frequency of 4.13 Hz, where complete energy transfer from the beam to the ball takes place.From Figure 16, one can see that the beam frequency (4.13 Hz) is twice the ball frequency (2.065 Hz) and the forcing frequency is equal to the beam frequency (4.13 Hz), which is the condition for the complete autoparametric interaction.Coupling between the modes of the beam and the ball can also be verified by looking to the phase plane curves of the system at steady state as shown in Figure 15.The loop shown in this figure is the evidence for one to two frequency relationships between the beam and the ball.For further proof of the energy transfer, one can refer to Figure 19, which shows the kinetic energies of the beam and the ball during the full simulation time.
From this figure, one can see that the beam kinetic energy diminishes as time passes; however, the ball kinetic energy still exists.Detailed system energy curves given in Figures 17 and 18 will be compared with ADAMS results found in Section 5.2 in order to validate the ANCF results.Referring to Figure 20, which compares the transverse displacement of the beam tip for active (ball unlocked) and passive (ball locked) cases, and Figure 21, which compares the rotational angle of the ball for unlocked and locked cases, one can claim the free pendulum as a suitable autoparametric vibration absorber under periodic excitation.
Referring to Figure 22, time history curves of the beam and the ball for the forcing frequency of 4.25 Hz show an interesting phenomenon called beating in vibration literature.Beating can be expressed as a special case in which the amplitude of the vibrations periodically varies when the forcing frequency is very close to the frequency of the free vibrations of the system.Beating is undesirable in absorber systems because energy can be transferred from the secondary system (ball) back to the primary system (beamtip mass).As it can be seen from the FFT of the ball in Figure 22, besides the peak value at 2.065 Hz that corresponds to the ball frequency, peaks with values around 1.9 Hz and The source of this peak can be attributed to the beating phenomenon.Therefore, a beating computation analysis will be performed to figure out and confirm this conclusion.

Shock and Vibration
Figure 23 shows the time history curve of the ball during a complete beating cycle.Referring to this figure, the beating period can be computed using the start and end time of the cycle as follows: beating period = 26.1 − 22.2 = 3.9 sec.(30) Similarly, the oscillation period can be computed as oscillation period = beating period 8 = 0.4875 sec, (31) where 8 is the number of oscillations during one complete beating cycle.Moreover, equations for the beating period and the oscillation period in terms of the beating frequencies  1 and  2 can be written as [54] beating period = 2  1 −  2 , Combining ( 30), (31), and (32) and solving for the unknown frequencies, one can obtain  1 = 13.82rad/s = 2.2 Hz,  2 = 11.94 rad/s = 1.9 Hz. (33) Therefore, it is confirmed that the additional peaks seen around 1.9 Hz and 2.2 Hz in the ball FFT are due to beating.Similarly, investigating the beam FFT curve in Figure 23, one can see two major peaks corresponding to the excitation frequency (4.25 Hz) and the beam natural frequency (4.13 Hz) and two minor peaks corresponding to integer multiples (×2) of the beating frequencies  1 and  2 .Figure 24 shows impact details of the system at the beating.

ADAMS Results.
Since the system has complete energy transfer at the forcing frequency of 4.13 Hz, ADAMS simulation is performed at this frequency, and the results are given in Figures 25-30.
Comparing the transverse displacement curves of the beam shown in Figures 16 and 27, one can see that both curves start from the static equilibrium position −0.016 m and have their maximum value approximately around 0.003 m before 5 seconds.After reaching the maximum value, they decrease gradually, and toward the end of the simulation they have small fluctuations around the static equilibrium position.Figures 26,27,and 28 show details of the beam kinetic, potential, and strain energies, respectively.Comparing Figure 28 with the beam strain energy curve given in Figure 18, maximum strain energy is observed before 5 seconds and its value is around 60 N⋅mm.Similarly, for potential energy and kinetic energy curves, maximum values are approximately 6 N⋅mm and 3.25 N⋅mm, respectively.In addition to the numerical matches, similarities can be seen between the decreasing and increasing trends of the strain energy, kinetic energy, and potential energy curves of the two models.Similar observations can be seen between the kinetic energy and the potential energy curves of the ball and the tip mass.
In conclusion, the results obtained from ADAMS and ANCF are observed to be in good quantitative and qualitative agreement even though two methods used different solution approaches.

Conclusion
This paper is concerned with the dynamics of a flexible beam with a tip mass-ball arrangement.The system is treated as a flexible multibody system interconnected with joints.The tip mass and the ball are assumed to be rigid, and the beam is treated as a flexible body.Connection between the tip mass and the free end of the beam is modeled using a fixed joint, and the contact between the ball and the tip mass is modeled using the geometry of the bodies.
The absolute nodal coordinate formulation (ANCF) is used to determine the mass matrix, stiffness matrix, and generalized forces of the system.Generalized elastic forces for the flexible beam are found using the continuum mechanics approach.Nonlinear equations of motion of the system are found using the Lagrangian Formulation, in which constraints are treated explicitly as extra equations by using Lagrange Multipliers.The resulting differential algebraic equations are solved using a two-loop sparse matrix numerical integration method, in which the kinematic constraint equations are satisfied at the position, velocity, and acceleration levels.
The detailed system dynamics including frequency response curves, time history curves, FFT curves, phase plane curves, and energy curves are plotted for various base excitation frequencies.Numerical results are compared with the results of previously studied similar systems and a good qualitative agreement is observed.Moreover, the same system with the same parameters is modeled using the mechanical analysis software, ADAMS, and the results are observed to be in good quantitative agreement, although the two  methods use different formulations.Therefore, in view of the numerical results, it is found that the free pendulum can  be considered a suitable autoparametric vibration absorber under periodic excitation.Futuristic structures will be made of materials like fiber reinforced polymers which are much lighter than steel and hence the vibration problem will be more acute.For example, traffic signal light structures, highway signs, and luminaires are observed to vibrate regularly at steady winds of 10 to 30 mph.The amplitude of vibration depends upon the characteristics of the wind like mean speed, mean direction, and gustiness; dynamic characteristics of the structures; and shape and size of the structure.Vortex shedding and buffeting are the two predominant wind-structure interaction phenomena which could cause vibrations in this class of structures, consisting mainly of a vertical pole and horizontal arm and lights or signs attached to the arm.
This study will provide useful information for designing passive vibration control devices and systems in an exposed environment.Also it will provide important information for the design-iteration process leading to an optimum passive vibration absorber for use in the real world.Results obtained from this study will generate knowledge to develop (a) better understanding of the working principles of control systems, (b) design guidelines and standards, and (c) practical approaches for design, fabrication, and field installations.

Figure 7 :
Figure 7: Computational algorithm for dynamic analysis.

Figure 8 :
Figure 8: Ball locked/detailed dynamics of the beam for the forcing frequency of 4.13 Hz.

Figure 9 :
Figure 9: Ball locked/phase plane curves for the forcing frequency of 4.13 Hz at steady state.

Figure 10 :Figure 11 :
Figure 10: Frequency response curves for an excitation amplitude of 1 mm peak-to-peak.

Figure 17 :
Figure 17: Energy curves of the ball and the tip mass for the forcing frequency of 4.13 Hz.

Figure 18 :
Figure 18: Beam energy curves for the forcing frequency of 4.13 Hz.

Figure 19 :Figure 20 :
Figure 19: Comparison of the beam kinetic energy and the ball kinetic energy for the forcing frequency of 4.13 Hz.

Figure 21 :
Figure 21: Comparison of ball rotation angle in active (ball unlocked) and passive (ball locked) cases.

Figure 22 :
Figure 22: FFT of the system for the forcing frequency of 4.25 Hz.

Figure 24 :Figure 25 :
Figure 24: Impact details for the forcing frequency of 4.25 Hz.

Figure 26 :
Figure 26: Beam kinetic energy curve for the forcing frequency of 4.13 Hz.

Figure 27 :Figure 28 :
Figure 27: Beam potential energy delta curve for the forcing frequency of 4.13 Hz.

Figure 29 :
Figure 29: Ball potential energy delta curve for the forcing frequency of 4.13 Hz.

Figure 30 :
Figure 30: Ball kinetic energy curve for the forcing frequency of 4.13 Hz.
1 ,  2 ,  5 , and  6 are the absolute displacement coordinates and  3 ,  4 ,  7 , and  8 are the global slopes of the nodes.