Performance enhancement of cantilever piezoelectric energy harvesters by sizing analysis

ABSTRACT

This paper presents the results of the performance of piezoelectric cantilever beams in relation to their size. The total produced power represents the main indicator of performance of a piezoelectric harvesting system while the area of the beams stays constant. Lightweight design is an important aspect in any industry, mainly in the aerospace. In this study, the effects of non-uniformity on the efficiency and power output are studied. Finite element method (FEM) with the application of superconvergent element (SCE) is adopted here to solve the equations. It is observed that the trapezoidal geometry (converging beam) provides a higher output power while the efficiency decreases. Moreover, in order to prove that the power enhancement is achievable while the amount of piezoelectric material consumed is constant the new configuration is proposed. In the configuration, an array of uniform beams connected in series is used instead of one single rectangular beam. The proposed setting generates an output power of 1.817 mW at a resonant frequency of 284.6 Hz when excited by an input acceleration of 1 g. The only challenge is the fundamental frequency difference which is met with the application of proof mass and thinner substrate and piezoelectric layers.

Abbreviations: $a\left(t\right)$: Generalized coordinates vector; $A_p$: The electrode area; A: The amplitude vector of generalized coordinates; $b\left(x\right)$: Width; $b_0$: Initial width; C: Damping matrix; $C_p$: Capacitance of the one piezoelectric layer; ${\bar{C}}_p$: Effective (equivalent) capacitance of the piezoelectric layers; $d_{31}$: Piezoelectric strain coefficient; $D_3$: The electric displacement; $E$: Young’s modulus; $E_3$: The electrical field along the thickness direction; F: Dynamic force vector; $g\left(t\right)$: The translation part of base motion; $G$: Shear modulus; $h_p$: Piezoelectric layers’ thickness; $h_0$: Initial thickness of substrate layer; $i$: Imaginary number; $I\left(t\right)$: Current output; $k_s$: Shear correction factor; K: Global stiffness matrix; $L$: Length of the beam; $l_e$: The length of one element; M: Global mass matrix; $n$: Polynomial’s degree; ${\mathit{N}}_{\mathit{w}},N_{\varphi }$: Matrices of shape functions (1×4 for one element); $P_{in}$: The input mechanical power; $P_{out}$: The output electrical power; $Q\left(t\right)$: Electric charge output; R: Displacement vectors of the generic point S (3×1); $R_l$: External load resistance; $\dot{R}$: Velocity vector of the generic point S (3×1); $t$: Time; $T$: Kinetic energy; $u$: Axial displacement; $U$: Potential energy; $v,w$: Transversal displacements; $V$: Volume; $v\left(t\right)$: The generated piezoelectric voltage$v_0,w_0$: Displacements on the middle-plane; $x,y,z$: Positions of a general point in the relative coordinates system; $X$: The amplitude of the effective displacement, $w_b$; ${\epsilon }_{33}^T$: The dielectric permittivity of piezoelectric layer at constant stress; ${\epsilon }_{33}^s$: The dielectric permittivity of piezoelectric layer at constant strain; $\mathit{\theta }$: The electromechanical coupling vector (3×1 for one element); ${\varphi }_x$: The binding rotation of the cross section; $\gamma$: Transverse shear strain; $\epsilon$: Axial strain; $\sigma$: Axial stress; $\tau$: Transverse shear stress${\psi }_j^e,{\phi }_j^e$: Shape functions of superconvergent element, $j=1,\dots ,4$; $W_i^e,{\varphi }_i^e,W_{\mathrm{i}}^{\prime e},\text{break}{\varphi }_{\mathrm{i}}^{\prime e}\left(t\right)$: Generalized coordinate elements, $i=1,2$; $\rho$: Mass density; $\mathrm{\Gamma }$: The efficiency of harvesting; $\omega$: Excitation frequency; $\mu$: Constant of mass proportionality; Constant of stiffness proportionality; $\zeta$: Damping ratio; ${\alpha }_b$: Width taper ratio; ${\alpha }_h$: Height taper ratio; $\delta$: The symbol of virtual work; $W_{IE}$: The internal electrical energy; $W_{nc}$: The non-conservative mechanical force; RD: Relative difference; PEHs: Piezoelectric energy harvesters; FEM: Finite element method; FRFs: Frequency response functions; DOF: Degree-of–freedom; MEMS: Micro-electromechanical system; MPGs: Micro-power generators; SCE: Superconvergent element; $({)}_p$: Piezoelectric layer properties; $({)}_s$: Substrate layer properties; ${()}^{\prime }$ : Partial differentiation with respect to x; ${()}^{\prime }$: Partial differentiation with respect to t.

GRAPHICAL ABSTRACT

KEYWORDS:

• Piezoelectric energy harvester
• superconvergent element
• tapered beam
• performance enhancement
• sizing

Introduction

Energy harvesting is a type of technology wherein the energy is captured from one or more external energy sources to be stored and eventually used for later applications. The subject has attracted considerable attention during the last decade due to the cost-effectiveness and efficiency of the system. The most important feature of power scavengers is their free energy sources which is collected from the wasted energy of the ambient surroundings. In comparison to the other common kinds of energy providers such as fuels, they are environmentally friendly and cost-effective. This rich source of energy is available in most of the industrial plants, machinery, vehicles of all types, and industrial machinery [1–3].

Nowadays, the energy scavenged from the harvesters have been known as an alternative source of energy. It has been considered as a support of the electric batteries supplying power for a wide variety of devices such as small wireless sensors and portable and wearable electronics owing to the noticeable advances in wireless technology [4,5]. In addition to the weight addition and the high cost of maintenance, short life of batteries represents their most important obstacle compared to the long working life of power scavengers. Although there is an opportunity to replace or recharge the batteries, this type of solution is inefficient and sometimes not feasible because of the poor access to their locations. Hence, local recharging from scavenging the free released energy represents a feasible alternative for solving this problem. As a result, power harvesters have been taken into consideration by many researchers as a renewable power source of portable or wireless devices in different areas [6]. Among the three types of converters including electromagnetic, electrostatic, and piezoelectric, piezoelectric converters are advantageous structures with high electromechanical coupling [7]. As an exclusive advantage, they do not need any external voltage source and are a promising choice for the micro-electromechanical systems (MEMS) application specially wireless sensor node [8]. The amount of generated power from these useful structures strongly depends on the different factors including the applied load, the frequency of vibration, the geometric features, and the boundary conditions [9].

Regarding the development of MEMS power generators, fixed or narrow operation frequency range and low power output are the most important challenges [10]. To address the issue, researchers always have been attempting to reach to higher performance using different practical methods such as modification of shape, size, material properties, and damping [11–15]. There were also other researchers who adopted different strategies to broaden the frequency range and improve the power output [16–22]. Optimization of the active electrode area in order to maximize output is another important strategy illustrated by [23–25]. Another effective adopted method is the application of an array of piezoelectric cantilevers instead of one single beam [26]. It not only covers a higher range of ambient frequencies, but also can help the researchers to improve the output power. To prove the efficiency of the strategy, Soliman et al. [27] presented a new architecture for wideband vibration-based micro-power generators (MPGs) with higher bandwidth compared to a traditional MPG. Also, Liu et al. [28] investigated a MEMS power generator array based on thick-film piezoelectric cantilever to improve both output power and frequency bandwidth. These cantilevers are tunable for the desired low frequencies in the range of 200–400 Hz. As an optimization method, a new architecture consisting of an array of trapezoidal cantilevers on a circular rim was proposed [29]. They could increase the output power and decrease the natural frequency using the diverging beam and tip mass [13,30–32]. The effectiveness of tapered beam configuration can also be seen in other researches [11,33–35]. Changing the thickness of the substrate layer or mounting a tip mass are two practical methods which have been exploited to change the resonance frequency of the cantilever energy harvesters [10]. Other studies were carried out to improve the output power and decrease the fundamental frequency through the application of a parallel or series array of piezoelectric energy harvesters (PEHs) [26,36]. Splitting the piezoelectric layer into smaller segments is another strategy which was implemented to change the resonant frequency of a cantilever energy harvester [37]. In order to cover a wide range of excitation frequencies and attain higher output power, a piezoelectric generator was fabricated by multiple circular diaphragm array [38]. They used four circular diaphragm piezoelectric harvesters with various tip masses which covered the resonance frequencies from 120 Hz to 225 Hz. Additionally, an AIN-based MEMS piezoelectric cantilever array including five piezoelectric cantilever beams and a single proof mass in series connection were presented to increase the output voltage and reduce the loss in a management circuit [39]. Although by the application of series and parallel connections there is possibility to improve the desired electrical output, utilization of more piezoelectric material is unavoidable. As a matter of fact, mass of the piezoelectric part plays an important role in the microfabrication cost of the energy harvesters. In fact, although conventional piezoelectric materials including PZT and PVDF are useful, their fabrication process is time-consuming and expensive. More importantly, they are not environmentally friendly materials. Thus, there is a need to improve the efficiency and electrical output without changing the used amount of piezoelectric material. In this paper, the effects of the tapered structure is investigated on the efficiency and power output of PEHS with the scope to figure out which structure provides the maximum efficiency or power. The beam is modeled with the assumptions of Timoshenko beam theory but it can be easily switched to the Euler–Bernoulli beam with the application of superconvergent element (SCE). Secondly, enhancement power is met with the application of an array of cantilever beams although the same amount of piezoelectric layer is consumed. This new design leads to 15.64% substrate weight reduction, as well. In order to tune the fundamental frequency the thinner substrate and piezoelectric layers are applied. Overall it is shown by sizing a beam, modeled using the SCE model, it is possible to improve the power produced by the structure without alteration of the operation natural frequency of main system and total structure weight.

Tapered energy harvester

In this section, the theoretical model of the tapered energy harvester is presented. The cantilever beam has a non-uniform width, b(x), and height, h(x), varying according to the polynomial function, as shown in Figure 1. It is assumed that perfect bonding is performed such that there is no slip between the piezoceramic layers and the substrate layer.

Figure 1. A 3D schematic model of piezoelectric-based cantilever energy harvester with non-uniform width and height.

Equations of motion for a tapered bimorph PEHs

By the assumption of a Timoshenko beam, the displacement of a generic point can be defined by

(1) $\mathit{R}=\left\{\begin{array}{c}{R}_x\\ R_y\\ R_z\end{array}\right\}=\left\{\begin{array}{c}u\left(x,y,z,t\right)\\ v\left(x,y,z,t\right)\\ w\left(x,y,z,t\right)\end{array}\right\}=\left\{\begin{array}{c}z{\varphi }_x\left(x,t\right)\\ 0\\ w_0\left(x,t\right)+g\left(t\right)\end{array}\right\}$(1)

where $u_0$, $v_0$, and $w_0$ are the displacements projected on the mid-plane, ${\varphi }_x$ is the rotation of the corresponding cross-section at point $x$ and time $t$ relative to the moving base, and $g\left(t\right)$ is the transverse translation coming from the base motion. Using the definition of Lagrange strain tensor and using the non-zero components of the flexible displacement, the normal and shear strains are obtained as

(2) $\begin{array}{l}{\in }_p=\frac{\partial u}{\partial x}=z\varphi {\text{'}}_x\frac{h_s}{2}\le z\le \frac{h_s}{2}+h_p\\ {\in }_{xx}=\frac{\partial u}{\partial x}=z\varphi {\text{'}}_x-\frac{h_s}{2}\le z\le \frac{h_s}{2}-\frac{h_s}{2}-h_p\le z\le -\frac{h_s}{2}\\ {\in }_p=\frac{\partial u}{\partial x}=z\varphi {\text{'}}_x\end{array}$(2)

(3) $\{\begin{array}{l}{\gamma }_p=\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}={\varphi }_x+w{\text{"}}_0\frac{h_s}{2}\le z\le \frac{h_s}{2}+h_p\\ {\gamma }_{xz}=\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}={\varphi }_x+w{\text{'}}_0-\frac{h_s}{2}\le z\le \frac{h_s}{2}\\ {\gamma }_p=\frac{\partial u}{\partial z}+\frac{\partial w}{\partial x}={\varphi }_x+w{\text{'}}_0-\frac{h_s}{2}-h_p\le z\le -\frac{h_s}{2}\end{array}$(3)

With the assumption of linear piezoelectric constitutive equations [40], the calculated piezoelectric constitutive equations are summarized in Equation 4(4) $\left\{\begin{array}{c}{ϵ}_p\\ {\gamma }_p\\ D_3\end{array}\right\}=\left[\begin{array}{ccc}\frac{1}{E_p}& 0& d_{31}\\ 0& \frac{1}{K_s{G}_p}& 0\\ d_{31}& 0& {\epsilon }_{33}^T\end{array}\right]\left\{\begin{array}{c}{\sigma }_p\\ {\tau }_p\\ E_3\end{array}\right\}$(4) [41].

(4) $\left\{\begin{array}{c}{ϵ}_p\\ {\gamma }_p\\ D_3\end{array}\right\}=\left[\begin{array}{ccc}\frac{1}{E_p}& 0& d_{31}\\ 0& \frac{1}{K_s{G}_p}& 0\\ d_{31}& 0& {\epsilon }_{33}^T\end{array}\right]\left\{\begin{array}{c}{\sigma }_p\\ {\tau }_p\\ E_3\end{array}\right\}$(4)

where ${\sigma }_p$ is the normal stress, ${\tau }_p$ is the shear stress, ${ϵ}_p$ is the normal strain, ${\gamma }_p$ is the shear strain, $G_p$ is the shear modulus, $E_p$ is the elastic modulus, $D_3$ is the electric displacement, d₃₁ is the piezoelectric strain coefficient, ${\epsilon }_{33}^T$ is the dielectric permittivity of the piezoelectric layer at constant stress and finally $K_s$ and $E_3$ are the shear correction factor and the applied electrical field, respectively.

The next step through the derivation of governing equations is the application of the Hamilton’s principle (Equation 5(5) $\left\{\begin{array}{c}{ϵ}_p\\ {\gamma }_p\\ D_3\end{array}\right\}=\left[\begin{array}{ccc}\frac{1}{E_p}& 0& d_{31}\\ 0& \frac{1}{K_s{G}_p}& 0\\ d_{31}& 0& {\epsilon }_{33}^T\end{array}\right]\left\{\begin{array}{c}{\sigma }_p\\ {\tau }_p\\ E_3\end{array}\right\}$(5) ).

(5) $\left\{\begin{array}{c}{ϵ}_p\\ {\gamma }_p\\ D_3\end{array}\right\}=\left[\begin{array}{ccc}\frac{1}{E_p}& 0& d_{31}\\ 0& \frac{1}{K_s{G}_p}& 0\\ d_{31}& 0& {\epsilon }_{33}^T\end{array}\right]\left\{\begin{array}{c}{\sigma }_p\\ {\tau }_p\\ E_3\end{array}\right\}$(5)

where $\delta T$, $\delta U$, $\delta W_{IE}$ and $\delta W_{nc}$ are the virtual works of the total kinetic energy, the strain energy, the internal electrical energy and the non-conservative mechanical force, respectively.

The kinetic and strain energies of the entire structure are described using the following equations.

(6) $U=\frac{1}{2}\int \left[{\sigma }_{xx}{ϵ}_{xx}+{\tau }_{xz}{\gamma }_{xz}\right]d{V}_s+\frac{1}{2}\int \left[{\sigma }_p{ϵ}_p+{\tau }_p{\gamma }_p\right]d{V}_p=\frac{1}{2}\underset{0}{\overset{L}{\int }}\left[I_2^E{{\stackrel{.}{\varphi }}_x}^2+I_0^G\left({{\stackrel{.}{w}}_0}^2+{{\varphi }_x}^2+2{\varphi }_x^{.}{\stackrel{.}{w}}_0\right)\right]dx+\frac{1}{2}\underset{0}{\overset{L}{\int }}\left[J_2^E{{\stackrel{.}{\varphi }}_x}^2+J_0^G\left({{\stackrel{.}{w}}_0}^2+{{\varphi }_x^{.}}^2+2{\varphi }_x{\stackrel{.}{w}}_0\right)+\frac{d_{31}}{h_p}{J}_1^{E}v\left(t\right){\stackrel{.}{\varphi }}_x\right]dx$(6)

(7) $\begin{array}{l}T=\frac{1}{2}\int \rho \stackrel{.}{\overrightarrow{R}}\overrightarrow{R}dV+\frac{1}{2}{M}_t{\left[\left({\dot{w}}_0+{\dot{w}}_b\right){|}_{x=l}\right]}^2+\frac{1}{2}{I}_t{\left[{\dot{\varphi }}_x{|}_{x=l}\right]}^2\\ =\frac{1}{2}{\displaystyle \underset{0}{\overset{L}{\int }}\left[I_0^{\rho }{\dot{w}}_{0}2+2I_0^{\rho }{\dot{w}}_0\dot{g}+I_0^{\rho }{\dot{g}}^2+I_2^{\rho }{\dot{\varphi }}_{x}2\right]}dx+\frac{1}{2}{\displaystyle \underset{0}{\overset{\text{L}}{\int }}\left[J_0^{\rho }{\dot{w}}_{0}2+2J_0^{\rho }{\dot{w}}_0\dot{g}+J_0^{\rho }{\dot{g}}^2+J_2^{\rho }{\dot{\varphi }}_{x}2\right]}dx\\ +\frac{1}{2}{M}_t{\left[\left({\dot{w}}_0+{\dot{w}}_b\right){|}_{x=l}\right]}^2+\frac{1}{2}{I}_t{\left[{\dot{\varphi }}_x{|}_{x=l}\right]}^2\end{array}$(7)

where $M_t$ and $I_t$ are the tip mass and its associated moment of inertia. Other coefficients are also defined as follows

(8) $I_i^E\stackrel{\mathrm{\Delta }}{=}\int \int z^i{E}_{s}dydz,I_0^G\stackrel{\mathrm{\Delta }}{=}\int \int K_s{z}^i{G}_{s}dydz,I_i^{\rho }\stackrel{\mathrm{\Delta }}{=}\int \int z^i{\rho }_{s}dydz,i=0,2$(8)

(9) $J_i^E\stackrel{\mathrm{\Delta }}{=}\int \int z^i{E}_{p}dydz,J_0^G\stackrel{\mathrm{\Delta }}{=}\int \int K_s{z}^i{G}_{p}dydz,J_i^{\rho }\stackrel{\mathrm{\Delta }}{=}\int \int z^i{\rho }_{p}dydz,i=0,1,2$(9)

Additionally, the internal electrical energy in the piezoelectric layer can be calculated as follows

(10) $\begin{array}{l}{\text{W}}_{IE}=\frac{1}{2}\int E_3{D}_{3}d{V}_p\\ =-\frac{1}{2}\frac{v\left(t\right)}{h_p}\int E_p{d}_{31}z\varphi {\text{'}}_{xx}-\frac{v\left(t\right)}{h_p}{\epsilon }_{33}^{S}d{V}_{p}c=-\frac{1}{2}\frac{v\left(t\right)}{h_p}{d}_{31}{\displaystyle \underset{0}{\overset{\text{L}}{\int }}{J}_1^E}{\varphi }^{\prime }xdx+\frac{1}{2}{C}_{p}v{\left(t\right)}^2\end{array}$(10)

where $C_p$ is the internal capacitance of the piezoceramic defined by Equation (11)(11) $C_p={\epsilon }_{33}^S\frac{A_p}{h_p}$(11) :

(11) $C_p={\epsilon }_{33}^S\frac{A_p}{h_p}$(11)

where ${\epsilon }_{33}^S$ is the dielectric permittivity at constant strain and $A_p$ is the electrode area.

In the above equations, dots and prime symbols identify derivatives with respect to time and space (i.e. x), respectively. The only non-conservative virtual work is owing to the electric charge output $Q_u\left(t\right)$ which is defined in Equation (12)(12) $\delta W_{nc}=Q_u\left(t\right)\delta v_u\left(t\right)$(12)

(12) $\delta W_{nc}=Q_u\left(t\right)\delta v_u\left(t\right)$(12)

To complete the modeling, the effects of damping should be taken into consideration. In this study, the proportional damping (Rayleigh damping) is applied [42].

Discretization of equations using superconvergent element

The final step is discretization of the coupled model. Finite element method (FEM) with the application of SCE is employed here. This type of element contains two nodes and four degrees-of-freedom (DOF) which include two DOF in the transverse direction and two DOF in the rotational direction (Figure 2). The element provides the precise response with the minimum number of elements since the shape functions are derived by using static equilibrium (Appendix A). One of the most important advantages of this practical element is the simple transferring method. In fact, it is possible to switch readily from the Timoshenko to the Euler–Bernoulli by changing the value of one coefficient ($\beta$). Hence, analysis based on both beam theories can be carried out without performing extra calculations [43]. Using the calculated shape functions available in Appendix A, the displacement parameters in terms of the nodal displacements are represented in Equation (13)(13) $\left[\begin{array}{c}{w}_0\left(x,t\right)\\ {\varphi }_x\left(x,t\right)\end{array}\right]=\left[\begin{array}{c}{N}_w\left(x\right)\\ N_{\varphi }\left(x\right)\end{array}\right]a\left(t\right)=\left[\begin{array}{cccc}{\psi }_1^e\left(x\right)& {\psi }_2^e\left(x\right)& {\psi }_3^e\left(x\right)& {\psi }_4^e\left(x\right)\\ {\phi }_1^e\left(x\right)& {\phi }_2^e\left(x\right)& {\phi }_3^e\left(x\right)& {\phi }_4^e\left(x\right)\end{array}\right]\left[\begin{array}{c}{W}_1^e\left(t\right)\\ {\varphi }_1^e\left(t\right)\\ W_2^e\left(t\right)\\ {\varphi }_2^e\left(t\right)\end{array}\right]$(13) .

(13) $\left[\begin{array}{c}{w}_0\left(x,t\right)\\ {\varphi }_x\left(x,t\right)\end{array}\right]=\left[\begin{array}{c}{N}_w\left(x\right)\\ N_{\varphi }\left(x\right)\end{array}\right]a\left(t\right)=\left[\begin{array}{cccc}{\psi }_1^e\left(x\right)& {\psi }_2^e\left(x\right)& {\psi }_3^e\left(x\right)& {\psi }_4^e\left(x\right)\\ {\phi }_1^e\left(x\right)& {\phi }_2^e\left(x\right)& {\phi }_3^e\left(x\right)& {\phi }_4^e\left(x\right)\end{array}\right]\left[\begin{array}{c}{W}_1^e\left(t\right)\\ {\varphi }_1^e\left(t\right)\\ W_2^e\left(t\right)\\ {\varphi }_2^e\left(t\right)\end{array}\right]$(13)

Figure 2. Two-node finite element with four DOF.

where $N_w\left(x\right)$ and $N_{\varphi }\left(x\right)$ are the shape functions and $a\left(t\right)$ is the unknown generalized coordinates. By substitution of Equation (13)(13) $\left[\begin{array}{c}{w}_0\left(x,t\right)\\ {\varphi }_x\left(x,t\right)\end{array}\right]=\left[\begin{array}{c}{N}_w\left(x\right)\\ N_{\varphi }\left(x\right)\end{array}\right]a\left(t\right)=\left[\begin{array}{cccc}{\psi }_1^e\left(x\right)& {\psi }_2^e\left(x\right)& {\psi }_3^e\left(x\right)& {\psi }_4^e\left(x\right)\\ {\phi }_1^e\left(x\right)& {\phi }_2^e\left(x\right)& {\phi }_3^e\left(x\right)& {\phi }_4^e\left(x\right)\end{array}\right]\left[\begin{array}{c}{W}_1^e\left(t\right)\\ {\varphi }_1^e\left(t\right)\\ W_2^e\left(t\right)\\ {\varphi }_2^e\left(t\right)\end{array}\right]$(13) into Equation (5)(5) $\left\{\begin{array}{c}{ϵ}_p\\ {\gamma }_p\\ D_3\end{array}\right\}=\left[\begin{array}{ccc}\frac{1}{E_p}& 0& d_{31}\\ 0& \frac{1}{K_s{G}_p}& 0\\ d_{31}& 0& {\epsilon }_{33}^T\end{array}\right]\left\{\begin{array}{c}{\sigma }_p\\ {\tau }_p\\ E_3\end{array}\right\}$(5) and by performing integration by parts, the electromechanical coupling equations for a typical finite element result in Equations (14)(14) $M^e\stackrel{\mathrm{..}}{a}\left(t\right)+C^e\stackrel{.}{a}\left(t\right)+K^{e}a\left(t\right)-{\theta }^{e}v\left(t\right)=f^e$(14) and (15(15) $\overline{C_p}{\dot{v}}_u\left(t\right)+\frac{v\left(t\right)}{R_l}+\theta u{e}^T\stackrel{.}{a}\left(t\right)=0$(15) ). Addition of the Rayleigh damping accounting for the mechanical dissipative effects is the last stage of modeling (Equation (16)(16) ${\mathit{C}}^{\mathit{e}}=\mu M^e+\gamma K^e$(16) ).

(14) $M^e\stackrel{\mathrm{..}}{a}\left(t\right)+C^e\stackrel{.}{a}\left(t\right)+K^{e}a\left(t\right)-{\theta }^{e}v\left(t\right)=f^e$(14)

(15) $\overline{C_p}{\dot{v}}_u\left(t\right)+\frac{v\left(t\right)}{R_l}+\theta u{e}^T\stackrel{.}{a}\left(t\right)=0$(15)

(16) ${\mathit{C}}^{\mathit{e}}=\mu M^e+\gamma K^e$(16)

where $\mu$ and $\gamma$ are the constants of mass and stiffness proportionality. Terms ${\mathit{M}}^{\mathit{e}}$, ${\mathit{C}}^{\mathit{e}}$, ${\mathit{K}}^{\mathit{e}},{\theta }^e$, ${\mathit{f}}^{\mathit{e}}$ are, respectively, the mass matrix; the damping matrix; the elastic stiffness matrix; the electromechanical coupling matrix and the vector of dynamical forces corresponding to one element of the piezoelectric cantilever beam, and ${\bar{C}}_p$ is the equivalent capacitance (Appendix B).

Solution of the eletromechanical coupling equations

By the assumption of the input in the harmonic form $g\left(t\right)=X{e}^{i\omega t}$ and assembling the elements of the cantilever energy harvester the components of the forcing vector become:

(17) ${\mathit{f}}_u^e=\left[{\omega }^{2}X\underset{0}{\overset{l_e}{\int }}\left(I_0^{\rho }+J_0^{\rho }\right){N_w}^{T}dx\right]e^{i\omega t}$(17)

The assembled steady-state forms of electromechanical coupling equations can be achieved as Equations (18)(18) $-M{\omega }^2+i\omega C+KA-\theta V=F$(18) and (19(19) $i\omega {\mathit{\theta }}^T\mathit{A}+\left(i\omega \overline{C_p}+\frac{1}{R_l}\right)V=0$(19) ) by the definition of the steady-state response and the voltage output. These two parameters are considered $\mathit{a}\left(t\right)=A{e}^{i\omega t}$ and $v=V{e}^{i\omega t}$, respectively.

(18) $-M{\omega }^2+i\omega C+KA-\theta V=F$(18)

(19) $i\omega {\mathit{\theta }}^T\mathit{A}+\left(i\omega \overline{C_p}+\frac{1}{R_l}\right)V=0$(19)

After solving the above equations the complex voltage is obtained as

(20) $v\left(t\right)=-i\omega {\left(i\omega \overline{C_p}+\frac{1}{R_l}\right)}^{-1}{\mathit{\theta }}^T{\left[-\mathit{M}{\omega }^2+i\omega \mathit{C}+\mathit{K}+i\omega {\left(i\omega \overline{C_p}+\frac{1}{R_l}\right)}^{-1}\theta {\theta }^T\right]}^{-1}\mathit{F}{e}^{i\omega t}$(20)

Expressions of power density and efficiency

To compare the performance of different systems and structures it is necessary to define one type of measurements. Regarding the PEHs, electrical output power extracted from the structure is one kind of calculation which has been widely used. Although this is a practical measurement, it cannot provide a precise and correct evaluation because the value strongly depends on the size, vibration environment and circuit conditions [44]. In fact, the application of a larger beam can result in higher output power. For the comparison of two different configurations the mechanical input power is one important factor which should not be neglected. In this regard, there are other measurement parameters including the power density (output power divided by the piezoelectric layer volume) and efficiency which have been of primary interest. In this paper, both parameters are taken into consideration for the comparison target and their results are analyzed and compared.

Efficiency is a sufficient comparison factor indicating the amount of the input mechanical energy ($P_{in}$) which can be converted into the output electrical power ($P_{out}$) (Equation (21)(21) $\mathrm{\Gamma }=\frac{P_{out}}{P_{in}}$(21) ).

(21) $\mathrm{\Gamma }=\frac{P_{out}}{P_{in}}$(21)

where $P_{in}$ is the energy coming from the shaker which creates the excitation. Thus, the amount of input energy can be written as [45]

(22) $P_{in}={\mathit{f}}^e\dot{\mathit{a}}\left(\mathit{t}\right)$(22)

By having the input and output energies, the harvesting efficiency is defined as

(23) $\mathrm{\Gamma }=\left|\frac{\mathit{v}\left(t\right)I\left(t\right)}{{\mathit{f}}^e\dot{\mathit{a}}\left(\mathit{t}\right)}\right|=\left|\frac{{\left[-i\omega {\left(i\omega \overline{C_p}+\frac{1}{R_l}\right)}^{-1}{\theta }^T{\left[-\mathit{M}{\omega }^2+i\omega \mathit{C}+\mathit{K}+i\omega {\left(i\omega \overline{C_p}+\frac{1}{R_l}\right)}^{-1}\theta {\theta }^T\right]}^{-1}\right]}^2{\mathit{F}}^2}{i{\left[-\mathit{M}{\omega }^2+i\omega \mathit{C}+\mathit{K}+i\omega {\left(i\omega \overline{C_p}+\frac{1}{R_l}\right)}^{-1}\theta {\theta }^T\right]}^{-1}{R}_l{\mathit{F}}^2\omega }\right|$(23)

Results and discussion

Before parametric study of the energy harvester is carried out, there is a need to validate the numerical solution. To satisfy the requirement, the output voltage and fundamental frequency of uniform and non-uniform beams are firstly compared with those available in literature. After that, the influence of non-uniformity on both the power and efficiency are studied to figure out which configuration (converging or diverging beams) can result in the maximum desired measurement parameters. Finally, the application of an array of PEHs (Figure 3) instead of having one single cantilever energy harvester is proposed as the array has the capability to harvest more power.

In the first example, the convergence of the SCE is verified by considering a bimorph cantilever beam investigated by Erturk and Inman [41] for both with and without the tip mass. The properties of the investigated bimorph piezoelectric beam which are connected in series are summarized in Table 1.

Table 1. Properties of the uniform bimorph piezoelectric cantilever beam.

Property	Piezoelectric layer	Substrate layer
Material	PZT-5 H	Brass
Density $\rho$ (kg/m^3)	7500	9000
Length (mm)	24.53	24.53
Thickness (mm)	0.265	0.14
Width (mm)	6.4	6.4
Young’s modulus E (GPa)	60.6	105
Piezoelectric strain coefficient d31 (pC/N)	−274	$-$
Dielectric permittivity at constant strain ${\epsilon }_{33}^S$ (nF/m)	25.55	$-$
Shear modulus G (GPa)	23	40

To be able to compare the results, the short-circuit ($R_l=470\mathrm{\Omega }$) and open-circuit ($R_l=995k\mathrm{\Omega }$) conditions are considered in agreement with the ones proposed in the open literature [41]. The fundamental resonance frequency of the PEH based on the Timoshenko beam theory for two types of considerations are provided in Table 2–3.

Table 2. Fundamental resonance frequency of piezoelectric energy harvester without tip mass.

First natural frequency	Experimental (Hz) [41]	SCE (Hz)	RD (%)
Short circuit ($R_l=470\mathrm{\Omega }$)	502.5	504.7	0.44
Open circuit ($R_l=995k\mathrm{\Omega }$)	524.7	526.9	0.42

Table 3. Fundamental resonance frequency of piezoelectric energy harvester with tip mass (0.239 g).

First natural frequency	Experimental (Hz) [41]	SCE (Hz)	RD (%)
Short circuit ($R_l=470\mathrm{\Omega }$)	338.4	338.7	0.089
Open circuit ($R_l=995k\mathrm{\Omega }$)	356.3	355.6	0.2

The excellent convergence rate of SCE is obvious from the results in which using only one SCE the difference between the experimental and theoretical values for both conditions is negligible (less than 0.5%). Hence, the represented SCE has a very good rate of convergence for the study of PEHs.

Figure 3. A 3D simulated structure of an array of diverging bimorph PEHs.

As the second part of validation, the output voltage of one non-uniform PEH is compared with the literature [10]. They assumed that the width of the beam is varying exponentially through the length by $b\left(x\right)=b_0{e}^{-cx}$. The properties of studied PEH are listed in Table 4. The output voltage for five different tapering parameters are calculated and provided in Table 5. The electrical resistance and damping ratio are taken 1000$\mathrm{\Omega }$ and 1%, respectively.

Table 4. Properties of the tapered bimorph piezoelectric cantilever beam.

Property	Piezoelectric layer	Substrate layer
Material	PZT-4	Steel
Density $\rho$ (kg/m^3)	7960	7800
Length (mm)	60	60
Thickness (mm)	0.6	0.7
Initial width $b_0$(mm)	10	10
Young’s modulus E (GPa)	81.3	200
Piezoelectric strain coefficient d31(pC/N)	−123	$-$
Dielectric permittivity at constant strain ${\epsilon }_{33}^S$ (nF/m)	13.059	$-$

Table 5. Validation of piezoelectric energy harvester’s voltage for different tapering values.

c	0	5	10	15	20
Analytical (V/g) [10]	0.6252	0.5479	0.4802	0.4214	0.3704
SCE element (V/g)	0.6252	0.5470	0.4800	0.4225	0.3729
RD (%)	0	0.16	0.042	0.26	0.67

It can be concluded again that SCE has a great rate of convergence so that with the application of only one SCE the relative error between the numerical and analytical results is less than 0.7%.

In the next example, the effects of tapering ratio on the power density and efficiency of the PEH is analyzed. The energy harvester has the same properties as described in the previous example. To optimize the system, the effects of tapering ratio ($\alpha$) and the degree of polynomial function are investigated and depicted in Figures 4–7. The results show that by increasing the tapering ratio (application of diverging beam) the power output increases while the efficiency decreases. Obviously, increasing the content of piezoelectric material enables the increasing pattern of the output power. The input absorbed energy from the environment will yield more power at the cost of reduced efficiency for the diverging beams. In other words, for the diverging beam the amount of input absorbed energy is higher in comparison to converging structure, leading to smaller efficiency. Additionally, by the comparison of graphs associated with ${\alpha }_b$ and ${\alpha }_h$ illustrated in Figures 4–7, it is concluded that the effects of tapering ratio in width direction (${\alpha }_b$) yields more power output and higher efficiency. It is necessary to point out that the variation of efficiency and the power output are investigated here for the fixed value of load resistance providing small values. In order to maximize the parameters the value of external load resistance should be optimized as well [41].

To study the influence of tapered beam on the performance, the variation of power FRF curves for both tapering parameters is shown separately in Figures 8 and 9. It can be found that the peak voltage of power increases dramatically by 2038% while the tapering ratio changes from ${\alpha }_b=-0.6$ to ${\alpha }_b=0.6$. Although the power output escalates, the fundamental natural frequency decreases significantly. As the tapering parameter varies from 0.6 for diverging beam to −0.6 for converging beam, the resonance frequency reduces by 622% from 917.81 Hz to 147.58 Hz. Consequently, we get variation in frequency range and the peak of output power when ${\alpha }_h$ is changing but the changes are not as dramatic as the ${\alpha }_b$ case. The degree of polynomial function is another determinant parameter for the output power as the influence is depicted in Figure 10. The power output is increased by 1744.5% for tapered beam with n = 5 with respect to uniform beam for which ${\alpha }_b={\alpha }_h=0$. The enhancement of electrical output is due to increasing the surface of the structure. In fact, since the piezoelectric layer is the part contributing to the generation of power output, increasing the polynomial degree results into the expansion of the surface of cantilever beam covered by the piezoelectric material and exposed to strain. This will be resulting in higher electrical output. Thus, the above analysis and investigation shows while the power output increases for the diverging beam the frequency reduces.

Figure 4. Power output variation versus polynomial degree under different tapering ratios in width direction (${\alpha }_b$).

Figure 5. Power output variation versus polynomial degree under different tapering ratios in height direction (${\alpha }_h$).

Figure 6. Efficiency variation versus polynomial degree under different tapering ratios in width direction (${\alpha }_b$).

Figure 7. Efficiency variation versus polynomial degree under different tapering ratios in height direction (${\alpha }_h$).

Figure 8. A power FRF curves for quintic tapered beam in width direction.

Figure 9. A power FRF curves for quintic tapered beam in height direction.

Figure 10. A power FRF curves for tapered beam (${\alpha }_b={\alpha }_h=0.5$) of different polynomial degrees.

In the next example, the amount of power output per volume of piezoelectric layers is investigated using the application of an array of PEHs instead of having one single beam. The properties of the PEHs are the same as in previous examples. Since the effects of sizing are important, the ratio of length over width of the beam is kept constant for comparison. Application of smaller beams will result in the increased fundamental frequency. This frequency difference for the beams with the smaller size should be surmounted since the target is maximization of output power without the change of operation frequency. Application of tip mass has been one of the simplest and most useful methods used to tune the frequency. Selection of an appropriate tip mass could change the mass of system to meet the desirable natural frequency. The tip mass attached to the tip of the cantilever PEH such that the center line of the cube is at the tip of the beam [41]. The variation of power for different arrays of beams with the same volume is shown in Figures 11 and 12. By analyzing the results, it can be concluded that the application of smaller cantilever beams and tip masses (for the frequency tuning) provides the possibility to increase the power noticeably. Decreasing the size of the beams resulted in increasing the output power. For example, by the application of 25 beams connecting in series, the corresponding power increases by 958% with respect to that of the uniform beam although the operation frequency remains constant, around 284.6 Hz. Although the smaller size has an increasing effect on the power output, the amount of required tip mass for the beam with this small size is not feasible. In fact, this tiny beam cannot withstand the $M_t=47.7g$. Moreover, the addition of this amount of mass to the structure will result in escalating the total mass. Thus, the modification should be made to decrease the amount of required proof mass or even remove it entirely. To satisfy the requirement, application of thinner beam method is employed.

Reducing the thickness of piezoelectric and substrate layers is the desired method considered here. From Figure 8 it appears that increasing the tapering ratio not only results in the reduction of resonance frequency, but also amplification of the electrical power output. This kind of useful feature is used to remove the tip mass completely. The new structure consists of an array of tapered beams (${\alpha }_b=0.6,{\alpha }_h=0,\mathrm{n}=5$). The thickness of the piezoelectric and substrate layers are considered to be the same, as well. In order to find a suitable value for the thickness, its effect on the natural frequency should be firstly analyzed. This analysis is necessary to figure out for which thickness one can reach the desired resonant frequency of 284.6 Hz. The effect of increasing thickness for the array of beams is shown in Figure 13. The active/passive layer thickness is varied from 0.1 mm to 0.5 mm. As it was expected, the graphs show that there is a linear relationship between the frequency and the thickness. Hence, for every array of the beams with specific size, there is an optimal thickness value. The comparison of power FRFs per unit mass of piezoelectric layers between the original rectangular beam and the modified tapered beams is represented in Figure 14. The results show that the new designed systems are capable of scavenging more power but the fundamental natural frequency is the same for all structures. The overlapping of the curves indicates that the analyzed configurations provide the same amount of power. The produced output power per unit mass is 88.32 mW/kg under 1.0 g which is 29.5% higher with respect to the one provided by a single uniform beam. Although the structure harvests more energy, there is still room to increase the output by the change of length-to-width ratio.

As the last step, the effects of ratio of the length over the width of the beam are investigated to find the optimum value. As is obvious from Figure 15, decreasing the ratio can result in increasing the maximum output power of the system. By reducing the ratio from 12 to 3 there is a possibility to enhance the output power by 267%. It is noticeable the width at the free end of the beam is the limiting factor. In other words, the ratio should be a feasible value in order to achieve an acceptable structure. The comparison of power FRFs per mass of the piezoelectric material between the primary uniform beam and the new designed ones with different lengths over width ratios is calculated and represented in Figure 16. As the completing part, a performance comparison between an array of 14 tapered bimorph cantilever beams with the geometry: ${\alpha }_b=0.6,{\alpha }_h=0,\mathrm{n}=5,l=15\mathrm{m}\mathrm{m},b_0=5\mathrm{m}\mathrm{m},\text{break}{h}_p=h_s=0.077\mathrm{m}\mathrm{m}$ and the primary uniform cantilever beam is analyzed and shown in Figure 17. The curves reveal that the output power dramatically increases from 0.391 mW to 0.940 mW (140% more power). In this enhanced design, not only the amount of piezoelectric martial consumed is constant, but the mass of the substrate layer is reduced by 15.64%.

By and large, it has been shown that by the adoption of sizing analysis and application of an array of beams in series connection instead of one single beam one can noticeably improve the power output of the system. It should be pointed out that the amount of used active material and fundamental natural frequency remain constant for the comparison. To offset the tip mass concern, the thinner layers can be used but the acceptable amount of thickness reduction depends on the fabrication limitation and the stability of the structure. In other words, it is not possible to decrease the thickness more below the limit of buckling or to the level that the beam could wrinkle. Another practical way which can help in the tuning target is the application of substrate layer with lower modulus of elasticity resulting in lower natural frequencies. This is for further analysis which can be considered to achieve more power output with less amount of required tip mass. Paper-based piezoelectric structures are among the potential structures which can be used to satisfy the requirements [46].

Figure 11. Comparison of power FRFs between one uniform beam and an array of multiple uniform beams.

Figure 12. Comparison of power FRFs per piezoelectric mass between different beam sizes.

Figure 13. Fundamental frequency as a function of piezoelectric/substrate layer thickness ($h_p=h_s$).

Figure 14. Comparison of power FRFs per piezoelectric mass between different beam sizes.

Figure 15. Comparison of maximum power output per piezoelectric mass between different beam ratios.

Figure 16. Comparison of power FRFs per piezoelectric mass between different length/width ratios.

Figure 17. Comparison of power FRFs between the primary uniform beam and an array of tapered beams.

Conclusions

In this study, the effects of tapered beam on both the efficiency and the power output of PEHs were investigated. The influence of sizing on the amount of scavenged power was also studied and an enhanced system consisting of an array of uniform beams was proposed. For the analysis of the structure, the superconvergent beam was taken into consideration. To validate the accuracy of the element, the results available in the literature were compared. Application of this proposed element provided the great possibility of switching readily from the Timoshenko to the Euler–Bernoulli by changing the value of one coefficient ($\beta$).

Analysis of the non-uniform beams showed that when the tapering ratio increases (diverging beam) the output power rises while the efficiency decreases, and vice versa. Expanding the degree of polynomial function also intensified the increasing/decreasing pattern. For example, changing the degree from $n=0$ to $n=5$ for the tapered beam with the tapering properties ${\alpha }_b={\alpha }_h=0.5$ resulted in the power output increase by 1644.5%. The mentioned paradox for the trends of efficiency and power is originated from different amounts of input energy. It should be pointed out that changing the geometry of the beam (application of converging or diverging beam) changes the fundamental resonance frequency of the structure, as well. It is feasible to improve the objective target (efficiency or power) by tapering according to a polynomial function. It can be concluded that although application of diverging beam results in the power output enhancement it increases the weight which is an important subject in different industries such as aerospace. In order to save more material the sizing analysis and application of an array of beams was targeted as the second feasible concept.

In order to improve the power output with the same amount of material, a new system consisting of an array of tapered beams (diverging beams) was proposed. In the new architecture although the same amount of piezoelectric material was used the harvested electrical power output increased significantly by 140%. Additionally, the weight of substrate layer is reduced by 15.64% by the application of thinner substrate layer. The results showed that in order to increase the power output it is not always necessary to employ thicker piezoelectric layers. The power output strongly depends on the properties of the system, especially the external load resistance. For the design of the new system, the operation frequency (fundamental natural frequency) was kept constant. To tune the natural frequency of the smaller arrays of tapered beams to the original uniform beam, the thinner layers were employed. Therefore, the new proposed architecture could provide a solution for the challenging task of recharging the batteries in wireless sensor nodes.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

The authors acknowledge the financial contribution of NSERC and Concordia University that made possible this research; Natural Sciences and Engineering Research Council of Canada; Concordia University.

AppendiesAppendix A. Shape functions of SCE

Lagrange linear shape functions

(A1) ${\chi }_1^e=\left(1-\eta \right)$(A1)

(A2) ${\chi }_2^e=\eta$(A2)

Cubic shape functions

(A3) ${\psi }_1^e=\frac{1}{1+\beta }\left[1+\beta \left(1-\eta \right)-\left(3-2\eta \right){\eta }^2\right]$(A3)

(A4) ${\psi }_2^e=\frac{-l_e}{1+\beta }\left[{\left(1-\eta \right)}^2\eta +\frac{\beta }{2}\left(1-\eta \right)\eta \right]$(A4)

(A5) ${\psi }_3^e=\frac{1}{1+\beta }\left[\left(3-2\eta \right){\eta }^2+\beta \eta \right]$(A5)

(A6) ${\psi }_4^e=\frac{l_e}{1+\beta }\left[\left(1-\eta \right){\eta }^2+\frac{\beta }{2}\left(1-\eta \right)\eta \right]$(A6)

Quadratic shape functions

(A7) ${\phi }_1^e=\frac{6}{l_e\left(1+\beta \right)}\left(1-\eta \right)\eta$(A7)

(A8) ${\phi }_2^e=\frac{1}{1+\beta }\left[1+\beta -\left(4+\beta \right)\eta +3{\eta }^2\right]$(A8)

(A9) ${\phi }_3^e=\frac{-6}{l_e\left(1+\beta \right)}\left(1-\eta \right)\eta$(A9)

(A10) ${\phi }_4^e=\frac{1}{1+\beta }\left[\left(\beta -2\right)\eta +3{\eta }^2\right]$(A10)

where

(A11) $\eta \stackrel{\mathrm{\Delta }}{=}\frac{x}{l_e}$(A11)

(A12) $\beta \stackrel{\mathrm{\Delta }}{=}\frac{12\left(I_2^E+J_2^E\right)}{{l_e}^2\left(I_0^G+J_0^G\right)}$(A12)

where $l_e$ is the length of a typical element of the beam. In the previous shape functions, if $\beta =0$ (or in other words $J_0^G$ → ∞ and $I_0^G$ → ∞, i.e. the assumption of no-shear-deformation theory), shape functions corresponding to the Euler–Bernoulli beam are obtained and the interpolating functions are reduced to cubic and quadratic Hermitian polynomials.

Appendix B. Matrices of piezoelectric bimorph beam in series state

(B1) ${\mathit{M}}^e=\underset{0}{\overset{l_e}{\int }}\left(I_0^{\rho }+2J_0^{\rho }\right){N_w}^T{N}_{w}dx+\underset{0}{\overset{l_e}{\int }}\left(I_2^{\rho }+2J_2^{\rho }\right){N_{\varphi }}^T{N}_{\varphi }dx+M_t{N_w}^T\left(L\right)N_w\left(L\right)+I_t{N_{\varphi }}^T\left(L\right)N_{\varphi }\left(L\right)$(B1)

(B2) $\begin{array}{l}{K}^e={\displaystyle \underset{0}{\overset{l_e}{\int }}\left(I_2^E+2J_2^E\right)}N{\varphi }^T{N}^{\prime }\varphi dx+{\displaystyle \underset{0}{\overset{l_e}{\int }}\left(I_0^G+2J_0^G\right)}{N}^{\prime }{w}^T{N}^{\prime }wdx+{\displaystyle \underset{0}{\overset{l_e}{\int }}\left(I_0^G+2J_0^G\right)}{N}_{\varphi }T{N}_{\varphi }dx\\ +{\displaystyle \underset{0}{\overset{l_e}{\int }}\left(I_0^G+2J_0^G\right)}N\varphi T{N}^{\prime }wdx+{\displaystyle \underset{0}{\overset{l_e}{\int }}\left(I_0^G+2J_0^G\right)}{N}^{\prime }{w}^T{N}_{\varphi }dx\end{array}$(B2)

(B3) $f^e=-{\displaystyle \underset{0}{\overset{l_e}{\int }}\left(I_0^{\rho }+2J_0^{\rho }\right)}{N}_{w}T\ddot{g}dx-M_t{N}_{w}T\left(L\right)\ddot{g}$(B3)

(B4) ${\theta }^e=-\underset{0}{\overset{l_e}{\int }}\frac{d_{31}}{h_p}{J}_1^E{N_{\varphi }^{\prime }}^{T}dxand{\bar{C}}_p=\frac{C_p}{2}$(B4)