Analytical model and energy harvesting analysis of a vibrating slender rod with added tip mass in three-dimensional space

In the paper, a new 3D energy harvesting system is provided. This work discussed the Lagrange approach to derive the differential equations of motion in the case of energy harvesting systems. An electromechanical system consists of a mechanical resonator, a piezoelectric transducer and electrical circuit with the load resistor. A flexible slender rod clamped at the bottom and loaded by the tip mass is proposed as the resonator. Moving in the 3D space, it enables the system to avoid the gravitational potential barrier of the straight vertical shape in case of buckling. This paper investigates the response of the rod deflection and the root mean square power output of selected vibration mode shapes with an attached tip mass.


Introduction
Vibration transducers for scavenging the ambient energy provide an emerging new technology to power small autonomous devices and/or to charge their batteries. The proper modeling of mechanical systems is useful in solutions, especially practical for systems that show nonlinear behaviors. In contrast to linear models where the harvester resonance is gauged to the ambient frequency [1][2][3], the nonlinear devices are characterized by the frequency broadband [4]. It originates from the existence of multiple solutions, and consequently additional resonances of rational or multiple frequencies and also from natural inclinations of the resonance curves in the nonlinear systems [4,5]. In other case [6], the extra degrees-of-freedom system are discussed. Induced by nonlinearities, the discrete real space models were applied to model the system resonator [7]. Such modeling is applied to the whole energy harvesting system including the power electronics in the electrical circuits. Consequently, the authors [8] derived the equations of the coupled tuning and standard resistive circuit systems for piezoelectric energy harvesters using the extended Hamiltonian principle. They analyzed the shunt circuit control responses to tune electromechanical piezoelectric vibration power harvesting structures with proof mass in 2D space. In systems where the influence of random excitation is considered [9,10], the nonlinear equations of motion of the electromechanical system are also in discrete form and describe the system displacement in the xy plane. In the paper [11], the bistable microelectromechanical system was devela e-mail: m.borowiec@pollub.pl (corresponding author) oped using a classical Newton approach. A nonlinear model was prepared with a classical spring damper system. Non-linearity was included by an additive nonlinear term described by a bistable double well potential energy function. Bistable systems of energy harvesting are useful for optimal electrical power production. In the paper [12], the authors presented such a system and discussed it including a bistable potential part in the differential equation of motion in 2D space. In the other paper [13], the authors obtained an analytical model of electromechanical systems by employing the Lagrange's equations. The model was analyzed also is 2D space. The first approaches to the above mathematical modeling were based on systems with linear resonators and piezoelectric (or electromagnetic) transducers [1,3,14]. Unfortunately, linear devices have an insufficient performance in the presence of variable (amplitude and frequency) vibration sources. The recent development on harvesters with nonlinear resonators showed new features for more efficient performance. They are characterized by broader frequency bands because of some additional solutions in nonlinear systems and consequently additional resonances. The frequency band is also broadened due to the inclination of the amplitudefrequency resonance curve [4,15]. The previous studies [16][17][18][19][20] investigated a system with a vertical piezoelectric planar beam. The beam was fixed at the bottom to a frame excited periodically or stochastically. For a large enough tip mass, nonlinear behavior was induced by gravity acting on it. In that case, the beam was clamped at the bottom end to the moving frame and its motion was limited to the plane spanned by frame in vertical direction of excitation. This paper presents the numerical model of a flexible slender rod clamped at the bottom to describe nonlinearities both from a large deflection of the rod and from coupling of mechanical and electrical subsystems. Moreover, despite the system is excited in one direction and the rod keeps symmetry along the vertical axis, its responses are demonstrated in 3D space. In contrast to previous considerations, the motion of the propose system includes rotation along the length of the rod which describes an additional degree of freedom. In the literature are not found mathematical models for such energy harvesting analysis which include all aspects presented in this paper.

Formulation of the problem
This section discussed the Lagrange approach to derive the differential equations of the 3D motion in the case of the flexible slender rod clamp. The assumed length of the rod is much greater than its diameter. Hence, the large deformations lead to geometric non-linearities because of curvature. It creates displacement-strain nonlinear relations. According to [21], the analyzed model of the rod is assumed to be a nonlinear elastic structure because its structure undergoes large displacement but insignificant strains occur. In this paper, the rod is modeled using the Euler-Bernoulli beam theory and the effects of shear deformation is neglected. This simplification is important for the slender rod because it allows a differential rod element to be considered as a rigid body. In the equations of motion, there are 6 degrees of freedom, three translational and three rotational. Moreover, it is easy to determine the displacement of any point of the rod, while the deformation of the neutral axis is known [22]. The slender rod is assumed to be an in-extensionality constraint, thus the number of independent variables is reduced to three. There are two transverse displacements denoted as v, w along the y and z directions and one rotational displacement φ around the x axis. The model of the vertical slender rod with a tip mass is presented in Fig. 1. The kinematic excitation is realized by the harmonic function y e = Q cos(2πf e t) in one horizontal direction y. The piezoelectric patches are connected to the rod surface parallel to the xz plane.
The Lagrange method is applied to derive the equations of motion.
where L = T − Π + W e is the Lagrangian, andq, q are generalized variables of velocities and displacements, respectively.
In the Lagrangian, the kinetic energies T = T 1 + T 2 + T M are: where T 1 and T 2 correspond to the translational and rotational energies of the rod, respectively, but T M corresponds to the translational energy of the tip mass. Note that the rotational energy of the tip is neglected due to the meaningless tip mass dimensions. The potential energies Π = Π 1 + Π 2 of the rod are: where Π 1 corresponds to the elastic strain energy of the rod and Π 2 is the gravity term of potential energy. The constants G and E in Eqs. The angular velocities around the corresponding axes ξ, η and ζ are defined as The curvature components can also be obtained using the definition of curvatures [21]: The work executed by the piezoelectric patches W e in moving and extracting the electrical charge is defined as [17,23]: where D 1 and D 2 are the electromechanical constants that correspond to the linear and nonlinear terms [17]: The piezoelectric patches are considered as a capacitor, and on their electrical side the charge is produced by Introducing Eqs. (2), (3) and (7) to the Lagrangian Eq. (1) and according to the derivations, the equations of motion are obtained and formulated in Eqs. (10)- (12). The equation of the first degree of freedom in the exciting direction is as follows: The equation of the second degree of freedom which is perpendicular to the exciting direction takes the form: Finally, the equation of the third degree of freedom with rotation around the x axis is as follows: The constants N i in the equations of motion are calculated from integers as below: These constants directly depend on the shape functions f ξ , f η and f ζ given in Eq. (14), which describe the mode shapes of the rod at resonance zones. where In this paper, the functions around the axes η and ζ are simultaneously assumed the same for the ith mode f η = f ζ ≡ f i . The successive values of λ i that correspond to the modes of the rod are presented in Table 1.  Figure 2 presents the normalized shape functions given in Eq. (14). In the case of the twisting mode shape, the tip mass does not influence their functions, (the tip mass is considered as dimensionless). Although the tip mass slightly influences the first bending mode shape, this influence on the second mode shape is clearly visible (Fig. 2).
The last equation that corresponds to the electric subsystem is: where R is the load resistance, C p is the capacitance of the piezoelectric patch and U is the voltage on the load resistor connected to the piezoelectric patches. The rod is equipped with the parallelly connected piezoelectric patches. The mechanical and electrical subsystems are coupled by the electromechanical constants D 1 and D 2 (see Table 2). Based on the simulated voltage U , the averaged power scavenged between the times t 1 and t 2 is obtained from: where t 2 − t 1 = Δt is a number of each simulation time period T .

Results and discussion
Section 3 presents the results of the model simulation.
The analysis includes the system responses of a simple rod and a rod with an attached tip mass to increase bending moments. The weight of the tip mass M is comparable to the rod mass m (M = 0.1 kg and m = 0.098 kg). All simulations are provided for the system parameters given in Table 2.
The characteristics reflect the rod behavior of the first two transverse mode shapes of vibration. Figure 3a, b present the results of the amplitude responses A out and root-mean-square P rms achieved by the system vibrating at the first resonance zone for both cases, i.e. without the tip mass (red-blue) and with it (greenblack), respectively. Note that the large deflection visible in Figs. 3 and 4 correspond to the buckle shape of the rod. The excitation corresponds to 0.2 g. The implemented tip mass caused the natural frequency to decrease from 2.8 to 1.2 Hz and the output amplitude of the end rod simultaneously increased almost four times. It generated the seven time increase in the output power which charged the piezoelectric element (Fig. 3b). The first resonance frequency of the rod with the attached tip mass is 5% higher than the frequency estimated by the analytical methods of the linear terms of Eq. (10). In the case without the tip mass, the simulated and analytical natural frequencies are consistent (Fig. 3a, b green-black characteristics). The conclusion is that while the system vibrates in the first mode shape, the nonlinear effect is much more visible when the tip mass M is attached. The tip mass effect the splitting curves of amplitude-frequency characteristic (Fig. 4a, b). Following the initial conditions approach, the excitation frequencies from forward-and backward-sweeps result in the different shape of the amplitude-frequency response in the vicinity of the first resonance zone. In Fig. 4a, the hardening effect of the output amplitude characteristic is confirmed and also clear on the output power characteristic (Fig. 4b).
The opposite effect of tip mass is observed for the second mode shape. Although the natural frequency also decreases, the output amplitude decreases three times (Fig. 5a) and the output power fourteen times (Fig. 5b). This means that energy harvesting efficiency significantly weakens and the nonlinear effect disappears. Reversing the excitation frequency direction does not affect the shape of the amplitude-frequency curve. Comparing the results of the first and the second modes, the attached tip mass enhances the efficiency   Table 3 of energy harvesting in the first case and weakens it in the second one.
To illustrate the influence of the attached tip mass on the dynamics of the system at consecutive resonances, the trajectory of the free end rod for the chosen initial conditions are plotted. The system operates under seven initial conditions: v, w, φ, their time derivatives and also the initial voltage of the piezoelement. They are chosen to emphasize the non-linear nature of the system. Their values are listed in Table 3. In Fig. 6a, b, the trajectories of the first resonance are presented for simple rod and with the attached tip mass, respectively. The system response for both initial conditions is plotted by red and blue lines. The next Fig. 6c-f show the time histories of the displacements v and w to better visualize the system behavior. It is clear that the case of the dynamic responses of the simple rod are convergent for both initial conditions (Fig. 6c-e), while the responses of the rod with tip mass significantly vary, depending on the initial conditions. Such behavior persists for more than a half time of the simulation and,  Table 3 finally, the v displacement converges, but the w visibly diverges (Fig. 6d, f). The results of the second resonance are shown in Fig. 7a-f. The solutions of the v displacement for both initial conditions are similar in this case (Fig. 7c, d) but the w displacement significant changes during the transient time (Fig. 7e, f).
Also, a tiny perturbation in the angle revolution of the rod with the tip mass (Fig. 8b) is visible, while for the simple rod no rotations occur (Fig. 8a). It means that the influence of the tip mass turned out to be important. Figures 9 and 10 depict the phase portraits for the steady state solution of the v direction. For both the simple rod and with the tip mass, the solutions reflect the harmonic system responses for the first and second resonance. The regular behavior is also confirmed by Poincare maps, shown by black dots in both Figs. 9 and 10 for all discussed cases.

Conclusions
The paper presents the analytical model of the rod with the attached tip mass and the piezoelectric transducer. Following the concept of the extension of degree-offreedom, the out-of-plane, (the xy plane Fig. 1) motion of the model was considered. Despite the excitation was in the one direction only, the out-of-plane motion could be induced by selection proper initial conditions. Note that the model includes the geometric nonlinearities resulting from deflections shape of the rod and the final splitting of equilibria (buckling) and the additional nonlinearities resulting from electro-mechanical coupling with the horizontal displacement of the tip mass. The influence of the tip mass on the system dynamics has been studied. It has been confirmed that the tip mass leads to the significant enhancement of the nonlinear effects observed in the dynamical response of the structure. As expected, the rod without the tip mass reaches smaller vibration amplitudes and behaves practically linearly. Consequently, its response does not depend on the initial conditions and out-of-plane vibrations response of the system without the tip mass were not observed. Attaching the tip mass to the rod system and following proper initial conditions, the tendency to out-of-plane responses is revealed. Although the tip mass decreases the first and second natural frequencies of the system response, a significant increase in the vibration amplitude is observed at the first resonance zone. Simultaneously, in the second resonance area vibration amplitude significantly decreases, comparing to the system without the tip mass. Presence of the tip mass changes the effective stiffness of the system and favors the first mode solution with respect to large amplitudes. Additionally, the tip mass leads to a more complex transient response of the system (3D like response), especially visible for the second vibration mode.
Acknowledgements This publication was supported by the program of the Ministry of Science and Higher Education in Poland under project DIALOG 0019/DLG/2019/10 for 2019-2021.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.